/[MITgcm]/MITgcm/model/src/cg3d.F
ViewVC logotype

Diff of /MITgcm/model/src/cg3d.F

Parent Directory Parent Directory | Revision Log Revision Log | View Revision Graph Revision Graph | View Patch Patch

revision 1.3 by adcroft, Mon May 24 14:15:15 1999 UTC revision 1.9 by heimbach, Mon May 14 21:33:30 2001 UTC
# Line 1  Line 1 
1  C $Header$  C $Header$
2    C $Name$
3    
4  #include "CPP_OPTIONS.h"  #include "CPP_OPTIONS.h"
5    
# Line 41  C     === Routine arguments === Line 42  C     === Routine arguments ===
42  C     myThid - Thread on which I am working.  C     myThid - Thread on which I am working.
43        INTEGER myThid        INTEGER myThid
44    
45    #ifdef ALLOW_NONHYDROSTATIC
46    
47  C     === Local variables ====  C     === Local variables ====
48  C     actualIts      - Number of iterations taken  C     actualIts      - Number of iterations taken
49  C     actualResidual - residual  C     actualResidual - residual
50  C     bi          - Block index in X and Y.  C     bi          - Block index in X and Y.
51  C     bj  C     bj
52  C     etaN        - Used in computing search directions  C     eta_qrN     - Used in computing search directions
53  C     etaNM1        suffix N and NM1 denote current and  C     eta_qrNM1     suffix N and NM1 denote current and
54  C     cgBeta        previous iterations respectively.  C     cgBeta        previous iterations respectively.
55  C     alpha    C     alpha  
56  C     sumRHS      - Sum of right-hand-side. Sometimes this is a  C     sumRHS      - Sum of right-hand-side. Sometimes this is a
# Line 62  C     I, J, N     - Loop counters ( N co Line 65  C     I, J, N     - Loop counters ( N co
65        INTEGER I, J, K, it3d        INTEGER I, J, K, it3d
66        INTEGER KM1, KP1        INTEGER KM1, KP1
67        _RL    err        _RL    err
68        _RL    etaN        _RL    eta_qrN
69        _RL    etaNM1        _RL    eta_qrNM1
70        _RL    cgBeta        _RL    cgBeta
71        _RL    alpha        _RL    alpha
72        _RL    sumRHS        _RL    sumRHS
# Line 77  C     I, J, N     - Loop counters ( N co Line 80  C     I, J, N     - Loop counters ( N co
80        INTEGER exchWidthX        INTEGER exchWidthX
81        INTEGER exchWidthY        INTEGER exchWidthY
82        INTEGER myNz        INTEGER myNz
83        _RL     topLevFac        _RL     topLevTerm
   
   
 CcnhDebugStarts  
       CHARACTER*(MAX_LEN_FNAM) suff  
 CcnhDebugEnds  
   
 #ifdef ALLOW_NONHYDROSTATIC  
84    
85  C--   Initialise inverter  C--   Initialise inverter
86        etaNM1              = 1. D0        eta_qrNM1 = 1. D0
87    
88  C--   Normalise RHS  C--   Normalise RHS
89        rhsMax = 0. _d 0        rhsMax = 0. _d 0
# Line 123  C--   Update overlaps Line 119  C--   Update overlaps
119        _EXCH_XYZ_R8( cg3d_b, myThid )        _EXCH_XYZ_R8( cg3d_b, myThid )
120        _EXCH_XYZ_R8( cg3d_x, myThid )        _EXCH_XYZ_R8( cg3d_x, myThid )
121    
122  #ifdef NONO  C--   Initial residual calculation (with free-Surface term)
 CcnhDebugStarts  
 C--   Initial residual calculation  
       err    = 0. _d 0  
       sumRHS = 0. _d 0  
       DO bj=myByLo(myThid),myByHi(myThid)  
        DO bi=myBxLo(myThid),myBxHi(myThid)  
          DO J=1,sNy  
           DO I=1,sNx  
           alpha = 0. _d 0  
           DO K=1,Nr  
            KM1 = K-1  
            IF ( KM1 .EQ. 0 )    KM1 = 1  
            KP1 = K+1  
            IF ( KP1 .EQ. Nr+1 ) KP1 = 1  
             cg3d_r(I,J,K,bi,bj) = cg3d_b(I,J,K,bi,bj) -( 0.  
      &     +aW3d(I  ,J  ,K  ,bi,bj)*cg3d_x(I-1,J  ,K  ,bi,bj)  
      &     +aW3d(I+1,J  ,K  ,bi,bj)*cg3d_x(I+1,J  ,K  ,bi,bj)  
      &     +aS3d(I  ,J  ,K  ,bi,bj)*cg3d_x(I  ,J-1,K  ,bi,bj)  
      &     +aS3d(I  ,J+1,K  ,bi,bj)*cg3d_x(I  ,J+1,K  ,bi,bj)  
      &     +aV3d(I  ,J  ,K  ,bi,bj)*cg3d_x(I  ,J  ,KM1,bi,bj)  
      &     +aV3d(I  ,J  ,KP1,bi,bj)*cg3d_x(I  ,J  ,KP1,bi,bj)  
      &     -aW3d(I  ,J  ,K  ,bi,bj)*cg3d_x(I  ,J  ,K  ,bi,bj)  
      &     -aW3d(I+1,J  ,K  ,bi,bj)*cg3d_x(I  ,J  ,K  ,bi,bj)  
      &     -aS3d(I  ,J  ,K  ,bi,bj)*cg3d_x(I  ,J  ,K  ,bi,bj)  
      &     -aS3d(I  ,J+1,K  ,bi,bj)*cg3d_x(I  ,J  ,K  ,bi,bj)  
      &     -aV3d(I  ,J  ,K  ,bi,bj)*cg3d_x(I  ,J  ,K  ,bi,bj)  
      &     -aV3d(I  ,J  ,KP1,bi,bj)*cg3d_x(I  ,J  ,K  ,bi,bj)  
      &     )  
            alpha = alpha  
      &     +cg3d_r(I,J,K,bi,bj)  
            sumRHS = sumRHS  
      &     +cg3d_b(I,J,K,bi,bj)  
           ENDDO  
           err = err + alpha*alpha  
          ENDDO  
         ENDDO  
        ENDDO  
       ENDDO  
       WRITE(6,*) 'DEBUG  mythid, err = ', mythid, SQRT(err)  
       _GLOBAL_SUM_R8( err   , myThid )  
       _GLOBAL_SUM_R8( sumRHS   , myThid )  
       _BEGIN_MASTER( myThid )  
       write(0,*) 'DEBUG cg3d: Sum(rhs) = ',sumRHS  
       _END_MASTER( )  
       actualIts      = 0  
       actualResidual = SQRT(err)  
 C     _BARRIER  
       _BEGIN_MASTER( myThid )  
        WRITE(0,'(A,I6,1PE30.14)') 'DEBUG  CG3D iters, err = ',  
      &  actualIts, actualResidual  
       _END_MASTER( )  
 CcnhDebugEnds  
 #endif  
   
 C--   Initial residual calculation  
123        err    = 0. _d 0        err    = 0. _d 0
124        sumRHS = 0. _d 0        sumRHS = 0. _d 0
125        DO bj=myByLo(myThid),myByHi(myThid)        DO bj=myByLo(myThid),myByHi(myThid)
# Line 188  C--   Initial residual calculation Line 129  C--   Initial residual calculation
129           IF ( K .EQ. 1 ) KM1 = 1           IF ( K .EQ. 1 ) KM1 = 1
130           KP1 = K+1           KP1 = K+1
131           IF ( K .EQ. Nr ) KP1 = 1           IF ( K .EQ. Nr ) KP1 = 1
132           topLevFac = 0.           topLevTerm = 0.
133           IF ( K .EQ. 1) topLevFac = 1.           IF ( K .EQ. 1) topLevTerm = freeSurfFac*cg3dNorm*
134         &      (horiVertRatio/gravity)/deltaTMom/deltaTMom
135           DO J=1,sNy           DO J=1,sNy
136            DO I=1,sNx            DO I=1,sNx
137             cg3d_s(I,J,K,bi,bj) = 0.             cg3d_s(I,J,K,bi,bj) = 0.
# Line 206  C--   Initial residual calculation Line 148  C--   Initial residual calculation
148       &     -aS3d(I  ,J+1,K  ,bi,bj)*cg3d_x(I  ,J  ,K  ,bi,bj)       &     -aS3d(I  ,J+1,K  ,bi,bj)*cg3d_x(I  ,J  ,K  ,bi,bj)
149       &     -aV3d(I  ,J  ,K  ,bi,bj)*cg3d_x(I  ,J  ,K  ,bi,bj)       &     -aV3d(I  ,J  ,K  ,bi,bj)*cg3d_x(I  ,J  ,K  ,bi,bj)
150       &     -aV3d(I  ,J  ,KP1,bi,bj)*cg3d_x(I  ,J  ,K  ,bi,bj)       &     -aV3d(I  ,J  ,KP1,bi,bj)*cg3d_x(I  ,J  ,K  ,bi,bj)
151       &     -freeSurfFac*_rA(i,j,bi,bj)* horiVertRatio*       &     -topLevTerm*_rA(I,J,bi,bj)*cg3d_x(I,J,K,bi,bj)
      &      cg3d_x(I  ,J  ,K,bi,bj)/deltaTMom/deltaTMom*cg3dNorm  
      &      *topLevFac  
152       &     )       &     )
153             err = err             err = err
154       &     +cg3d_r(I,J,K,bi,bj)*cg3d_r(I,J,K,bi,bj)       &     +cg3d_r(I,J,K,bi,bj)*cg3d_r(I,J,K,bi,bj)
# Line 228  C     _EXCH_XYZ_R8( cg3d_r, myThid ) Line 168  C     _EXCH_XYZ_R8( cg3d_r, myThid )
168         exchWidthY = 1         exchWidthY = 1
169         myNz       = Nr         myNz       = Nr
170         CALL EXCH_RL( cg3d_r,         CALL EXCH_RL( cg3d_r,
171       I            OLw, OLe, OLs, OLn, myNz,       I            OLw, OLe, OLs, OLn, sNx, sNy, myNz,
172       I            exchWidthX, exchWidthY,       I            exchWidthX, exchWidthY,
173       I            FORWARD_SIMULATION, EXCH_IGNORE_CORNERS, myThid )       I            FORWARD_SIMULATION, EXCH_IGNORE_CORNERS, myThid )
174  C     _EXCH_XYZ_R8( cg3d_s, myThid )  C     _EXCH_XYZ_R8( cg3d_s, myThid )
# Line 240  C     _EXCH_XYZ_R8( cg3d_s, myThid ) Line 180  C     _EXCH_XYZ_R8( cg3d_s, myThid )
180         exchWidthY = 1         exchWidthY = 1
181         myNz       = Nr         myNz       = Nr
182         CALL EXCH_RL( cg3d_s,         CALL EXCH_RL( cg3d_s,
183       I            OLw, OLe, OLs, OLn, myNz,       I            OLw, OLe, OLs, OLn, sNx, sNy, myNz,
184       I            exchWidthX, exchWidthY,       I            exchWidthX, exchWidthY,
185       I            FORWARD_SIMULATION, EXCH_IGNORE_CORNERS, myThid )       I            FORWARD_SIMULATION, EXCH_IGNORE_CORNERS, myThid )
186        _GLOBAL_SUM_R8( sumRHS, myThid )        _GLOBAL_SUM_R8( sumRHS, myThid )
187        _GLOBAL_SUM_R8( err   , myThid )        _GLOBAL_SUM_R8( err   , myThid )
188                
189        _BEGIN_MASTER( myThid )        _BEGIN_MASTER( myThid )
190        write(0,'(A,1PE30.14)') ' cg3d: Sum(rhs) = ',sumRHS        write(*,'(A,1PE30.14)') ' cg3d: Sum(rhs) = ',sumRHS
191        _END_MASTER( )        _END_MASTER( )
192    
193        actualIts      = 0        actualIts      = 0
194        actualResidual = SQRT(err)        actualResidual = SQRT(err)
195  C     _BARRIER  C     _BARRIER
196        _BEGIN_MASTER( myThid )        _BEGIN_MASTER( myThid )
197         WRITE(0,'(A,I6,1PE30.14)') ' CG3D iters, err = ',         WRITE(*,'(A,I6,1PE30.14)') ' CG3D iters, err = ',
198       &  actualIts, actualResidual       &  actualIts, actualResidual
199        _END_MASTER( )        _END_MASTER( )
200    
# Line 264  C     >>>>>>>>>>>>>>> BEGIN SOLVER <<<<< Line 204  C     >>>>>>>>>>>>>>> BEGIN SOLVER <<<<<
204  CcnhDebugStarts  CcnhDebugStarts
205  #ifdef VERBOSE  #ifdef VERBOSE
206         IF ( mod(it3d-1,10).EQ.0)         IF ( mod(it3d-1,10).EQ.0)
207       &   WRITE(0,*) ' CG3D: Iteration ',it3d-1,       &   WRITE(*,*) ' CG3D: Iteration ',it3d-1,
208       &      ' residual = ',actualResidual       &      ' residual = ',actualResidual
209  #endif  #endif
210  CcnhDebugEnds  CcnhDebugEnds
# Line 274  C--    conjugate direction vector "s". Line 214  C--    conjugate direction vector "s".
214  C      Note. On the next to loops over all tiles the inner loop ranges  C      Note. On the next to loops over all tiles the inner loop ranges
215  C            in sNx and sNy are expanded by 1 to avoid a communication  C            in sNx and sNy are expanded by 1 to avoid a communication
216  C            step. However this entails a bit of gynamastics because we only  C            step. However this entails a bit of gynamastics because we only
217  C            want etaN for the interior points.  C            want eta_qrN for the interior points.
218         etaN = 0. _d 0         eta_qrN = 0. _d 0
219         DO bj=myByLo(myThid),myByHi(myThid)         DO bj=myByLo(myThid),myByHi(myThid)
220          DO bi=myBxLo(myThid),myBxHi(myThid)          DO bi=myBxLo(myThid),myBxHi(myThid)
221           DO K=1,1           DO K=1,1
# Line 307  caja       ENDDO Line 247  caja       ENDDO
247  caja      ENDIF  caja      ENDIF
248            DO J=1,sNy            DO J=1,sNy
249             DO I=1,sNx             DO I=1,sNx
250              etaN = etaN              eta_qrN = eta_qrN
251       &      +cg3d_q(I,J,K,bi,bj)*cg3d_r(I,J,K,bi,bj)       &      +cg3d_q(I,J,K,bi,bj)*cg3d_r(I,J,K,bi,bj)
252             ENDDO             ENDDO
253            ENDDO            ENDDO
# Line 322  caja      ENDIF Line 262  caja      ENDIF
262            ENDDO            ENDDO
263            DO J=1,sNy            DO J=1,sNy
264             DO I=1,sNx             DO I=1,sNx
265              etaN = etaN              eta_qrN = eta_qrN
266       &      +cg3d_q(I,J,K,bi,bj)*cg3d_r(I,J,K,bi,bj)       &      +cg3d_q(I,J,K,bi,bj)*cg3d_r(I,J,K,bi,bj)
267             ENDDO             ENDDO
268            ENDDO            ENDDO
# Line 330  caja      ENDIF Line 270  caja      ENDIF
270          ENDDO          ENDDO
271         ENDDO         ENDDO
272  caja  caja
273  caja  etaN=0.  caja  eta_qrN=0.
274  caja   DO bj=myByLo(myThid),myByHi(myThid)  caja   DO bj=myByLo(myThid),myByHi(myThid)
275  caja    DO bi=myBxLo(myThid),myBxHi(myThid)  caja    DO bi=myBxLo(myThid),myBxHi(myThid)
276  caja     DO K=1,Nr  caja     DO K=1,Nr
277  caja      DO J=1,sNy  caja      DO J=1,sNy
278  caja       DO I=1,sNx  caja       DO I=1,sNx
279  caja        etaN = etaN  caja        eta_qrN = eta_qrN
280  caja &      +cg3d_q(I,J,K,bi,bj)*cg3d_r(I,J,K,bi,bj)  caja &      +cg3d_q(I,J,K,bi,bj)*cg3d_r(I,J,K,bi,bj)
281  caja       ENDDO  caja       ENDDO
282  caja      ENDDO  caja      ENDDO
# Line 345  caja    ENDDO Line 285  caja    ENDDO
285  caja   ENDDO  caja   ENDDO
286  caja  caja
287    
288         _GLOBAL_SUM_R8(etaN, myThid)         _GLOBAL_SUM_R8(eta_qrN, myThid)
289  CcnhDebugStarts  CcnhDebugStarts
290  C      WRITE(0,*) ' CG3D: Iteration ',it3d-1,' etaN = ',etaN  C      WRITE(*,*) ' CG3D: Iteration ',it3d-1,' eta_qrN = ',eta_qrN
291  CcnhDebugEnds  CcnhDebugEnds
292         cgBeta   = etaN/etaNM1         cgBeta   = eta_qrN/eta_qrNM1
293  CcnhDebugStarts  CcnhDebugStarts
294  C      WRITE(0,*) ' CG3D: Iteration ',it3d-1,' beta = ',cgBeta  C      WRITE(*,*) ' CG3D: Iteration ',it3d-1,' beta = ',cgBeta
295  CcnhDebugEnds  CcnhDebugEnds
296         etaNM1 = etaN         eta_qrNM1 = eta_qrN
297    
298         DO bj=myByLo(myThid),myByHi(myThid)         DO bj=myByLo(myThid),myByHi(myThid)
299          DO bi=myBxLo(myThid),myBxHi(myThid)          DO bi=myBxLo(myThid),myBxHi(myThid)
# Line 371  CcnhDebugEnds Line 311  CcnhDebugEnds
311  C==    Evaluate laplace operator on conjugate gradient vector  C==    Evaluate laplace operator on conjugate gradient vector
312  C==    q = A.s  C==    q = A.s
313         alpha = 0. _d 0         alpha = 0. _d 0
314           topLevTerm = freeSurfFac*cg3dNorm*
315         &      (horiVertRatio/gravity)/deltaTMom/deltaTMom
316         DO bj=myByLo(myThid),myByHi(myThid)         DO bj=myByLo(myThid),myByHi(myThid)
317          DO bi=myBxLo(myThid),myBxHi(myThid)          DO bi=myBxLo(myThid),myBxHi(myThid)
318           IF ( Nr .GT. 1 ) THEN           IF ( Nr .GT. 1 ) THEN
# Line 388  C==    q = A.s Line 330  C==    q = A.s
330       &      -aS3d(I  ,J  ,K  ,bi,bj)*cg3d_s(I  ,J  ,K  ,bi,bj)       &      -aS3d(I  ,J  ,K  ,bi,bj)*cg3d_s(I  ,J  ,K  ,bi,bj)
331       &      -aS3d(I  ,J+1,K  ,bi,bj)*cg3d_s(I  ,J  ,K  ,bi,bj)       &      -aS3d(I  ,J+1,K  ,bi,bj)*cg3d_s(I  ,J  ,K  ,bi,bj)
332       &      -aV3d(I  ,J  ,K+1,bi,bj)*cg3d_s(I  ,J  ,K  ,bi,bj)       &      -aV3d(I  ,J  ,K+1,bi,bj)*cg3d_s(I  ,J  ,K  ,bi,bj)
333       &      -freeSurfFac*_rA(i,j,bi,bj)* horiVertRatio*       &      -topLevTerm*_rA(I,J,bi,bj)*cg3d_s(I,J,K,bi,bj)
      &       cg3d_s(I  ,J  ,K,bi,bj)/deltaTMom/deltaTMom*cg3dNorm  
334               alpha = alpha+cg3d_s(I,J,K,bi,bj)*cg3d_q(I,J,K,bi,bj)               alpha = alpha+cg3d_s(I,J,K,bi,bj)*cg3d_q(I,J,K,bi,bj)
335              ENDDO              ENDDO
336             ENDDO             ENDDO
# Line 407  C==    q = A.s Line 348  C==    q = A.s
348       &      -aW3d(I+1,J  ,K  ,bi,bj)*cg3d_s(I  ,J  ,K  ,bi,bj)       &      -aW3d(I+1,J  ,K  ,bi,bj)*cg3d_s(I  ,J  ,K  ,bi,bj)
349       &      -aS3d(I  ,J  ,K  ,bi,bj)*cg3d_s(I  ,J  ,K  ,bi,bj)       &      -aS3d(I  ,J  ,K  ,bi,bj)*cg3d_s(I  ,J  ,K  ,bi,bj)
350       &      -aS3d(I  ,J+1,K  ,bi,bj)*cg3d_s(I  ,J  ,K  ,bi,bj)       &      -aS3d(I  ,J+1,K  ,bi,bj)*cg3d_s(I  ,J  ,K  ,bi,bj)
351       &      -freeSurfFac*_rA(i,j,bi,bj)* horiVertRatio*       &      -topLevTerm*_rA(I,J,bi,bj)*cg3d_s(I,J,K,bi,bj)
      &       cg3d_s(I  ,J  ,K,bi,bj)/deltaTMom/deltaTMom*cg3dNorm  
