| 73 |
_RL JFNKresidual |
_RL JFNKresidual |
| 74 |
_RL JFNKresidualKm1 |
_RL JFNKresidualKm1 |
| 75 |
C parameters to compute convergence criterion |
C parameters to compute convergence criterion |
| 76 |
_RL phi_e, alp_e, JFNKgamma_lin |
_RL JFNKgamma_lin |
| 77 |
_RL FGMRESeps |
_RL FGMRESeps |
| 78 |
_RL JFNKtol |
_RL JFNKtol |
| 79 |
|
C backward differences extrapolation factors |
| 80 |
|
_RL bdfFac, bdfAlpha |
| 81 |
|
C |
| 82 |
_RL recip_deltaT |
_RL recip_deltaT |
| 83 |
LOGICAL JFNKconverged, krylovConverged |
LOGICAL JFNKconverged, krylovConverged |
| 84 |
LOGICAL writeNow |
LOGICAL writeNow |
| 87 |
C u/vIceRes :: residual of sea-ice momentum equations |
C u/vIceRes :: residual of sea-ice momentum equations |
| 88 |
_RL uIceRes(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
_RL uIceRes(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
| 89 |
_RL vIceRes(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
_RL vIceRes(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
| 90 |
|
C extra time level required for backward difference time stepping |
| 91 |
|
_RL duIcNm1(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
| 92 |
|
_RL dvIcNm1(1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
| 93 |
C du/vIce :: ice velocity increment to be added to u/vIce |
C du/vIce :: ice velocity increment to be added to u/vIce |
| 94 |
_RL duIce (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
_RL duIce (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
| 95 |
_RL dvIce (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
_RL dvIce (1-OLx:sNx+OLx,1-OLy:sNy+OLy,nSx,nSy) |
| 118 |
& DIFFERENT_MULTIPLE( SEAICE_monFreq, myTime, deltaTClock ) ) |
& DIFFERENT_MULTIPLE( SEAICE_monFreq, myTime, deltaTClock ) ) |
| 119 |
& iOutFGMRES=1 |
& iOutFGMRES=1 |
| 120 |
|
|
| 121 |
|
C backward difference extrapolation factors |
| 122 |
|
bdfFac = 0. _d 0 |
| 123 |
|
IF ( SEAICEuseBDF2 ) THEN |
| 124 |
|
IF ( myIter.EQ.nIter0 .AND. SEAICEmomStartBDF.EQ.0 ) THEN |
| 125 |
|
bdfFac = 0. _d 0 |
| 126 |
|
ELSE |
| 127 |
|
bdfFac = 0.5 _d 0 |
| 128 |
|
ENDIF |
| 129 |
|
ENDIF |
| 130 |
|
bdfAlpha = 1. _d 0 + bdfFac |
| 131 |
|
|
| 132 |
DO bj=myByLo(myThid),myByHi(myThid) |
DO bj=myByLo(myThid),myByHi(myThid) |
| 133 |
DO bi=myBxLo(myThid),myBxHi(myThid) |
DO bi=myBxLo(myThid),myBxHi(myThid) |
| 134 |
DO J=1-OLy,sNy+OLy |
DO J=1-OLy,sNy+OLy |
| 137 |
vIceRes(I,J,bi,bj) = 0. _d 0 |
vIceRes(I,J,bi,bj) = 0. _d 0 |
| 138 |
duIce (I,J,bi,bj) = 0. _d 0 |
duIce (I,J,bi,bj) = 0. _d 0 |
| 139 |
dvIce (I,J,bi,bj) = 0. _d 0 |
dvIce (I,J,bi,bj) = 0. _d 0 |
| 140 |
|
ENDDO |
| 141 |
|
ENDDO |
| 142 |
|
C cycle ice velocities |
| 143 |
|
DO J=1-OLy,sNy+OLy |
| 144 |
|
DO I=1-OLx,sNx+OLx |
| 145 |
|
duIcNm1(I,J,bi,bj) = uIce(I,J,bi,bj) * bdfAlpha |
| 146 |
|
& + ( uIce(I,J,bi,bj) - uIceNm1(I,J,bi,bj) ) * bdfFac |
| 147 |
|
dvIcNm1(I,J,bi,bj) = vIce(I,J,bi,bj) * bdfAlpha |
| 148 |
|
& + ( vIce(I,J,bi,bj) - vIceNm1(I,J,bi,bj) ) * bdfFac |
| 149 |
uIceNm1(I,J,bi,bj) = uIce(I,J,bi,bj) |
uIceNm1(I,J,bi,bj) = uIce(I,J,bi,bj) |
| 150 |
vIceNm1(I,J,bi,bj) = vIce(I,J,bi,bj) |
vIceNm1(I,J,bi,bj) = vIce(I,J,bi,bj) |
| 151 |
ENDDO |
ENDDO |
| 152 |
ENDDO |
ENDDO |
| 153 |
|
C As long as IMEX is not properly implemented leave this commented out |
| 154 |
|
CML IF ( .NOT.SEAICEuseIMEX ) THEN |
| 155 |
C Compute things that do no change during the Newton iteration: |
C Compute things that do no change during the Newton iteration: |
| 156 |
C sea-surface tilt and wind stress: |
C sea-surface tilt and wind stress: |
| 157 |
C FORCEX/Y0 - mass*(u/vIceNm1)/deltaT |
C FORCEX/Y0 - mass*(1.5*u/vIceNm1+0.5*(u/vIceNm1-u/vIceNm2))/deltaT |
| 158 |
DO J=1-OLy,sNy+OLy |
DO J=1-OLy,sNy+OLy |
| 159 |
DO I=1-OLx,sNx+OLx |
DO I=1-OLx,sNx+OLx |
| 160 |
FORCEX(I,J,bi,bj) = FORCEX0(I,J,bi,bj) |
FORCEX(I,J,bi,bj) = FORCEX0(I,J,bi,bj) |
| 161 |
& + seaiceMassU(I,J,bi,bj)*uIceNm1(I,J,bi,bj)*recip_deltaT |
& + seaiceMassU(I,J,bi,bj)*duIcNm1(I,J,bi,bj)*recip_deltaT |
| 162 |
FORCEY(I,J,bi,bj) = FORCEY0(I,J,bi,bj) |
FORCEY(I,J,bi,bj) = FORCEY0(I,J,bi,bj) |
| 163 |
& + seaiceMassV(I,J,bi,bj)*vIceNm1(I,J,bi,bj)*recip_deltaT |
& + seaiceMassV(I,J,bi,bj)*dvIcNm1(I,J,bi,bj)*recip_deltaT |
| 164 |
ENDDO |
ENDDO |
| 165 |
ENDDO |
ENDDO |
| 166 |
|
CML ENDIF |
| 167 |
ENDDO |
ENDDO |
| 168 |
ENDDO |
ENDDO |
| 169 |
C Start nonlinear Newton iteration: outer loop iteration |
C Start nonlinear Newton iteration: outer loop iteration |
| 195 |
JFNKgamma_lin = JFNKgamma_lin_max |
JFNKgamma_lin = JFNKgamma_lin_max |
| 196 |
IF ( newtonIter.GT.1.AND.newtonIter.LE.SEAICE_JFNK_tolIter |
IF ( newtonIter.GT.1.AND.newtonIter.LE.SEAICE_JFNK_tolIter |
| 197 |
& .AND.JFNKresidual.LT.JFNKres_t ) THEN |
& .AND.JFNKresidual.LT.JFNKres_t ) THEN |
| 198 |
C Eisenstat, 1996, equ.(2.6) |
C Eisenstat and Walker (1996), eq.(2.6) |
| 199 |
phi_e = 1. _d 0 |
JFNKgamma_lin = SEAICE_JFNKphi |
| 200 |
alp_e = 1. _d 0 |
& *( JFNKresidual/JFNKresidualKm1 )**SEAICE_JFNKalpha |
|
JFNKgamma_lin = phi_e*( JFNKresidual/JFNKresidualKm1 )**alp_e |
|
| 201 |
JFNKgamma_lin = min(JFNKgamma_lin_max, JFNKgamma_lin) |
JFNKgamma_lin = min(JFNKgamma_lin_max, JFNKgamma_lin) |
| 202 |
JFNKgamma_lin = max(JFNKgamma_lin_min, JFNKgamma_lin) |
JFNKgamma_lin = max(JFNKgamma_lin_min, JFNKgamma_lin) |
| 203 |
ENDIF |
ENDIF |