/[MITgcm]/MITgcm/pkg/seaice/seaice_jfnk.F
ViewVC logotype

Diff of /MITgcm/pkg/seaice/seaice_jfnk.F

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

revision 1.29 by mlosch, Wed Jan 27 14:03:34 2016 UTC revision 1.31 by mlosch, Thu Jun 8 15:33:50 2017 UTC
# Line 63  C     !FUNCTIONS: Line 63  C     !FUNCTIONS:
63    
64  C     !LOCAL VARIABLES:  C     !LOCAL VARIABLES:
65  C     === Local variables ===  C     === Local variables ===
66  C     i,j,k,bi,bj :: loop indices  C     i,j,bi,bj :: loop indices
67        INTEGER i,j,k,bi,bj        INTEGER i,j,bi,bj
68  C     loop indices  C     loop indices
69        INTEGER newtonIter        INTEGER newtonIter
70        INTEGER krylovIter, krylovFails        INTEGER krylovIter, krylovFails
# Line 181  CML        ENDIF Line 181  CML        ENDIF
181         ENDDO         ENDDO
182        ENDDO        ENDDO
183  C     Start nonlinear Newton iteration: outer loop iteration  C     Start nonlinear Newton iteration: outer loop iteration
184        DO WHILE ( newtonIter.LT.SEAICEnewtonIterMax .AND.        DO WHILE ( newtonIter.LT.SEAICEnonLinIterMax .AND.
185       &     .NOT.JFNKconverged )       &     .NOT.JFNKconverged )
186         newtonIter = newtonIter + 1         newtonIter = newtonIter + 1
187  C     Compute initial residual F(u), (includes computation of global  C     Compute initial residual F(u), (includes computation of global
# Line 228  C     in that routine Line 228  C     in that routine
228         iCode         = 0         iCode         = 0
229    
230         JFNKconverged = JFNKresidual.LT.JFNKtol         JFNKconverged = JFNKresidual.LT.JFNKtol
231         &      .OR.JFNKresidual.EQ.0. _d 0
232    
233  C     do Krylov loop only if convergence is not reached  C     do Krylov loop only if convergence is not reached
234    
# Line 240  C     start Krylov iteration (FGMRES) Line 241  C     start Krylov iteration (FGMRES)
241  C     map first guess sol; it is zero because the solution is a correction  C     map first guess sol; it is zero because the solution is a correction
242         CALL SEAICE_MAP2VEC(nVec,duIce,dvIce,sol,.TRUE.,myThid)         CALL SEAICE_MAP2VEC(nVec,duIce,dvIce,sol,.TRUE.,myThid)
243  C     map rhs and change its sign because we are solving J*u = -F  C     map rhs and change its sign because we are solving J*u = -F
244          CALL SEAICE_MAP2VEC(nVec,-uIceRes,-vIceRes,rhs,.TRUE.,myThid)          CALL SEAICE_MAP2VEC(nVec,uIceRes,vIceRes,rhs,.TRUE.,myThid)
245            DO bj=myByLo(myThid),myByHi(myThid)
246             DO bi=myBxLo(myThid),myBxHi(myThid)
247              DO j=1,nVec
248               rhs(j,bi,bj) = - rhs(j,bi,bj)
249              ENDDO
250             ENDDO
251            ENDDO
252          DO WHILE ( .NOT.krylovConverged )          DO WHILE ( .NOT.krylovConverged )
253  C     solution vector sol = du/vIce  C     solution vector sol = du/vIce
254  C     residual vector (rhs) Fu = u/vIceRes  C     residual vector (rhs) Fu = u/vIceRes
# Line 253  C     because du/vIce = 0 Line 261  C     because du/vIce = 0
261  C  C
262           CALL SEAICE_FGMRES (nVec,im,rhs,sol,ifgmres,krylovIter,           CALL SEAICE_FGMRES (nVec,im,rhs,sol,ifgmres,krylovIter,
263       U        vv,w,wk1,wk2,       U        vv,w,wk1,wk2,
264       I        FGMRESeps,SEAICEkrylovIterMax,iOutFGMRES,       I        FGMRESeps,SEAICElinearIterMax,iOutFGMRES,
265       U        iCode,       U        iCode,
266       I        myThid)       I        myThid)
267  C      C    
# Line 272  C     or product of matrix (Jacobian) ti Line 280  C     or product of matrix (Jacobian) ti
280  C     iteration  C     iteration
281           IF (iCode.EQ.1) THEN           IF (iCode.EQ.1) THEN
282  C     Call preconditioner  C     Call preconditioner
283            IF ( SOLV_MAX_ITERS .GT. 0 )            IF ( SEAICEpreconLinIter .GT. 0 )
284       &         CALL SEAICE_PRECONDITIONER(       &         CALL SEAICE_PRECONDITIONER(
285       U         duIce, dvIce,       U         duIce, dvIce,
286       I         zetaPre, etaPre, etaZpre, zetaZpre, dwatPre,       I         zetaPre, etaPre, etaZpre, zetaZpre, dwatPre,
# Line 292  C     some output diagnostics Line 300  C     some output diagnostics
300          IF ( debugLevel.GE.debLevA ) THEN          IF ( debugLevel.GE.debLevA ) THEN
301           _BEGIN_MASTER( myThid )           _BEGIN_MASTER( myThid )
302           totalNewtonItersLoc =           totalNewtonItersLoc =
303       &        SEAICEnewtonIterMax*(myIter-nIter0)+newtonIter       &        SEAICEnonLinIterMax*(myIter-nIter0)+newtonIter
304           WRITE(msgBuf,'(2A,2(1XI6),2E12.5)')           WRITE(msgBuf,'(2A,2(1XI6),2E12.5)')
305       &        ' S/R SEAICE_JFNK: Newton iterate / total, ',       &        ' S/R SEAICE_JFNK: Newton iterate / total, ',
306       &        'JFNKgamma_lin, initial norm = ',       &        'JFNKgamma_lin, initial norm = ',
# Line 308  C     some output diagnostics Line 316  C     some output diagnostics
316       &        SQUEEZE_RIGHT, myThid )       &        SQUEEZE_RIGHT, myThid )
317           _END_MASTER( myThid )           _END_MASTER( myThid )
318          ENDIF          ENDIF
319          IF ( krylovIter.EQ.SEAICEkrylovIterMax ) THEN          IF ( krylovIter.EQ.SEAICElinearIterMax ) THEN
320           krylovFails = krylovFails + 1           krylovFails = krylovFails + 1
321          ENDIF          ENDIF
322  C     Set the stopping criterion for the Newton iteration and the  C     Set the stopping criterion for the Newton iteration and the
323  C     criterion for the transition from accurate to approximate FGMRES  C     criterion for the transition from accurate to approximate FGMRES
324          IF ( newtonIter .EQ. 1 ) THEN          IF ( newtonIter .EQ. 1 ) THEN
325           JFNKtol=JFNKgamma_nonlin*JFNKresidual           JFNKtol=SEAICEnonLinTol*JFNKresidual
326           IF ( JFNKres_tFac .NE. UNSET_RL )           IF ( JFNKres_tFac .NE. UNSET_RL )
327       &        JFNKres_t = JFNKresidual * JFNKres_tFac       &        JFNKres_t = JFNKresidual * JFNKres_tFac
328          ENDIF          ENDIF
# Line 353  C     Count iterations Line 361  C     Count iterations
361         totalKrylovIters   = totalKrylovIters + totalKrylovItersLoc         totalKrylovIters   = totalKrylovIters + totalKrylovItersLoc
362  C     Record failure  C     Record failure
363         totalKrylovFails   = totalKrylovFails + krylovFails         totalKrylovFails   = totalKrylovFails + krylovFails
364         IF ( newtonIter .EQ. SEAICEnewtonIterMax ) THEN         IF ( newtonIter .EQ. SEAICEnonLinIterMax ) THEN
365          totalNewtonFails = totalNewtonFails + 1          totalNewtonFails = totalNewtonFails + 1
366         ENDIF         ENDIF
367        ENDIF        ENDIF
# Line 427  C     reset and start again Line 435  C     reset and start again
435    
436  C     Print more debugging information  C     Print more debugging information
437        IF ( debugLevel.GE.debLevA ) THEN        IF ( debugLevel.GE.debLevA ) THEN
438         IF ( newtonIter .EQ. SEAICEnewtonIterMax ) THEN         IF ( newtonIter .EQ. SEAICEnonLinIterMax ) THEN
439          _BEGIN_MASTER( myThid )          _BEGIN_MASTER( myThid )
440          WRITE(msgBuf,'(A,I10)')          WRITE(msgBuf,'(A,I10)')
441       &       ' S/R SEAICE_JFNK: JFNK did not converge in timestep ',       &       ' S/R SEAICE_JFNK: JFNK did not converge in timestep ',
Legend:
Removed from v.1.29  
changed lines
  Added in v.1.31
  ViewVC Help
Powered by ViewVC 1.1.22