Changeset 16 for trunk/NEMO/OPA_SRC/SOL/solsor.F90
- Timestamp:
- 2004-02-17T09:06:15+01:00 (20 years ago)
- File:
-
- 1 edited
-
trunk/NEMO/OPA_SRC/SOL/solsor.F90 (modified) (5 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/OPA_SRC/SOL/solsor.F90
r3 r16 22 22 !! * Routine accessibility 23 23 PUBLIC sol_sor ! ??? 24 24 25 !!---------------------------------------------------------------------- 25 26 !! OPA 9.0 , LODYC-IPSL (2003) … … 28 29 CONTAINS 29 30 30 SUBROUTINE sol_sor( k t, kindic )31 SUBROUTINE sol_sor( kindic ) 31 32 !!---------------------------------------------------------------------- 32 33 !! *** ROUTINE sol_sor *** 33 34 !! 34 35 !! ** Purpose : Solve the ellipic equation for the barotropic stream 35 !! function system ( default option) or the transport divergence36 !! system ( key_dynspg_fsc =T) using a successive-over-relaxation36 !! function system (lk_dynspg_rl=T) or the transport divergence 37 !! system (lk_dynspg_fsc=T) using a successive-over-relaxation 37 38 !! method. 38 39 !! In the former case, the barotropic stream function trend has a … … 59 60 !!---------------------------------------------------------------------- 60 61 !! * Arguments 61 INTEGER, INTENT( in ) :: kt ! ocean time-step62 62 INTEGER, INTENT( inout ) :: kindic ! solver indicator, < 0 if the conver- 63 63 ! ! gence is not reached: the model is … … 67 67 !! * Local declarations 68 68 INTEGER :: ji, jj, jn ! dummy loop indices 69 REAL(wp) :: zgwgt ! temporary scalar70 69 !!---------------------------------------------------------------------- 71 70 72 73 ! Iterative loop 74 ! ============== 75 76 IF( kt == nit000 ) gccd(:,:) = sor * gcp(:,:,2) 71 ! ! ============== 72 DO jn = 1, nmax ! Iterative loop 73 ! ! ============== 77 74 78 79 DO jn = 1, nmax 80 81 !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,, 75 CALL lbc_lnk( gcx, c_solver_pt, 1. ) ! applied the lateral boubary conditions 82 76 83 ! boundary conditions (at each sor iteration) only cyclic b. c. are required 84 #if defined key_dynspg_fsc 85 # if defined key_mpp 86 ! Mpp: export boundary values to neighbouring processors 87 CALL lbc_lnk( gcx, 'S', 1. ) ! S=T with special staff ??? 88 # else 89 CALL lbc_lnk( gcx, 'T', 1. ) 90 # endif 91 #else 92 # if defined key_mpp 93 ! Mpp: export boundary values to neighbouring processors 94 CALL lbc_lnk( gcx, 'G', 1. ) ! G= F with special staff ??? 95 # else 96 CALL lbc_lnk( gcx, 'F', 1. ) 97 # endif 98 #endif 99 100 !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,, 101 102 ! 1. Residus 103 ! ---------- 104 77 ! Residus 78 ! ------- 105 79 DO jj = 2, jpjm1 106 80 DO ji = 2, jpim1 107 gcr(ji,jj) = gcb(ji,jj ) - gcx(ji ,jj ) &108 - gcp(ji,jj,1)*gcx(ji ,jj-1) &109 - gcp(ji,jj,2)*gcx(ji-1,jj ) &110 - gcp(ji,jj,3)*gcx(ji+1,jj ) &111 - gcp(ji,jj,4)*gcx(ji ,jj+1)81 gcr(ji,jj) = gcb(ji,jj ) - gcx(ji ,jj ) & 82 - gcp(ji,jj,1) * gcx(ji ,jj-1) & 83 - gcp(ji,jj,2) * gcx(ji-1,jj ) & 84 - gcp(ji,jj,3) * gcx(ji+1,jj ) & 85 - gcp(ji,jj,4) * gcx(ji ,jj+1) 112 86 END DO 113 87 END DO 88 CALL lbc_lnk( gcr, c_solver_pt, 1. ) ! applied the lateral boubary conditions 114 89 115 !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,, 116 117 ! 1.1 Boundary conditions (at each sor iteration) only cyclic b. c. are required 118 #if defined key_dynspg_fsc 119 # if defined key_mpp 120 ! Mpp: export boundary values to neighbouring processors 121 CALL lbc_lnk( gcr, 'S', 1. ) 122 # else 123 ! mono- or macro-tasking: W-point, >0, 2D array, no slab 124 CALL lbc_lnk( gcr, 'T', 1. ) 125 # endif 126 #else 127 # if defined key_mpp 128 ! Mpp: export boundary values to neighbouring processors 129 CALL lbc_lnk( gcr, 'G', 1. ) 130 # else 131 ! mono- or macro-tasking: W-point, >0, 2D array, no slab 132 CALL lbc_lnk( gcr, 'F', 1. ) 133 # endif 134 #endif 135 136 ! 1.2 Successive over relaxation 137 90 ! Successive over relaxation 138 91 DO jj = 2, jpj 139 92 DO ji = 1, jpi 140 gcr(ji,jj) = gcr(ji,jj) - sor *gcp(ji,jj,1)*gcr(ji,jj-1)93 gcr(ji,jj) = gcr(ji,jj) - sor * gcp(ji,jj,1) * gcr(ji,jj-1) 141 94 END DO 142 95 DO ji = 2, jpi 143 gcr(ji,jj) = gcr(ji,jj) - sor *gcp(ji,jj,2)*gcr(ji-1,jj)96 gcr(ji,jj) = gcr(ji,jj) - sor * gcp(ji,jj,2) * gcr(ji-1,jj) 144 97 END DO 145 98 END DO 146 99 147 !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,148 149 100 ! gcx guess 150 151 101 DO jj = 2, jpjm1 152 102 DO ji = 1, jpi 153 gcx(ji,jj) = ( gcx(ji,jj)+sor*gcr(ji,jj))*bmask(ji,jj)103 gcx(ji,jj) = ( gcx(ji,jj) + sor * gcr(ji,jj) ) * bmask(ji,jj) 154 104 END DO 155 105 END DO 156 157 !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,, 158 159 ! boundary conditions (at each sor iteration) only cyclic b. c. are required 160 #if defined key_dynspg_fsc 161 # if defined key_mpp 162 ! Mpp: export boundary values to neighbouring processors 163 CALL lbc_lnk( gcx, 'S', 1. ) 164 # else 165 ! mono- or macro-tasking: W-point, >0, 2D array, no slab 166 CALL lbc_lnk( gcx, 'T', 1. ) 167 # endif 168 #else 169 # if defined key_mpp 170 ! Mpp: export boundary values to neighbouring processors 171 CALL lbc_lnk( gcx, 'G', 1. ) 172 # else 173 ! mono- or macro-tasking: W-point, >0, 2D array, no slab 174 CALL lbc_lnk( gcx, 'F', 1. ) 175 # endif 176 #endif 177 178 ! maximal residu (old exit test on the maximum value of residus) 179 ! 180 ! imax = isamax( jpi*jpj, gcr, 1 ) 181 182 ! avoid an out of bound in no bounds compilation 183 184 ! iimax1 = mod( imax, jpi ) 185 ! ijmax1 = int( float(imax) / float(jpi)) + 1 186 ! resmax = abs( gcr(iimax1,ijmax1) ) 106 CALL lbc_lnk( gcx, c_solver_pt, 1. ) 187 107 188 108 ! relative precision 189 190 rnorme = 0. 191 DO jj = 1, jpj 192 DO ji = 1, jpi 193 zgwgt = gcdmat(ji,jj) * gcr(ji,jj) 194 rnorme= rnorme + gcr(ji,jj)*zgwgt 195 END DO 196 END DO 197 198 #if defined key_mpp 199 ! mpp sum over all the global domain 200 CALL mpp_sum( rnorme ) 201 #endif 109 rnorme = SUM( gcr(:,:) * gcdmat(:,:) * gcr(:,:) ) 110 IF( lk_mpp ) CALL mpp_sum( rnorme ) ! sum over the global domain 202 111 203 112 ! test of convergence … … 216 125 !**** 217 126 218 !,,,,,,,,,,,,,,,,,,,,,,,,,,,,,synchro,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,,219 220 127 ! indicator of non-convergence or explosion 221 222 128 IF( jn == nmax .OR. SQRT(epsr)/eps > 1.e+20 ) kindic = -2 223 129 IF( ncut == 999 ) GOTO 999 224 130 225 226 ! END of iterative loop 227 ! ===================== 228 229 END DO 230 131 ! ! ===================== 132 END DO ! END of iterative loop 133 ! ! ===================== 231 134 232 135 999 CONTINUE 233 136 234 137 235 ! 2.Output in gcx236 ! ------------- ----138 ! Output in gcx 139 ! ------------- 237 140 238 ! boundary conditions (est-ce necessaire? je ne crois pas!!!!) 239 240 #if defined key_dynspg_fsc 241 # if defined key_mpp 242 ! Mpp: export boundary values to neighbouring processors 243 CALL lbc_lnk( gcx, 'S', 1. ) 244 # else 245 IF( nperio /= 0 ) THEN 246 CALL lbc_lnk( gcx, 'T', 1. ) 247 ENDIF 248 # endif 249 #else 250 # if defined key_mpp 251 ! Mpp: export boundary values to neighbouring processors 252 CALL lbc_lnk( gcx, 'G', 1. ) 253 # else 254 IF( nperio /= 0 ) THEN 255 CALL lbc_lnk( gcx, 'F', 1. ) 256 ENDIF 257 # endif 258 #endif 141 CALL lbc_lnk( gcx, c_solver_pt, 1. ) ! boundary conditions (est-ce necessaire? je ne crois pas!!!!) 259 142 260 143 END SUBROUTINE sol_sor
Note: See TracChangeset
for help on using the changeset viewer.