352               alpha = alpha+cg3d_s(I,J,K,bi,bj)*cg3d_q(I,J,K,bi,bj)               alpha = alpha+cg3d_s(I,J,K,bi,bj)*cg3d_q(I,J,K,bi,bj)
353              ENDDO              ENDDO
354             ENDDO             ENDDO
# Line 458  C==    q = A.s Line 398  C==    q = A.s
398         ENDDO         ENDDO
399         _GLOBAL_SUM_R8(alpha,myThid)         _GLOBAL_SUM_R8(alpha,myThid)
400  CcnhDebugStarts  CcnhDebugStarts
401  C      WRITE(0,*) ' CG3D: Iteration ',it3d-1,' SUM(s*q)= ',alpha  C      WRITE(*,*) ' CG3D: Iteration ',it3d-1,' SUM(s*q)= ',alpha
402  CcnhDebugEnds  CcnhDebugEnds
403         alpha = etaN/alpha         alpha = eta_qrN/alpha
404  CcnhDebugStarts  CcnhDebugStarts
405  C      WRITE(0,*) ' CG3D: Iteration ',it3d-1,' alpha= ',alpha  C      WRITE(*,*) ' CG3D: Iteration ',it3d-1,' alpha= ',alpha
406  CcnhDebugEnds  CcnhDebugEnds
407            
408  C==    Update solution and residual vectors  C==    Update solution and residual vectors
# Line 498  C      _EXCH_XYZ_R8(cg3d_r, myThid ) Line 438  C      _EXCH_XYZ_R8(cg3d_r, myThid )
438         exchWidthY = 1         exchWidthY = 1
439         myNz       = Nr         myNz       = Nr
440         CALL EXCH_RL( cg3d_r,         CALL EXCH_RL( cg3d_r,
441       I             OLw, OLe, OLs, OLn, myNz,       I             OLw, OLe, OLs, OLn, sNx, sNy, myNz,
442       I             exchWidthX, exchWidthY,       I             exchWidthX, exchWidthY,
443       I             FORWARD_SIMULATION, EXCH_IGNORE_CORNERS, myThid )       I             FORWARD_SIMULATION, EXCH_IGNORE_CORNERS, myThid )
444    
# Line 520  C--   Un-normalise the answer Line 460  C--   Un-normalise the answer
460    
461  Cadj  _EXCH_XYZ_R8(cg3d_x, myThid )  Cadj  _EXCH_XYZ_R8(cg3d_x, myThid )
462        _BEGIN_MASTER( myThid )        _BEGIN_MASTER( myThid )
463         WRITE(0,'(A,I6,1PE30.14)') ' CG3D iters, err = ',         WRITE(*,'(A,I6,1PE30.14)') ' CG3D iters, err = ',
464       &  actualIts, actualResidual       &  actualIts, actualResidual
465        _END_MASTER( )        _END_MASTER( )
466    
 CcnhDebugStarts  
 C     CALL PLOT_FIELD_XYRL( cg2d_x, 'CALC_MOM_RHS CG2D_X' , 1, myThid )  
 C     err    = 0. _d 0  
 C     DO bj=myByLo(myThid),myByHi(myThid)  
 C      DO bi=myBxLo(myThid),myBxHi(myThid)  
 C       DO J=1,sNy  
 C        DO I=1,sNx  
 C         cg2d_r(I,J,bi,bj) = cg2d_b(I,J,bi,bj) -  
 C    &    (aW2d(I  ,J  ,bi,bj)*cg2d_x(I-1,J  ,bi,bj)  
 C    &    +aW2d(I+1,J  ,bi,bj)*cg2d_x(I+1,J  ,bi,bj)  
 C    &    +aS2d(I  ,J  ,bi,bj)*cg2d_x(I  ,J-1,bi,bj)  
 C    &    +aS2d(I  ,J+1,bi,bj)*cg2d_x(I  ,J+1,bi,bj)  
 C    &    -aW2d(I  ,J  ,bi,bj)*cg2d_x(I  ,J  ,bi,bj)  
 C    &    -aW2d(I+1,J  ,bi,bj)*cg2d_x(I  ,J  ,bi,bj)  
 C    &    -aS2d(I  ,J  ,bi,bj)*cg2d_x(I  ,J  ,bi,bj)  
 C    &    -aS2d(I  ,J+1,bi,bj)*cg2d_x(I  ,J  ,bi,bj))  
 C         err            = err            +  
 C    &     cg2d_r(I,J,bi,bj)*cg2d_r(I,J,bi,bj)  
 C        ENDDO  
 C       ENDDO  
 C      ENDDO  
 C     ENDDO  
 C     _GLOBAL_SUM_R8( err   , myThid )  
 C     write(0,*) 'cg2d: Ax - b = ',SQRT(err)  
 CcnhDebugEnds  
   
467  #endif /* ALLOW_NONHYDROSTATIC */  #endif /* ALLOW_NONHYDROSTATIC */
468    
469        RETURN        RETURN
Legend:
Removed from v.1.3  
changed lines
  Added in v.1.9
  ViewVC Help
Powered by ViewVC 1.1.22