Actual source code: frelax.F
1: !
2: ! "$Id: frelax.F,v 1.2 2001/08/21 21:05:00 bsmith Exp $";
3: !
4: ! Fortran kernels for SOR relaxations
5: !
6: #include include/finclude/petscdef.h
7: !
8: subroutine FortranRelaxAIJForwardZero(n,omega,x,ai,aj,adiag,aa,b)
9: implicit none
10: PetscScalar x(0:*),aa(0:*),b(0:*)
11: PetscReal omega
12: integer n,ai(0:*),aj(0:*),adiag(0:*)
14: integer i,j,jstart,jend
15: PetscScalar sum
16: !
17: ! Forward Solve with zero initial guess
18: !
19: x(0) = omega*b(0) / aa(adiag(0))
20: do 20 i=1,n-1
21: jstart = ai(i)
22: jend = adiag(i) - 1
23: sum = b(i)
24: do 30 j=jstart,jend
25: sum = sum - aa(j) * x(aj(j))
26: 30 continue
27: x(i) = omega*sum / aa(adiag(i))
28: 20 continue
29:
30: ! return
31: end
32: !
33: !-------------------------------------------------------------------
34: !
35: subroutine FortranRelaxAIJBackwardZero(n,omega,x,ai,aj,adiag,aa,b)
36: implicit none
37: PetscScalar x(0:*),aa(0:*),b(0:*)
38: PetscReal omega
39: integer n,ai(0:*),aj(0:*),adiag(0:*)
41: integer 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 / aa(adiag(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: integer n,ai(0:*),aj(0:*),adiag(0:*)
66: integer 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: integer n,ai(0:*),aj(0:*),adiag(0:*)
93: integer 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: