Actual source code: frelax.F
1: !
2: ! Fortran kernels for SOR relaxations
3: !
4: #include include/finclude/petscdef.h
5: !
6: subroutine FortranRelaxAIJForwardZero(n,omega,x,ai,aj, &
7: & adiag,idiag,aa,b)
8: implicit none
9: PetscScalar x(0:*),aa(0:*),b(0:*),idiag(0:*)
10: PetscReal omega
11: PetscInt n,ai(0:*),aj(0:*),adiag(0:*)
13: PetscInt i,j,jstart,jend
14: PetscScalar sum
15: !
16: ! Forward Solve with zero initial guess
17: !
18: x(0) = omega*b(0)*idiag(0)
19: do 20 i=1,n-1
20: jstart = ai(i)
21: jend = adiag(i) - 1
22: sum = b(i)
23: do 30 j=jstart,jend
24: sum = sum - aa(j) * x(aj(j))
25: 30 continue
26: x(i) = omega*sum*idiag(i)
27: 20 continue
28:
29: ! return
30: end
31: !
32: !-------------------------------------------------------------------
33: !
34: subroutine FortranRelaxAIJBackwardZero(n,omega,x,ai,aj, &
35: & adiag,idiag,aa,b)
36: implicit none
37: PetscScalar x(0:*),aa(0:*),b(0:*),idiag(0:*)
38: PetscReal omega
39: PetscInt n,ai(0:*),aj(0:*),adiag(0:*)
41: PetscInt i,j,jstart,jend
42: PetscScalar sum
43: !
44: ! Backward solve with zero initial guess
45: !
46: do 40 i=n-1,0,-1
47: jstart = adiag(i) + 1
48: jend = ai(i+1) - 1
49: sum = b(i)
50: do 50 j=jstart,jend
51: sum = sum - aa(j)* x(aj(j))
52: 50 continue
53: x(i) = omega*sum*idiag(i)
54: 40 continue
55: return
56: end
57:
58: !-------------------------------------------------------------------
59: !
60: subroutine FortranRelaxAIJForward(n,omega,x,ai,aj,adiag,aa,b)
61: implicit none
62: PetscScalar x(0:*),aa(0:*),b(0:*)
63: PetscReal omega
64: PetscInt n,ai(0:*),aj(0:*),adiag(0:*)
66: PetscInt i,j,jstart,jend
67: PetscScalar sum
68: !
69: ! Forward solve with non-zero initial guess
70: !
71: do 40 i=0,n-1
72: sum = b(i)
74: jstart = ai(i)
75: jend = ai(i+1) - 1
76: do 50 j=jstart,jend
77: sum = sum - aa(j)* x(aj(j))
78: 50 continue
79: x(i) = (1.0 - omega)*x(i) + &
80: & omega*(sum + aa(adiag(i))*x(i))/ aa(adiag(i))
81: 40 continue
82: return
83: end
84:
85: !-------------------------------------------------------------------
86: !
87: subroutine FortranRelaxAIJBackward(n,omega,x,ai,aj,adiag,aa,b)
88: implicit none
89: PetscScalar x(0:*),aa(0:*),b(0:*)
90: PetscReal omega
91: PetscInt n,ai(0:*),aj(0:*),adiag(0:*)
93: PetscInt i,j,jstart,jend
94: PetscScalar sum
95: !
96: ! Backward solve with non-zero initial guess
97: !
98: do 40 i=n-1,0,-1
99: sum = b(i)
101: jstart = ai(i)
102: jend = ai(i+1) - 1
103: do 50 j=jstart,jend
104: sum = sum - aa(j)* x(aj(j))
105: 50 continue
106: x(i) = (1.0 - omega)*x(i) + &
107: & omega*(sum + aa(adiag(i))*x(i)) / aa(adiag(i))
108: 40 continue
109: return
110: end
111: