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

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

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

revision 1.18 by adcroft, Wed Dec 9 16:11:51 1998 UTC revision 1.33 by adcroft, Tue May 29 14:01:36 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    
6        SUBROUTINE CG2D(          SUBROUTINE CG2D(  
7       I                cg2d_b,       I                cg2d_b,
8       U                cg2d_x,       U                cg2d_x,
9         O                firstResidual,
10         O                lastResidual,
11         U                numIters,
12       I                myThid )       I                myThid )
13  C     /==========================================================\  C     /==========================================================\
14  C     | SUBROUTINE CG2D                                          |  C     | SUBROUTINE CG2D                                          |
# Line 35  C     === Global data === Line 39  C     === Global data ===
39  #include "EEPARAMS.h"  #include "EEPARAMS.h"
40  #include "PARAMS.h"  #include "PARAMS.h"
41  #include "GRID.h"  #include "GRID.h"
42  #include "CG2D_INTERNAL.h"  #include "CG2D.h"
43    #include "SURFACE.h"
44    
45  C     === Routine arguments ===  C     === Routine arguments ===
46  C     myThid - Thread on which I am working.  C     myThid    - Thread on which I am working.
47        INTEGER myThid  C     cg2d_b    - The source term or "right hand side"
48        _RL  cg2d_x(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)  C     cg2d_x    - The solution
49    C     firstResidual - the initial residual before any iterations
50    C     lastResidual  - the actual residual reached
51    C     numIters  - Entry: the maximum number of iterations allowed
52    C                 Exit:  the actual number of iterations used
53        _RL  cg2d_b(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)        _RL  cg2d_b(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
54          _RL  cg2d_x(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy)
55          _RL  firstResidual
56          _RL  lastResidual
57          INTEGER numIters
58          INTEGER myThid
59    
60  C     === Local variables ====  C     === Local variables ====
61  C     actualIts      - Number of iterations taken  C     actualIts      - Number of iterations taken
62  C     actualResidual - residual  C     actualResidual - residual
63  C     bi          - Block index in X and Y.  C     bi          - Block index in X and Y.
64  C     bj  C     bj
65  C     etaN        - Used in computing search directions  C     eta_qrN     - Used in computing search directions
66  C     etaNM1        suffix N and NM1 denote current and  C     eta_qrNM1     suffix N and NM1 denote current and
67  C     cgBeta        previous iterations respectively.  C     cgBeta        previous iterations respectively.
68  C     alpha    C     alpha  
69  C     sumRHS      - Sum of right-hand-side. Sometimes this is a  C     sumRHS      - Sum of right-hand-side. Sometimes this is a
# Line 64  C     I, J, N     - Loop counters ( N co Line 77  C     I, J, N     - Loop counters ( N co
77        INTEGER bi, bj                      INTEGER bi, bj              
78        INTEGER I, J, it2d        INTEGER I, J, it2d
79        _RL    err        _RL    err
80        _RL    etaN        _RL    eta_qrN
81        _RL    etaNM1        _RL    eta_qrNM1
82        _RL    cgBeta        _RL    cgBeta
83        _RL    alpha        _RL    alpha
84        _RL    sumRHS        _RL    sumRHS
# Line 82  C     I, J, N     - Loop counters ( N co Line 95  C     I, J, N     - Loop counters ( N co
95    
96    
97  CcnhDebugStarts  CcnhDebugStarts
98        CHARACTER*(MAX_LEN_FNAM) suff  C     CHARACTER*(MAX_LEN_FNAM) suff
99  CcnhDebugEnds  CcnhDebugEnds
100    
101    
102  C--   Initialise inverter  C--   Initialise inverter
103        etaNBuf(1,myThid)   = 0. _d 0        eta_qrNM1 = 1. _d 0
       errBuf(1,myThid)    = 0. _d 0  
       sumRHSBuf(1,myThid) = 0. _d 0  
       etaNM1              = 1. _d 0  
104    
105  CcnhDebugStarts  CcnhDebugStarts
106  C     _EXCH_XY_R8( cg2d_b, myThid )  C     _EXCH_XY_R8( cg2d_b, myThid )
# Line 112  C--   Normalise RHS Line 122  C--   Normalise RHS
122          ENDDO          ENDDO
123         ENDDO         ENDDO
124        ENDDO        ENDDO
125        rhsMaxBuf(1,myThid) = rhsMax  
126        _GLOBAL_MAX_R8( rhsMaxbuf, rhsMax, myThid )        IF (cg2dNormaliseRHS) THEN
127        rhsMax =  rhsMaxBuf(1,1)  C-  Normalise RHS :
128    #ifdef LETS_MAKE_JAM
129    C     _GLOBAL_MAX_R8( rhsMax, myThid )
130          rhsMax=1.
131    #else
132          _GLOBAL_MAX_R8( rhsMax, myThid )
133    Catm  rhsMax=1.
134    #endif
135        rhsNorm = 1. _d 0        rhsNorm = 1. _d 0
136        IF ( rhsMax .NE. 0. ) rhsNorm = 1. _d 0 / rhsMax        IF ( rhsMax .NE. 0. ) rhsNorm = 1. _d 0 / rhsMax
137        DO bj=myByLo(myThid),myByHi(myThid)        DO bj=myByLo(myThid),myByHi(myThid)
# Line 127  C--   Normalise RHS Line 144  C--   Normalise RHS
144          ENDDO          ENDDO
145         ENDDO         ENDDO
146        ENDDO        ENDDO
147    C- end Normalise RHS
148          ENDIF
149    
150  C--   Update overlaps  C--   Update overlaps
151        _EXCH_XY_R8( cg2d_b, myThid )        _EXCH_XY_R8( cg2d_b, myThid )
# Line 154  C--   Initial residual calculation Line 173  C--   Initial residual calculation
173       &    -aW2d(I+1,J  ,bi,bj)*cg2d_x(I  ,J  ,bi,bj)       &    -aW2d(I+1,J  ,bi,bj)*cg2d_x(I  ,J  ,bi,bj)
174       &    -aS2d(I  ,J  ,bi,bj)*cg2d_x(I  ,J  ,bi,bj)       &    -aS2d(I  ,J  ,bi,bj)*cg2d_x(I  ,J  ,bi,bj)
175       &    -aS2d(I  ,J+1,bi,bj)*cg2d_x(I  ,J  ,bi,bj)       &    -aS2d(I  ,J+1,bi,bj)*cg2d_x(I  ,J  ,bi,bj)
176       &    -freeSurfFac*_rA(i,j,bi,bj)* horiVertRatio*       &    -freeSurfFac*_rA(i,j,bi,bj)*recip_Bo(i,j,bi,bj)*
177       &     cg2d_x(I  ,J  ,bi,bj)/deltaTMom/deltaTMom*cg2dNorm       &     cg2d_x(I  ,J  ,bi,bj)/deltaTMom/deltaTMom*cg2dNorm
178       &    )       &    )
179            err            = err            +            err            = err            +
# Line 166  C--   Initial residual calculation Line 185  C--   Initial residual calculation
185         ENDDO         ENDDO
186        ENDDO        ENDDO
187  C     _EXCH_XY_R8( cg2d_r, myThid )  C     _EXCH_XY_R8( cg2d_r, myThid )
188    #ifdef LETS_MAKE_JAM
189          CALL EXCH_XY_O1_R8_JAM( cg2d_r )
190    #else
191         OLw        = 1         OLw        = 1
192         OLe        = 1         OLe        = 1
193         OLn        = 1         OLn        = 1
# Line 173  C     _EXCH_XY_R8( cg2d_r, myThid ) Line 195  C     _EXCH_XY_R8( cg2d_r, myThid )
195         exchWidthX = 1         exchWidthX = 1
196         exchWidthY = 1         exchWidthY = 1
197         myNz       = 1         myNz       = 1
198          IF (useCubedSphereExchange) THEN
199           CALL EXCH_RL_CUBE( cg2d_r,
200         I            OLw, OLe, OLs, OLn, myNz,
201         I            exchWidthX, exchWidthY,
202         I            FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )
203          ELSE
204         CALL EXCH_RL( cg2d_r,         CALL EXCH_RL( cg2d_r,
205       I            OLw, OLe, OLs, OLn, myNz,       I            OLw, OLe, OLs, OLn, myNz,
206       I            exchWidthX, exchWidthY,       I            exchWidthX, exchWidthY,
207       I            FORWARD_SIMULATION, EXCH_IGNORE_CORNERS, myThid )       I            FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )
208          ENDIF
209    #endif
210  C     _EXCH_XY_R8( cg2d_s, myThid )  C     _EXCH_XY_R8( cg2d_s, myThid )
211    #ifdef LETS_MAKE_JAM
212          CALL EXCH_XY_O1_R8_JAM( cg2d_s )
213    #else
214         OLw        = 1         OLw        = 1
215         OLe        = 1         OLe        = 1
216         OLn        = 1         OLn        = 1
# Line 185  C     _EXCH_XY_R8( cg2d_s, myThid ) Line 218  C     _EXCH_XY_R8( cg2d_s, myThid )
218         exchWidthX = 1         exchWidthX = 1
219         exchWidthY = 1         exchWidthY = 1
220         myNz       = 1         myNz       = 1
221          IF (useCubedSphereExchange) THEN
222           CALL EXCH_RL_CUBE( cg2d_s,
223         I            OLw, OLe, OLs, OLn, myNz,
224         I            exchWidthX, exchWidthY,
225         I            FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )
226          ELSE
227         CALL EXCH_RL( cg2d_s,         CALL EXCH_RL( cg2d_s,
228       I            OLw, OLe, OLs, OLn, myNz,       I            OLw, OLe, OLs, OLn, myNz,
229       I            exchWidthX, exchWidthY,       I            exchWidthX, exchWidthY,
230       I            FORWARD_SIMULATION, EXCH_IGNORE_CORNERS, myThid )       I            FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )
231        sumRHSBuf(1,myThid) = sumRHS        ENDIF
232        _GLOBAL_SUM_R8( sumRHSBuf , sumRHS, myThid )  #endif
233        sumRHS = sumRHSBuf(1,1)         _GLOBAL_SUM_R8( sumRHS, myThid )
234        errBuf(1,myThid)    = err         _GLOBAL_SUM_R8( err   , myThid )
235  C     WRITE(6,*) ' mythid, err = ', mythid, SQRT(err)         err = SQRT(err)
236        _GLOBAL_SUM_R8( errBuf    , err   , myThid )         actualIts      = 0
237        err    = errBuf(1,1)         actualResidual = err
238                
239        _BEGIN_MASTER( myThid )        _BEGIN_MASTER( myThid )
240        write(0,'(A,1PE30.14)') 'cg2d: Sum(rhs) = ',sumRHS        write(*,'(A,1P2E22.14)')' cg2d: Sum(rhs),rhsMax = ',
241         &                                  sumRHS,rhsMax
242        _END_MASTER( )        _END_MASTER( )
   
       actualIts      = 0  
       actualResidual = SQRT(err)  
243  C     _BARRIER  C     _BARRIER
244        _BEGIN_MASTER( myThid )  c     _BEGIN_MASTER( myThid )
245         WRITE(0,'(A,I6,1PE30.14)') ' CG2D iters, err = ',  c      WRITE(*,'(A,I6,1PE30.14)') ' CG2D iters, err = ',
246       & actualIts, actualResidual  c    & actualIts, actualResidual
247        _END_MASTER( )  c     _END_MASTER( )
248          firstResidual=actualResidual
249    
250          IF ( err .LT. cg2dTolerance ) GOTO 11
251    
252  C     >>>>>>>>>>>>>>> BEGIN SOLVER <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<  C     >>>>>>>>>>>>>>> BEGIN SOLVER <<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<<
253        DO 10 it2d=1, cg2dMaxIters        DO 10 it2d=1, numIters
254    
255  CcnhDebugStarts  CcnhDebugStarts
256  C      WRITE(0,*) ' CG2D: Iteration ',it2d-1,' residual = ',  C      WRITE(*,*) ' CG2D: Iteration ',it2d-1,' residual = ',
257  C    &  actualResidual  C    &  actualResidual
258  CcnhDebugEnds  CcnhDebugEnds
        IF ( err .LT. cg2dTargetResidual ) GOTO 11  
259  C--    Solve preconditioning equation and update  C--    Solve preconditioning equation and update
260  C--    conjugate direction vector "s".  C--    conjugate direction vector "s".
261         etaN = 0. _d 0         eta_qrN = 0. _d 0
262         DO bj=myByLo(myThid),myByHi(myThid)         DO bj=myByLo(myThid),myByHi(myThid)
263          DO bi=myBxLo(myThid),myBxHi(myThid)          DO bi=myBxLo(myThid),myBxHi(myThid)
264           DO J=1,sNy           DO J=1,sNy
# Line 233  C--    conjugate direction vector "s". Line 272  C--    conjugate direction vector "s".
272  CcnhDebugStarts  CcnhDebugStarts
273  C          cg2d_q(I,J,bi,bj) = cg2d_r(I  ,J  ,bi,bj)  C          cg2d_q(I,J,bi,bj) = cg2d_r(I  ,J  ,bi,bj)
274  CcnhDebugEnds  CcnhDebugEnds
275             etaN = etaN             eta_qrN = eta_qrN
276       &     +cg2d_q(I,J,bi,bj)*cg2d_r(I,J,bi,bj)       &     +cg2d_q(I,J,bi,bj)*cg2d_r(I,J,bi,bj)
277            ENDDO            ENDDO
278           ENDDO           ENDDO
279          ENDDO          ENDDO
280         ENDDO         ENDDO
281    
282         etanBuf(1,myThid) = etaN         _GLOBAL_SUM_R8(eta_qrN, myThid)
        _GLOBAL_SUM_R8(etaNbuf,etaN, myThid)  
        etaN = etaNBuf(1,1)  
283  CcnhDebugStarts  CcnhDebugStarts
284  C      WRITE(0,*) ' CG2D: Iteration ',it2d-1,' etaN = ',etaN  C      WRITE(*,*) ' CG2D: Iteration ',it2d-1,' eta_qrN = ',eta_qrN
285  CcnhDebugEnds  CcnhDebugEnds
286         cgBeta   = etaN/etaNM1         cgBeta   = eta_qrN/eta_qrNM1
287  CcnhDebugStarts  CcnhDebugStarts
288  C      WRITE(0,*) ' CG2D: Iteration ',it2d-1,' beta = ',cgBeta  C      WRITE(*,*) ' CG2D: Iteration ',it2d-1,' beta = ',cgBeta
289  CcnhDebugEnds  CcnhDebugEnds
290         etaNM1 = etaN         eta_qrNM1 = eta_qrN
291    
292         DO bj=myByLo(myThid),myByHi(myThid)         DO bj=myByLo(myThid),myByHi(myThid)
293          DO bi=myBxLo(myThid),myBxHi(myThid)          DO bi=myBxLo(myThid),myBxHi(myThid)
# Line 266  CcnhDebugEnds Line 303  CcnhDebugEnds
303  C--    Do exchanges that require messages i.e. between  C--    Do exchanges that require messages i.e. between
304  C--    processes.  C--    processes.
305  C      _EXCH_XY_R8( cg2d_s, myThid )  C      _EXCH_XY_R8( cg2d_s, myThid )
306    #ifdef LETS_MAKE_JAM
307          CALL EXCH_XY_O1_R8_JAM( cg2d_s )
308    #else
309         OLw        = 1         OLw        = 1
310         OLe        = 1         OLe        = 1
311         OLn        = 1         OLn        = 1
# Line 273  C      _EXCH_XY_R8( cg2d_s, myThid ) Line 313  C      _EXCH_XY_R8( cg2d_s, myThid )
313         exchWidthX = 1         exchWidthX = 1
314         exchWidthY = 1         exchWidthY = 1
315         myNz       = 1         myNz       = 1
316          IF (useCubedSphereExchange) THEN
317           CALL EXCH_RL_CUBE( cg2d_s,
318         I            OLw, OLe, OLs, OLn, myNz,
319         I            exchWidthX, exchWidthY,
320         I            FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )
321          ELSE
322         CALL EXCH_RL( cg2d_s,         CALL EXCH_RL( cg2d_s,
323       I            OLw, OLe, OLs, OLn, myNz,       I            OLw, OLe, OLs, OLn, myNz,
324       I            exchWidthX, exchWidthY,       I            exchWidthX, exchWidthY,
325       I            FORWARD_SIMULATION, EXCH_IGNORE_CORNERS, myThid )       I            FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )
326          ENDIF
327    #endif
328    
329  C==    Evaluate laplace operator on conjugate gradient vector  C==    Evaluate laplace operator on conjugate gradient vector
330  C==    q = A.s  C==    q = A.s
# Line 295  C==    q = A.s Line 342  C==    q = A.s
342       &    -aW2d(I+1,J  ,bi,bj)*cg2d_s(I  ,J  ,bi,bj)       &    -aW2d(I+1,J  ,bi,bj)*cg2d_s(I  ,J  ,bi,bj)
343       &    -aS2d(I  ,J  ,bi,bj)*cg2d_s(I  ,J  ,bi,bj)       &    -aS2d(I  ,J  ,bi,bj)*cg2d_s(I  ,J  ,bi,bj)
344       &    -aS2d(I  ,J+1,bi,bj)*cg2d_s(I  ,J  ,bi,bj)       &    -aS2d(I  ,J+1,bi,bj)*cg2d_s(I  ,J  ,bi,bj)
345       &    -freeSurfFac*_rA(i,j,bi,bj)* horiVertRatio*       &    -freeSurfFac*_rA(i,j,bi,bj)*recip_Bo(i,j,bi,bj)*
346       &     cg2d_s(I  ,J  ,bi,bj)/deltaTMom/deltaTMom*cg2dNorm       &     cg2d_s(I  ,J  ,bi,bj)/deltaTMom/deltaTMom*cg2dNorm
347            alpha = alpha+cg2d_s(I,J,bi,bj)*cg2d_q(I,J,bi,bj)            alpha = alpha+cg2d_s(I,J,bi,bj)*cg2d_q(I,J,bi,bj)
348            ENDDO            ENDDO
349           ENDDO           ENDDO
350          ENDDO          ENDDO
351         ENDDO         ENDDO
352         alphaBuf(1,myThid) = alpha         _GLOBAL_SUM_R8(alpha,myThid)
        _GLOBAL_SUM_R8(alphaBuf,alpha,myThid)  
        alpha = alphaBuf(1,1)  
