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,            &
  9:      &           adiag,idiag,aa,b)
 10:       implicit none
 11:       PetscScalar      x(0:*),aa(0:*),b(0:*),idiag(0:*)
 12:       PetscReal        omega
 13:       integer          n,ai(0:*),aj(0:*),adiag(0:*)

 15:       integer          i,j,jstart,jend
 16:       PetscScalar      sum
 17: !
 18: !     Forward Solve with zero initial guess
 19: !
 20:       x(0) = omega*b(0)*idiag(0)
 21:       do 20 i=1,n-1
 22:          jstart = ai(i)
 23:          jend   = adiag(i) - 1
 24:          sum    = b(i)
 25:          do 30 j=jstart,jend
 26:             sum  = sum -  aa(j) * x(aj(j))
 27:  30      continue
 28:          x(i) = omega*sum*idiag(i)
 29:  20   continue
 30: 
 31: !     return
 32:       end
 33: !
 34: !-------------------------------------------------------------------
 35: !
 36:       subroutine FortranRelaxAIJBackwardZero(n,omega,x,ai,aj,           &
 37:      &                                  adiag,idiag,aa,b)
 38:       implicit none
 39:       PetscScalar      x(0:*),aa(0:*),b(0:*),idiag(0:*)
 40:       PetscReal        omega
 41:       integer          n,ai(0:*),aj(0:*),adiag(0:*)

 43:       integer          i,j,jstart,jend
 44:       PetscScalar      sum
 45: !
 46: !     Backward solve with zero initial guess
 47: !
 48:       do 40 i=n-1,0,-1
 49:          jstart  = adiag(i) + 1
 50:          jend    = ai(i+1) - 1
 51:          sum     = b(i)
 52:          do 50 j=jstart,jend
 53:             sum = sum - aa(j)* x(aj(j))
 54:  50      continue
 55:          x(i)    = omega*sum*idiag(i)
 56:  40   continue
 57:       return
 58:       end
 59: 
 60: !-------------------------------------------------------------------
 61: !
 62:       subroutine FortranRelaxAIJForward(n,omega,x,ai,aj,adiag,aa,b)
 63:       implicit none
 64:       PetscScalar      x(0:*),aa(0:*),b(0:*)
 65:       PetscReal        omega
 66:       integer          n,ai(0:*),aj(0:*),adiag(0:*)

 68:       integer          i,j,jstart,jend
 69:       PetscScalar      sum
 70: !
 71: !     Forward solve with non-zero initial guess
 72: !
 73:       do 40 i=0,n-1
 74:          sum    = b(i)

 76:          jstart = ai(i)
 77:          jend    = ai(i+1) - 1
 78:          do 50 j=jstart,jend
 79:             sum = sum - aa(j)* x(aj(j))
 80:  50      continue
 81:          x(i)    = (1.0 - omega)*x(i) +                                 &
 82:      &              omega*(sum + aa(adiag(i))*x(i))/ aa(adiag(i))
 83:  40   continue
 84:       return
 85:       end
 86: 
 87: !-------------------------------------------------------------------
 88: !
 89:       subroutine FortranRelaxAIJBackward(n,omega,x,ai,aj,adiag,aa,b)
 90:       implicit none
 91:       PetscScalar      x(0:*),aa(0:*),b(0:*)
 92:       PetscReal        omega
 93:       integer          n,ai(0:*),aj(0:*),adiag(0:*)

 95:       integer          i,j,jstart,jend
 96:       PetscScalar      sum
 97: !
 98: !     Backward solve with non-zero initial guess
 99: !
100:       do 40 i=n-1,0,-1
101:          sum    = b(i)

103:          jstart = ai(i)
104:          jend   = ai(i+1) - 1
105:          do 50 j=jstart,jend
106:             sum = sum - aa(j)* x(aj(j))
107:  50      continue
108:          x(i)    = (1.0 - omega)*x(i) +                                 &
109:      &              omega*(sum + aa(adiag(i))*x(i)) / aa(adiag(i))
110:  40   continue
111:       return
112:       end
113: