Changeset 868
- Timestamp:
- 2008-03-14T19:53:00+01:00 (16 years ago)
- Location:
- trunk/NEMO/LIM_SRC_3
- Files:
-
- 8 edited
-
limadv.F90 (modified) (11 diffs)
-
limdyn.F90 (modified) (4 diffs)
-
limhdf.F90 (modified) (2 diffs)
-
limitd_me.F90 (modified) (18 diffs)
-
limitd_th.F90 (modified) (1 diff)
-
limrhg.F90 (modified) (23 diffs)
-
limtrp.F90 (modified) (1 diff)
-
limvar.F90 (modified) (11 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/LIM_SRC_3/limadv.F90
r834 r868 25 25 PUBLIC lim_adv_x ! called by lim_trp 26 26 PUBLIC lim_adv_y ! called by lim_trp 27 28 !! * Substitutions 29 # include "vectopt_loop_substitute.h90" 27 30 28 31 !! * Module variables … … 110 113 111 114 ! Calculate fluxes and moments between boxes i<-->i+1 112 DO jj = 2, jpjm1! Flux from i to i+1 WHEN u GT 0115 DO jj = 1, jpj ! Flux from i to i+1 WHEN u GT 0 113 116 !i bug DO ji = 1, jpim1 114 117 !i DO jj = 1, jpj ! Flux from i to i+1 WHEN u GT 0 … … 138 141 END DO 139 142 140 DO jj = 2, jpjm1 ! Flux from i+1 to i when u LT 0.141 !i DO jj = 1, jpjm1! Flux from i+1 to i when u LT 0.142 DO ji = 1, jpim1143 DO jj = 1, jpjm1 ! Flux from i+1 to i when u LT 0. 144 !i DO jj = 1, fs_jpjm1 ! Flux from i+1 to i when u LT 0. 145 DO ji = 1, fs_jpim1 143 146 zalf = MAX( rzero, -put(ji,jj) ) * zrdt * e2u(ji,jj) / psm(ji+1,jj) 144 147 zalg (ji,jj) = zalf … … 159 162 160 163 DO jj = 2, jpjm1 ! Readjust moments remaining in the box. 161 DO ji = 2, jpim1164 DO ji = fs_2, jpi 162 165 zbt = zbet(ji-1,jj) 163 166 zbt1 = 1.0 - zbet(ji-1,jj) … … 174 177 ! Put the temporary moments into appropriate neighboring boxes. 175 178 DO jj = 2, jpjm1 ! Flux from i to i+1 IF u GT 0. 176 DO ji = 2,jpim1179 DO ji = fs_2, fs_jpim1 177 180 zbt = zbet(ji-1,jj) 178 181 zbt1 = 1.0 - zbet(ji-1,jj) … … 195 198 196 199 DO jj = 2, jpjm1 ! Flux from i+1 to i IF u LT 0. 197 DO ji = 2,jpim1200 DO ji = fs_2, fs_jpim1 198 201 zbt = zbet(ji,jj) 199 202 zbt1 = 1.0 - zbet(ji,jj) … … 305 308 !!bug DO jj = 2, jpjm1 306 309 DO jj = 1, jpj 307 DO ji = 2, jpim1310 DO ji = 1, jpi 308 311 !!bug DO ji = 1, jpim1 309 312 ! Flux from j to j+1 WHEN v GT 0 … … 333 336 334 337 DO jj = 1, jpjm1 ! Flux from j+1 to j when v LT 0. 335 DO ji = 2, jpim1338 DO ji = 1, jpi 336 339 !i DO jj = 1, jpjm1 ! Flux from j+1 to j when v LT 0. 337 340 !i DO ji = 2, jpim1 … … 354 357 355 358 ! Readjust moments remaining in the box. 356 DO jj = 2, jpj m1357 DO ji = 2, jpim1359 DO jj = 2, jpj 360 DO ji = 1, jpi 358 361 zbt = zbet(ji,jj-1) 359 362 zbt1 = ( 1.0 - zbet(ji,jj-1) ) … … 370 373 ! Put the temporary moments into appropriate neighboring boxes. 371 374 DO jj = 2, jpjm1 ! Flux from j to j+1 IF v GT 0. 372 DO ji = 2, jpim1375 DO ji = 1, jpi 373 376 zbt = zbet(ji,jj-1) 374 377 zbt1 = ( 1.0 - zbet(ji,jj-1) ) … … 395 398 396 399 DO jj = 2, jpjm1 ! Flux from j+1 to j IF v LT 0. 397 DO ji = 2, jpim1400 DO ji = 1, jpi 398 401 zbt = zbet(ji,jj) 399 402 zbt1 = ( 1.0 - zbet(ji,jj) ) -
trunk/NEMO/LIM_SRC_3/limdyn.F90
r867 r868 33 33 PUBLIC lim_dyn ! routine called by ice_step 34 34 35 !! * Substitutions 36 # include "vectopt_loop_substitute.h90" 37 35 38 !! * Module variables 36 39 REAL(wp) :: rone = 1.e0 ! constant value … … 172 175 zsang = SIGN(1.e0, gphif(1,jj-1) ) * sangvg 173 176 174 DO ji = 2,jpim1177 DO ji = fs_2, fs_jpim1 175 178 ! computation of wind stress over ocean in X and Y direction 176 179 #if defined key_coupled && defined key_lim_cp1 … … 216 219 ! computation of friction velocity 217 220 DO jj = 2, jpjm1 218 DO ji = 2,jpim1221 DO ji = fs_2, fs_jpim1 219 222 220 223 zu_ice = u_ice(ji,jj) - u_io(ji,jj) … … 245 248 ! virer ca (key_lim_cp1) 246 249 DO jj = 2, jpjm1 247 DO ji = 2,jpim1250 DO ji = fs_2, fs_jpim1 248 251 #if defined key_coupled && defined key_lim_cp1 249 252 tio_u(ji,jj) = - ( gtaux(ji ,jj ) + gtaux(ji-1,jj ) & -
trunk/NEMO/LIM_SRC_3/limhdf.F90
r834 r868 147 147 END DO 148 148 149 ! lateral boundary condition on zrlx 150 CALL lbc_lnk( zrlx, 'T', 1. ) 151 149 152 ! convergence test 150 153 zconv = 0.e0 151 154 DO jj = 2, jpjm1 152 DO ji = 2,jpim1155 DO ji = fs_2, fs_jpim1 153 156 zconv = MAX( zconv, ABS( zrlx(ji,jj) - ptab(ji,jj) ) ) 154 157 END DO … … 156 159 IF( lk_mpp ) CALL mpp_max( zconv ) ! max over the global domain 157 160 158 DO jj = 2, jpjm1159 DO ji = 2 , jpim1161 DO jj = 1, jpj 162 DO ji = 1 , jpi 160 163 ptab(ji,jj) = zrlx(ji,jj) 161 164 END DO 162 165 END DO 163 166 164 ! lateral boundary condition on ptab165 CALL lbc_lnk( ptab, 'T', 1. )166 167 ! !========================== 167 168 END DO ! end of sub-time step loop 168 169 ! !========================== 169 170 ptab(:,:) = ptab(:,:)171 170 172 171 IF(ln_ctl) THEN -
trunk/NEMO/LIM_SRC_3/limitd_me.F90
r867 r868 977 977 IF ( raftswi .EQ. 1 ) THEN 978 978 979 DO jl = 1, jpl 980 DO jj = 1, jpj 981 DO ji = 1, jpi 982 IF ( aridge(ji,jj,jl) + araft(ji,jj,jl) - athorn(ji,jj,jl) .GT. & 983 epsi11 ) THEN 984 WRITE(numout,*) ' ALERTE 96 : wrong participation function ... ' 985 WRITE(numout,*) ' ji, jj, jl : ', ji, jj, jl 986 WRITE(numout,*) ' lat, lon : ', gphit(ji,jj), glamt(ji,jj) 987 WRITE(numout,*) ' aridge : ', aridge(ji,jj,1:jpl) 988 WRITE(numout,*) ' araft : ', araft(ji,jj,1:jpl) 989 WRITE(numout,*) ' athorn : ', athorn(ji,jj,1:jpl) 990 ENDIF 979 IF( MAXVAL(aridge + araft - athorn(:,:,1:jpl)) .GT. epsi11 ) THEN 980 DO jl = 1, jpl 981 DO jj = 1, jpj 982 DO ji = 1, jpi 983 IF ( aridge(ji,jj,jl) + araft(ji,jj,jl) - athorn(ji,jj,jl) .GT. & 984 epsi11 ) THEN 985 WRITE(numout,*) ' ALERTE 96 : wrong participation function ... ' 986 WRITE(numout,*) ' ji, jj, jl : ', ji, jj, jl 987 WRITE(numout,*) ' lat, lon : ', gphit(ji,jj), glamt(ji,jj) 988 WRITE(numout,*) ' aridge : ', aridge(ji,jj,1:jpl) 989 WRITE(numout,*) ' araft : ', araft(ji,jj,1:jpl) 990 WRITE(numout,*) ' athorn : ', athorn(ji,jj,1:jpl) 991 ENDIF 992 END DO 991 993 END DO 992 994 END DO 993 END DO995 ENDIF 994 996 995 997 ENDIF … … 1238 1240 vsnon_init(ji,jj,jl) = v_s(ji,jj,jl) 1239 1241 1240 esnon_init(ji,jj,jl) = e_s(ji,jj,1,jl)1241 1242 smv_i_init(ji,jj,jl) = smv_i(ji,jj,jl) 1242 1243 oa_i_init (ji,jj,jl) = oa_i(ji,jj,jl) … … 1244 1245 END DO ! jj 1245 1246 END DO !jl 1247 1248 esnon_init(:,:,:) = e_s(:,:,1,:) 1246 1249 1247 1250 DO jl = 1, jpl … … 1283 1286 large_afrft = .false. 1284 1287 1288 !CDIR NODEP 1285 1289 DO ij = 1, icells 1286 1290 ji = indxi(ij) … … 1426 1430 !-------------------------------------------------------------------- 1427 1431 DO jk = 1, nlay_i 1432 !CDIR NODEP 1428 1433 DO ij = 1, icells 1429 1434 ji = indxi(ij) … … 1468 1473 IF ( con_i ) THEN 1469 1474 DO jk = 1, nlay_i 1475 !CDIR NODEP 1470 1476 DO ij = 1, icells 1471 1477 ji = indxi(ij) … … 1478 1484 1479 1485 IF (large_afrac) THEN ! there is a bug 1486 !CDIR NODEP 1480 1487 DO ij = 1, icells 1481 1488 ji = indxi(ij) … … 1490 1497 ENDIF ! large_afrac 1491 1498 IF (large_afrft) THEN ! there is a bug 1499 !CDIR NODEP 1492 1500 DO ij = 1, icells 1493 1501 ji = indxi(ij) … … 1508 1516 DO jl2 = ice_cat_bounds(1,1), ice_cat_bounds(1,2) 1509 1517 ! over categories to which ridged ice is transferred 1518 !CDIR NODEP 1510 1519 DO ij = 1, icells 1511 1520 ji = indxi(ij) … … 1542 1551 ! Transfer ice energy to category jl2 by ridging 1543 1552 DO jk = 1, nlay_i 1553 !CDIR NODEP 1544 1554 DO ij = 1, icells 1545 1555 ji = indxi(ij) … … 1555 1565 DO jl2 = ice_cat_bounds(1,1), ice_cat_bounds(1,2) 1556 1566 1567 !CDIR NODEP 1557 1568 DO ij = 1, icells 1558 1569 ji = indxi(ij) … … 1579 1590 ! Transfer rafted ice energy to category jl2 1580 1591 DO jk = 1, nlay_i 1592 !CDIR NODEP 1581 1593 DO ij = 1, icells 1582 1594 ji = indxi(ij) … … 1729 1741 jl, & ! ice category index 1730 1742 jk, & ! ice layer index 1731 icells, & ! number of cells with ice to zap 1732 ij ! combined i/j horizontal index 1733 1734 INTEGER, DIMENSION(1:(jpi+1)*(jpj+1)) :: & 1735 indxi, & ! compressed indices for i/j directions 1736 indxj 1743 ! ij, & ! combined i/j horizontal index 1744 icells ! number of cells with ice to zap 1745 1746 ! INTEGER, DIMENSION(1:(jpi+1)*(jpj+1)) :: & 1747 ! indxi, & ! compressed indices for i/j directions 1748 ! indxj 1749 1750 INTEGER, DIMENSION(jpi,jpj) :: zmask 1751 1737 1752 1738 1753 REAL(wp) :: & … … 1745 1760 ! Abort model in case of negative area. 1746 1761 !----------------------------------------------------------------- 1747 1748 icells = 0 1749 DO jj = 1, jpj 1762 IF( MAXVAL(a_i(:,:,jl)) .LT. -epsi11 ) THEN 1763 DO jj = 1, jpj 1764 DO ji = 1, jpi 1765 IF ( a_i(ji,jj,jl) .LT. -epsi11 ) THEN 1766 WRITE (numout,*) ' ALERTE 98 ' 1767 WRITE (numout,*) ' Negative ice area: ji, jj, jl: ', ji, jj,jl 1768 WRITE (numout,*) ' a_i ', a_i(ji,jj,jl) 1769 ENDIF 1770 END DO 1771 END DO 1772 ENDIF 1773 1774 icells = 0 1775 zmask = 0.e0 1776 DO jj = 1, jpj 1750 1777 DO ji = 1, jpi 1751 IF ( a_i(ji,jj,jl) .LT. -1.0e-11 ) THEN 1752 WRITE (numout,*) ' ALERTE 98 ' 1753 WRITE (numout,*) ' Negative ice area: ji, jj, jl: ', ji, jj,jl 1754 WRITE (numout,*) ' a_i ', a_i(ji,jj,jl) 1755 ELSEIF ( ( a_i(ji,jj,jl) .GE. -epsi11 .AND. a_i(ji,jj,jl) .LT. 0.0) & 1778 IF ( ( a_i(ji,jj,jl) .GE. -epsi11 .AND. a_i(ji,jj,jl) .LT. 0.0) & 1756 1779 .OR. & 1757 1780 ( a_i(ji,jj,jl) .GT. 0.0 .AND. a_i(ji,jj,jl) .LE. 1.0e-11 ) & … … 1761 1784 .OR. & 1762 1785 ( v_i(ji,jj,jl) .GT. 0.0 .AND. v_i(ji,jj,jl) .LT. 1.e-12 ) ) THEN 1763 icells = icells + 1 1764 indxi(icells) = ji 1765 indxj(icells) = jj 1786 zmask(ji,jj) = 1 1766 1787 ENDIF 1767 1788 END DO 1768 1789 END DO 1769 WRITE(numout,*) icells, ' cells of ice zapped in the ocean '1790 WRITE(numout,*) SUM(zmask), ' cells of ice zapped in the ocean ' 1770 1791 1771 1792 !----------------------------------------------------------------- … … 1774 1795 1775 1796 DO jk = 1, nlay_i 1776 DO ij = 1, icells 1777 ji = indxi(ij) 1778 jj = indxj(ij) 1797 DO jj = 1 , jpj 1798 DO ji = 1 , jpi 1779 1799 1780 1800 xtmp = e_i(ji,jj,jk,jl) / area(ji,jj) / rdt_ice 1781 1801 xtmp = xtmp * unit_fac 1782 1802 ! fheat_res(ji,jj) = fheat_res(ji,jj) - xtmp 1783 e_i(ji,jj,jk,jl) = 0.01784 1785 END DO ! ij1803 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * ( 1 - zmask(ji,jj) ) 1804 END DO ! ji 1805 END DO ! jj 1786 1806 END DO ! jk 1787 1807 1788 DO ij = 1, icells 1789 ji = indxi(ij) 1790 jj = indxj(ij) 1808 DO jj = 1 , jpj 1809 DO ji = 1 , jpi 1791 1810 1792 1811 !----------------------------------------------------------------- … … 1803 1822 ! fheat_res(ji,jj) = fheat_res(ji,jj) - xtmp 1804 1823 1805 t_s(ji,jj,1,jl) = rtt 1824 xtmp = ( rhosn*cpic*( rtt-t_s(ji,jj,1,jl) ) + rhosn*lfus ) / rdt_ice !RB ??????? 1825 1826 t_s(ji,jj,1,jl) = rtt * zmask(ji,jj) + t_s(ji,jj,1,jl) * ( 1 - zmask(ji,jj) ) 1806 1827 1807 1828 !----------------------------------------------------------------- … … 1822 1843 ! fsalt_hist(i,j) = fsalt_hist(i,j) + xtmp 1823 1844 1824 ato_i(ji,jj) = ato_i(ji,jj) + a_i(ji,jj,jl) 1825 a_i(ji,jj,jl) = 0.0 1826 v_i(ji,jj,jl) = 0.0 1827 v_s(ji,jj,jl) = 0.0 1828 t_su(ji,jj,jl) = t_bo(ji,jj) 1829 oa_i(ji,jj,jl) = 0.0 1830 smv_i(ji,jj,jl) = 0.0 1831 1832 END DO ! ij 1845 ato_i(ji,jj) = a_i(ji,jj,jl) * zmask(ji,jj) + ato_i(ji,jj) 1846 a_i(ji,jj,jl) = a_i(ji,jj,jl) * ( 1 - zmask(ji,jj) ) 1847 v_i(ji,jj,jl) = v_i(ji,jj,jl) * ( 1 - zmask(ji,jj) ) 1848 v_s(ji,jj,jl) = v_s(ji,jj,jl) * ( 1 - zmask(ji,jj) ) 1849 t_su(ji,jj,jl) = t_su(ji,jj,jl) * (1 -zmask(ji,jj) ) + t_bo(ji,jj) * zmask(ji,jj) 1850 oa_i(ji,jj,jl) = oa_i(ji,jj,jl) * ( 1 - zmask(ji,jj) ) 1851 smv_i(ji,jj,jl) = smv_i(ji,jj,jl) * ( 1 - zmask(ji,jj) ) 1852 1853 END DO ! ji 1854 END DO ! jj 1833 1855 1834 1856 END DO ! jl -
trunk/NEMO/LIM_SRC_3/limitd_th.F90
r867 r868 985 985 986 986 DO jk = 1, nlay_i 987 !CDIR NODEP 987 988 DO ji = 1, nbrem 988 989 zji = nind_i(ji) -
trunk/NEMO/LIM_SRC_3/limrhg.F90
r866 r868 31 31 PUBLIC lim_rhg ! routine called by lim_dyn 32 32 33 !! * Substitutions 34 # include "vectopt_loop_substitute.h90" 35 33 36 !! * Module variables 34 37 REAL(wp) :: & ! constant values … … 111 114 112 115 INTEGER :: & 113 iter, jter!: temporary integers116 jter !: temporary integers 114 117 115 118 CHARACTER (len=50) :: charout … … 182 185 ! 183 186 ! Put every vector to 0 184 zpresh(:,:) = 0.0 ; zpreshc(:,:) = 0.0 ; zfrld1(:,:) = 0.0 185 zmass1(:,:) = 0.0 ; zcorl1(:,:) = 0.0 ; zcorl2(:,:) = 0.0 186 za1ct(:,:) = 0.0 ; za2ct(:,:) = 0.0 187 zc1(:,:) = 0.0 ; zusw(:,:) = 0.0 188 u_oce1(:,:) = 0.0 ; v_oce1(:,:) = 0.0 ; u_oce2(:,:) = 0.0 189 v_oce2(:,:) = 0.0 ; u_ice2(:,:) = 0.0 ; v_ice1(:,:) = 0.0 190 zf1(:,:) = 0.0 ; zf2(:,:) = 0.0 187 zpresh(:,:) = 0.0 ; zc1(:,:) = 0.0 188 zpreshc(:,:) = 0.0 189 u_ice2(:,:) = 0.0 ; v_ice1(:,:) = 0.0 191 190 zdd(:,:) = 0.0 ; zdt(:,:) = 0.0 ; zds(:,:) = 0.0 192 deltat(:,:) = 0.0 ; deltac(:,:) = 0.0193 zpresh(:,:) = 0.0194 191 195 192 ! Ice strength on T-points … … 197 194 198 195 ! Ice mass and temp variables 199 DO jj = k_j1 , k_jpj-1 196 !CDIR NOVERRCHK 197 DO jj = k_j1 , k_jpj 198 !CDIR NOVERRCHK 200 199 DO ji = 1 , jpi 201 200 zc1(ji,jj) = tms(ji,jj) * ( rhosn * vt_s(ji,jj) + rhoic * vt_i(ji,jj) ) … … 207 206 END DO 208 207 209 CALL lbc_lnk( zc1(:,:), 'T', 1. )210 CALL lbc_lnk( zpresh(:,:), 'T', 1. )211 212 CALL lbc_lnk( tmi(:,:), 'T', 1. )213 214 208 ! Ice strength on grid cell corners (zpreshc) 215 209 ! needed for calculation of shear stress 210 !CDIR NOVERRCHK 216 211 DO jj = k_j1+1, k_jpj-1 217 DO ji = 2, jpim1 212 !CDIR NOVERRCHK 213 DO ji = 2, jpim1 !RB caution no fs_ (ji+1,jj+1) 218 214 zstms = tms(ji+1,jj+1) * wght(ji+1,jj+1,2,2) + & 219 215 & tms(ji,jj+1) * wght(ji+1,jj+1,1,2) + & … … 249 245 250 246 DO jj = k_j1+1, k_jpj-1 251 DO ji = 2,jpim1247 DO ji = fs_2, fs_jpim1 252 248 253 249 zt11 = tms(ji,jj)*e1t(ji,jj) … … 276 272 277 273 ! Ocean has no slip boundary condition 278 ! GG bug279 ! v_oce1(ji,jj) = 0.5*( (v_io(ji,jj)+v_io(ji,jj-1))*e1t(ji+1,jj) &280 ! & +(v_io(ji+1,jj)+v_io(ji+1,jj-1))*e1t(ji,jj)) &281 ! & /(e1t(ji+1,jj)+e1t(ji,jj)) * tmu(ji,jj)282 274 v_oce1(ji,jj) = 0.5*( (v_io(ji,jj)+v_io(ji,jj-1))*e1t(ji,jj) & 283 275 & +(v_io(ji+1,jj)+v_io(ji+1,jj-1))*e1t(ji+1,jj)) & 284 276 & /(e1t(ji+1,jj)+e1t(ji,jj)) * tmu(ji,jj) 285 277 286 ! GG bug287 ! u_oce2(ji,jj) = 0.5*((u_io(ji,jj)+u_io(ji-1,jj))*e2t(ji,jj+1) &288 ! & +(u_io(ji,jj+1)+u_io(ji-1,jj+1))*e2t(ji,jj)) &289 ! & / (e2t(ji,jj+1)+e2t(ji,jj)) * tmv(ji,jj)290 278 u_oce2(ji,jj) = 0.5*((u_io(ji,jj)+u_io(ji-1,jj))*e2t(ji,jj) & 291 279 & +(u_io(ji,jj+1)+u_io(ji-1,jj+1))*e2t(ji,jj+1)) & … … 330 318 zs12(:,:) = stress12_i(:,:) 331 319 332 v_ice1(:,:) = 0.0333 u_ice2(:,:) = 0.0334 335 zdd(:,:) = 0.0336 zdt(:,:) = 0.0337 zds(:,:) = 0.0338 deltat(:,:) = 0.0339 320 !----------------------! 340 DO iter = 1 , nevp ! loop over iter !321 DO jter = 1 , nevp ! loop over jter ! 341 322 !----------------------! 342 323 DO jj = k_j1, k_jpj-1 … … 346 327 347 328 DO jj = k_j1+1, k_jpj-1 348 DO ji = 2,jpim1329 DO ji = fs_2, fs_jpim1 349 330 350 331 ! … … 407 388 END DO 408 389 409 CALL lbc_lnk( zdd(:,:), 'T', 1. ) 410 CALL lbc_lnk( zdt(:,:), 'T', 1. ) 411 CALL lbc_lnk( zds(:,:), 'F', 1. ) 412 390 !CDIR NOVERRCHK 413 391 DO jj = k_j1+1, k_jpj-1 414 DO ji = 2, jpim1 392 !CDIR NOVERRCHK 393 DO ji = fs_2, fs_jpim1 415 394 416 395 !- Calculate Delta at centre of grid cells … … 446 425 CALL lbc_lnk( zs2(:,:), 'T', 1. ) 447 426 427 !CDIR NOVERRCHK 448 428 DO jj = k_j1+1, k_jpj-1 449 DO ji = 2, jpim1 429 !CDIR NOVERRCHK 430 DO ji = fs_2, fs_jpim1 450 431 !- Calculate Delta on corners 451 432 zddc = ( ( v_ice1(ji,jj+1)/e1u(ji,jj+1) & … … 482 463 ! Ice internal stresses (Appendix C of Hunke and Dukowicz, 2002) 483 464 DO jj = k_j1+1, k_jpj-1 484 DO ji = 2,jpim1465 DO ji = fs_2, fs_jpim1 485 466 !- contribution of zs1, zs2 and zs12 to zf1 486 467 zf1(ji,jj) = 0.5*( (zs1(ji+1,jj)-zs1(ji,jj))*e2u(ji,jj) & … … 501 482 ! Both the Coriolis term and the ice-ocean drag are solved semi-implicitly. 502 483 ! 503 IF (MOD(iter,2).eq.0) THEN 504 484 IF (MOD(jter,2).eq.0) THEN 485 486 !CDIR NOVERRCHK 505 487 DO jj = k_j1+1, k_jpj-1 506 DO ji = 2, jpim1 488 !CDIR NOVERRCHK 489 DO ji = fs_2, fs_jpim1 507 490 zmask = (1.0-MAX(rzero,SIGN(rone,-zmass1(ji,jj))))*tmu(ji,jj) 508 491 zsang = SIGN ( 1.0 , fcor(ji,jj) ) * sangvg … … 510 493 511 494 ! SB modif because ocean has no slip boundary condition 512 ! GG bug513 ! zv_ice1 = 0.5*( (v_ice(ji,jj)+v_ice(ji,jj-1))*e1t(ji+1,jj) &514 ! & +(v_ice(ji+1,jj)+v_ice(ji+1,jj-1))*e1t(ji,jj)) &515 ! & /(e1t(ji+1,jj)+e1t(ji,jj)) * tmu(ji,jj)516 495 zv_ice1 = 0.5*( (v_ice(ji,jj)+v_ice(ji,jj-1))*e1t(ji,jj) & 517 496 & +(v_ice(ji+1,jj)+v_ice(ji+1,jj-1))*e1t(ji+1,jj)) & … … 530 509 CALL lbc_lnk( u_ice(:,:), 'U', -1. ) 531 510 511 !CDIR NOVERRCHK 532 512 DO jj = k_j1+1, k_jpj-1 533 DO ji = 2, jpim1 534 535 zmask = (1.0-max(rzero,sign(rone,-zmass2(ji,jj))))*tmv(ji,jj) 536 zsang = sign(1.0,fcor(ji,jj))*sangvg 513 !CDIR NOVERRCHK 514 DO ji = fs_2, fs_jpim1 515 516 zmask = (1.0-MAX(rzero,SIGN(rone,-zmass2(ji,jj))))*tmv(ji,jj) 517 zsang = SIGN(1.0,fcor(ji,jj))*sangvg 537 518 z0 = zmass2(ji,jj)/dtevp 538 519 ! SB modif because ocean has no slip boundary condition 539 ! GG bug540 ! zu_ice2 = 0.5*( (u_ice(ji,jj)+u_ice(ji-1,jj))*e2t(ji,jj+1) &541 ! & + (u_ice(ji,jj+1)+u_ice(ji-1,jj+1))*e2t(ji,jj)) &542 ! & /(e2t(ji,jj+1)+e2t(ji,jj)) * tmv(ji,jj)543 520 zu_ice2 = 0.5*( (u_ice(ji,jj)+u_ice(ji-1,jj))*e2t(ji,jj) & 544 521 & + (u_ice(ji,jj+1)+u_ice(ji-1,jj+1))*e2t(ji,jj+1)) & … … 558 535 559 536 ELSE 537 !CDIR NOVERRCHK 560 538 DO jj = k_j1+1, k_jpj-1 561 DO ji = 2, jpim1 562 zmask = (1.0-max(rzero,sign(rone,-zmass2(ji,jj))))*tmv(ji,jj) 563 zsang = sign(1.0,fcor(ji,jj))*sangvg 539 !CDIR NOVERRCHK 540 DO ji = fs_2, fs_jpim1 541 zmask = (1.0-MAX(rzero,SIGN(rone,-zmass2(ji,jj))))*tmv(ji,jj) 542 zsang = SIGN(1.0,fcor(ji,jj))*sangvg 564 543 z0 = zmass2(ji,jj)/dtevp 565 544 ! SB modif because ocean has no slip boundary condition 566 ! GG Bug567 ! zu_ice2 = 0.5*( (u_ice(ji,jj)+u_ice(ji-1,jj))*e2t(ji,jj+1) &568 ! & +(u_ice(ji,jj+1)+u_ice(ji-1,jj+1))*e2t(ji,jj)) &569 ! & /(e2t(ji,jj+1)+e2t(ji,jj)) * tmv(ji,jj)570 545 zu_ice2 = 0.5*( (u_ice(ji,jj)+u_ice(ji-1,jj))*e2t(ji,jj) & 571 546 & +(u_ice(ji,jj+1)+u_ice(ji-1,jj+1))*e2t(ji,jj+1)) & … … 585 560 CALL lbc_lnk( v_ice(:,:), 'V', -1. ) 586 561 562 !CDIR NOVERRCHK 587 563 DO jj = k_j1+1, k_jpj-1 588 DO ji = 2, jpim1 564 !CDIR NOVERRCHK 565 DO ji = fs_2, fs_jpim1 589 566 zmask = (1.0-MAX(rzero,SIGN(rone,-zmass1(ji,jj))))*tmu(ji,jj) 590 567 zsang = SIGN(1.0,fcor(ji,jj))*sangvg … … 613 590 ENDIF 614 591 615 !--- Convergence test. 616 DO jj = k_j1+1 , k_jpj-1 617 zresr(:,jj) = MAX( ABS( u_ice(:,jj) - zu_ice(:,jj) ) , & 618 ABS( v_ice(:,jj) - zv_ice(:,jj) ) ) 619 END DO 620 zresm = MAXVAL( zresr( 1:jpi , k_j1+1:k_jpj-1 ) ) 592 IF(ln_ctl) THEN 593 !--- Convergence test. 594 DO jj = k_j1+1 , k_jpj-1 595 zresr(:,jj) = MAX( ABS( u_ice(:,jj) - zu_ice(:,jj) ) , & 596 ABS( v_ice(:,jj) - zv_ice(:,jj) ) ) 597 END DO 598 zresm = MAXVAL( zresr( 1:jpi , k_j1+1:k_jpj-1 ) ) 599 IF( lk_mpp ) CALL mpp_max( zresm ) ! max over the global domain 600 ENDIF 621 601 622 602 ! ! ==================== ! 623 END DO ! end loop over iter !603 END DO ! end loop over jter ! 624 604 ! ! ==================== ! 625 605 … … 632 612 ! ocean velocity, 633 613 ! This prevents high velocity when ice is thin 614 !CDIR NOVERRCHK 634 615 DO jj = k_j1+1, k_jpj-1 635 DO ji = 2, jpim1 616 !CDIR NOVERRCHK 617 DO ji = fs_2, fs_jpim1 636 618 zindb = MAX( 0.0, SIGN( 1.0, at_i(ji,jj) - 1.0e-6 ) ) 637 619 zdummy = zindb * vt_i(ji,jj) / MAX(at_i(ji,jj) , 1.0e-06 ) … … 643 625 END DO 644 626 627 DO jj = k_j1+1, k_jpj-1 628 DO ji = fs_2, fs_jpim1 629 zindb = MAX( 0.0, SIGN( 1.0, at_i(ji,jj) - 1.0e-6 ) ) 630 zdummy = zindb * vt_i(ji,jj) / MAX(at_i(ji,jj) , 1.0e-06 ) 631 IF ( zdummy .LE. 5.0e-2 ) THEN 632 v_ice1(ji,jj) = 0.5*( (v_ice(ji,jj)+v_ice(ji,jj-1))*e1t(ji+1,jj) & 633 & +(v_ice(ji+1,jj)+v_ice(ji+1,jj-1))*e1t(ji,jj)) & 634 & /(e1t(ji+1,jj)+e1t(ji,jj)) * tmu(ji,jj) 635 636 u_ice2(ji,jj) = 0.5*( (u_ice(ji,jj)+u_ice(ji-1,jj))*e2t(ji,jj+1) & 637 & +(u_ice(ji,jj+1)+u_ice(ji-1,jj+1))*e2t(ji,jj)) & 638 & /(e2t(ji,jj+1)+e2t(ji,jj)) * tmv(ji,jj) 639 ENDIF ! zdummy 640 END DO 641 END DO 642 643 CALL lbc_lnk( u_ice2(:,:), 'U', -1. ) 644 CALL lbc_lnk( v_ice1(:,:), 'V', -1. ) 645 645 646 ! Recompute delta, shear and div, inputs for mechanical redistribution 647 !CDIR NOVERRCHK 646 648 DO jj = k_j1+1, k_jpj-1 647 DO ji = 2, jpim1 649 !CDIR NOVERRCHK 650 DO ji = fs_2, fs_jpim1 648 651 !- zdd(:,:), zdt(:,:): divergence and tension at centre 649 652 !- zds(:,:): shear on northeast corner of grid cells … … 680 683 & * tmi(ji+1,jj) * tmi(ji+1,jj+1) 681 684 682 !- Calculate Delta at centre of grid cells683 v_ice1(ji,jj) = 0.5*( (v_ice(ji,jj)+v_ice(ji,jj-1))*e1t(ji+1,jj) &684 & +(v_ice(ji+1,jj)+v_ice(ji+1,jj-1))*e1t(ji,jj)) &685 & /(e1t(ji+1,jj)+e1t(ji,jj)) * tmu(ji,jj)686 687 u_ice2(ji,jj) = 0.5*( (u_ice(ji,jj)+u_ice(ji-1,jj))*e2t(ji,jj+1) &688 & +(u_ice(ji,jj+1)+u_ice(ji-1,jj+1))*e2t(ji,jj)) &689 & /(e2t(ji,jj+1)+e2t(ji,jj)) * tmv(ji,jj)690 691 685 zdst = ( e2u( ji , jj ) * v_ice1(ji,jj) & 692 686 & - e2u( ji-1, jj ) * v_ice1(ji-1,jj) & … … 710 704 ! 711 705 ! * Invariants of the stress tensor are required for limitd_me 712 ! * Store the stress tensor for the next time step713 706 ! accelerates convergence and improves stability 714 707 DO jj = k_j1+1, k_jpj-1 715 DO ji = 2, jpim1 716 divu_i(ji,jj) = zdd(ji,jj) 717 delta_i(ji,jj) = deltat (ji,jj) 718 shear_i(ji,jj) = zds (ji,jj) 719 stress1_i(ji,jj) = zs1(ji,jj) 720 stress2_i(ji,jj) = zs2(ji,jj) 721 stress12_i(ji,jj) = zs12(ji,jj) 708 DO ji = fs_2, fs_jpim1 709 divu_i (ji,jj) = zdd (ji,jj) 710 delta_i(ji,jj) = deltat(ji,jj) 711 shear_i(ji,jj) = zds (ji,jj) 722 712 END DO 723 713 END DO 724 714 725 715 ! Lateral boundary condition 726 CALL lbc_lnk( divu_i (:,:), 'T', 1. )716 CALL lbc_lnk( divu_i (:,:), 'T', 1. ) 727 717 CALL lbc_lnk( delta_i(:,:), 'T', 1. ) 728 718 CALL lbc_lnk( shear_i(:,:), 'F', 1. ) 729 CALL lbc_lnk( stress1_i(:,:), 'T', 1. ) 730 CALL lbc_lnk( stress2_i(:,:), 'T', 1. ) 731 CALL lbc_lnk( stress12_i(:,:), 'F', 1. ) 719 720 ! * Store the stress tensor for the next time step 721 stress1_i (:,:) = zs1 (:,:) 722 stress2_i (:,:) = zs2 (:,:) 723 stress12_i(:,:) = zs12(:,:) 732 724 733 725 ! … … 738 730 ! print the residual for convergence 739 731 IF(ln_ctl) THEN 740 WRITE(charout,FMT="('lim_rhg : res =',D23.16, ' iter =',I4)") zresm, iter732 WRITE(charout,FMT="('lim_rhg : res =',D23.16, ' iter =',I4)") zresm, jter 741 733 CALL prt_ctl_info(charout) 742 734 CALL prt_ctl(tab2d_1=u_ice, clinfo1=' lim_rhg : u_ice :', tab2d_2=v_ice, clinfo2=' v_ice :') -
trunk/NEMO/LIM_SRC_3/limtrp.F90
r867 r868 127 127 zvbord = 1.0 + ( 1.0 - bound ) 128 128 DO jj = 1, jpjm1 129 DO ji = 1, jpim1129 DO ji = 1, fs_jpim1 130 130 zui_u(ji,jj) = u_ice(ji,jj) 131 131 zvi_v(ji,jj) = v_ice(ji,jj) -
trunk/NEMO/LIM_SRC_3/limvar.F90
r834 r868 236 236 ! Ice thickness, snow thickness, ice salinity, ice age 237 237 !------------------------------------------------------- 238 !CDIR NOVERRCHK 238 239 DO jl = 1, jpl 240 !CDIR NOVERRCHK 239 241 DO jj = 1, jpj 242 !CDIR NOVERRCHK 240 243 DO ji = 1, jpi 241 244 zindb = 1.0-MAX(0.0,SIGN(1.0,- a_i(ji,jj,jl))) !0 if no ice and 1 if yes … … 249 252 IF ( ( num_sal .EQ. 2 ) .OR. ( num_sal .EQ. 4 ) )THEN 250 253 254 !CDIR NOVERRCHK 251 255 DO jl = 1, jpl 256 !CDIR NOVERRCHK 252 257 DO jj = 1, jpj 258 !CDIR NOVERRCHK 253 259 DO ji = 1, jpi 254 260 zindb = 1.0-MAX(0.0,SIGN(1.0,-a_i(ji,jj,jl))) !0 if no ice and 1 if yes … … 266 272 ! Ice temperatures 267 273 !------------------- 274 !CDIR NOVERRCHK 268 275 DO jl = 1, jpl 276 !CDIR NOVERRCHK 269 277 DO jk = 1, nlay_i 278 !CDIR NOVERRCHK 270 279 DO jj = 1, jpj 280 !CDIR NOVERRCHK 271 281 DO ji = 1, jpi 272 282 !Energy of melting q(S,T) [J.m-3] … … 298 308 zfac1 = 1. / ( rhosn * cpic ) 299 309 zfac2 = lfus / cpic 310 !CDIR NOVERRCHK 300 311 DO jl = 1, jpl 312 !CDIR NOVERRCHK 301 313 DO jk = 1, nlay_s 314 !CDIR NOVERRCHK 302 315 DO jj = 1, jpj 316 !CDIR NOVERRCHK 303 317 DO ji = 1, jpi 304 318 !Energy of melting q(S,T) [J.m-3] … … 321 335 !------------------- 322 336 tm_i(:,:) = 0.0 337 !CDIR NOVERRCHK 323 338 DO jl = 1, jpl 339 !CDIR NOVERRCHK 324 340 DO jk = 1, nlay_i 341 !CDIR NOVERRCHK 325 342 DO jj = 1, jpj 343 !CDIR NOVERRCHK 326 344 DO ji = 1, jpi 327 345 zindb = 1.0-MAX(0.0,SIGN(1.0,-a_i(ji,jj,jl))) … … 462 480 zalpha(:,:,:) = 0.0 463 481 482 !CDIR NOVERRCHK 464 483 DO jl = 1, jpl 484 !CDIR NOVERRCHK 465 485 DO jj = 1, jpj 486 !CDIR NOVERRCHK 466 487 DO ji = 1, jpi 467 488 ! zind0 = 1 if sm_i le s_i_0 and 0 otherwise … … 511 532 sm_i(:,:,:) = 2.30 512 533 534 !CDIR NOVERRCHK 513 535 DO jl = 1, jpl 536 !CDIR NOVERRCHK 514 537 DO jk = 1, nlay_i 538 !CDIR NOVERRCHK 515 539 DO jj = 1, jpj 540 !CDIR NOVERRCHK 516 541 DO ji = 1, jpi 517 542 zargtemp = ( jk - 0.5 ) / nlay_i … … 567 592 zeps = 1.0e-13 568 593 bv_i(:,:) = 0.0 594 !CDIR NOVERRCHK 569 595 DO jl = 1, jpl 596 !CDIR NOVERRCHK 570 597 DO jk = 1, nlay_i 598 !CDIR NOVERRCHK 571 599 DO jj = 1, jpj 600 !CDIR NOVERRCHK 572 601 DO ji = 1, jpi 573 602 zindb = 1.0-MAX(0.0,SIGN(1.0,-a_i(ji,jj,jl))) !0 if no ice and 1 if yes … … 639 668 ! Slope of the linear profile zs_zero 640 669 !------------------------------------- 670 !CDIR NOVERRCHK 641 671 DO ji = kideb, kiut 642 672 z_slope_s(ji) = 2.0 * sm_i_b(ji) / MAX( 0.01 & … … 650 680 dummy_fac2 = 1. / nlay_i 651 681 682 !CDIR NOVERRCHK 652 683 DO jk = 1, nlay_i 684 !CDIR NOVERRCHK 653 685 DO ji = kideb, kiut 654 686 zji = MOD( npb(ji) - 1, jpi ) + 1 … … 688 720 sm_i_b(:) = 2.30 689 721 722 !CDIR NOVERRCHK 690 723 DO ji = kideb, kiut 724 !CDIR NOVERRCHK 691 725 DO jk = 1, nlay_i 692 726 zargtemp = ( jk - 0.5 ) / nlay_i
Note: See TracChangeset
for help on using the changeset viewer.