353  CcnhDebugStarts  CcnhDebugStarts
354  C      WRITE(0,*) ' CG2D: Iteration ',it2d-1,' SUM(s*q)= ',alpha  C      WRITE(*,*) ' CG2D: Iteration ',it2d-1,' SUM(s*q)= ',alpha
355  CcnhDebugEnds  CcnhDebugEnds
356         alpha = etaN/alpha         alpha = eta_qrN/alpha
357  CcnhDebugStarts  CcnhDebugStarts
358  C      WRITE(0,*) ' CG2D: Iteration ',it2d-1,' alpha= ',alpha  C      WRITE(*,*) ' CG2D: Iteration ',it2d-1,' alpha= ',alpha
359  CcnhDebugEnds  CcnhDebugEnds
360            
361  C==    Update solution and residual vectors  C==    Update solution and residual vectors
# Line 328  C      Now compute "interior" points. Line 373  C      Now compute "interior" points.
373          ENDDO          ENDDO
374         ENDDO         ENDDO
375    
376         errBuf(1,myThid) = err         _GLOBAL_SUM_R8( err   , myThid )
        _GLOBAL_SUM_R8( errBuf    , err   , myThid )  
        err = errBuf(1,1)  
377         err = SQRT(err)         err = SQRT(err)
378         actualIts      = it2d         actualIts      = it2d
379         actualResidual = err         actualResidual = err
380         IF ( err .LT. cg2dTargetResidual ) GOTO 11         IF ( err .LT. cg2dTolerance ) GOTO 11
381  C      _EXCH_XY_R8(cg2d_r, myThid )  C      _EXCH_XY_R8(cg2d_r, myThid )
382    #ifdef LETS_MAKE_JAM
383          CALL EXCH_XY_O1_R8_JAM( cg2d_r )
384    #else
385         OLw        = 1         OLw        = 1
386         OLe        = 1         OLe        = 1
387         OLn        = 1         OLn        = 1
# Line 343  C      _EXCH_XY_R8(cg2d_r, myThid ) Line 389  C      _EXCH_XY_R8(cg2d_r, myThid )
389         exchWidthX = 1         exchWidthX = 1
390         exchWidthY = 1         exchWidthY = 1
391         myNz       = 1         myNz       = 1
392          IF (useCubedSphereExchange) THEN
393           CALL EXCH_RL_CUBE( cg2d_r,
394         I            OLw, OLe, OLs, OLn, myNz,
395         I            exchWidthX, exchWidthY,
396         I            FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )
397          ELSE
398         CALL EXCH_RL( cg2d_r,         CALL EXCH_RL( cg2d_r,
399       I             OLw, OLe, OLs, OLn, myNz,       I            OLw, OLe, OLs, OLn, myNz,
400       I             exchWidthX, exchWidthY,       I            exchWidthX, exchWidthY,
401       I             FORWARD_SIMULATION, EXCH_IGNORE_CORNERS, myThid )       I            FORWARD_SIMULATION, EXCH_UPDATE_CORNERS, myThid )
402          ENDIF
403    #endif
404    
405     10 CONTINUE     10 CONTINUE
406     11 CONTINUE     11 CONTINUE
407    
408          IF (cg2dNormaliseRHS) THEN
409  C--   Un-normalise the answer  C--   Un-normalise the answer
410        DO bj=myByLo(myThid),myByHi(myThid)          DO bj=myByLo(myThid),myByHi(myThid)
411         DO bi=myBxLo(myThid),myBxHi(myThid)           DO bi=myBxLo(myThid),myBxHi(myThid)
412          DO J=1,sNy            DO J=1,sNy
413           DO I=1,sNx             DO I=1,sNx
414            cg2d_x(I  ,J  ,bi,bj) = cg2d_x(I  ,J  ,bi,bj)/rhsNorm              cg2d_x(I  ,J  ,bi,bj) = cg2d_x(I  ,J  ,bi,bj)/rhsNorm
415               ENDDO
416              ENDDO
417           ENDDO           ENDDO
418          ENDDO          ENDDO
419         ENDDO        ENDIF
       ENDDO  
420    
421        _EXCH_XY_R8(cg2d_x, myThid )  C     The following exchange was moved up to solve_for_pressure
422        _BEGIN_MASTER( myThid )  C     for compatibility with TAMC.
423         WRITE(0,'(A,I6,1PE30.14)') ' CG2D iters, err = ',  C     _EXCH_XY_R8(cg2d_x, myThid )
424       & actualIts, actualResidual  c     _BEGIN_MASTER( myThid )
425        _END_MASTER( )  c      WRITE(*,'(A,I6,1PE30.14)') ' CG2D iters, err = ',
426    c    & actualIts, actualResidual
427    c     _END_MASTER( )
428    
429    C--   Return parameters to caller
430          lastResidual=actualResidual
431          numIters=actualIts
432    
433  CcnhDebugStarts  CcnhDebugStarts
434  C     CALL PLOT_FIELD_XYRL( cg2d_x, 'CALC_MOM_RHS CG2D_X' , 1, myThid )  C     CALL PLOT_FIELD_XYRL( cg2d_x, 'CALC_MOM_RHS CG2D_X' , 1, myThid )
# Line 390  C        ENDDO Line 452  C        ENDDO
452  C       ENDDO  C       ENDDO
453  C      ENDDO  C      ENDDO
454  C     ENDDO  C     ENDDO
455  C     errBuf(1,myThid)    = err  C     _GLOBAL_SUM_R8( err   , myThid )
456  C     _GLOBAL_SUM_R8( errBuf    , err   , myThid )  C     write(*,*) 'cg2d: Ax - b = ',SQRT(err)
 C     err    = errBuf(1,1)  
 C     write(0,*) 'cg2d: Ax - b = ',SQRT(err)  
457  CcnhDebugEnds  CcnhDebugEnds
458    
459          RETURN
460        END        END
Legend:
Removed from v.1.18  
changed lines
  Added in v.1.33
  ViewVC Help
Powered by ViewVC 1.1.22