Changeset 10994
- Timestamp:
- 2019-05-17T15:08:34+02:00 (6 years ago)
- Location:
- NEMO/trunk/src/ICE
- Files:
-
- 9 edited
-
icecor.F90 (modified) (2 diffs)
-
icectl.F90 (modified) (3 diffs)
-
icedyn.F90 (modified) (2 diffs)
-
icedyn_rdgrft.F90 (modified) (13 diffs)
-
iceitd.F90 (modified) (14 diffs)
-
icestp.F90 (modified) (1 diff)
-
icethd.F90 (modified) (3 diffs)
-
icethd_do.F90 (modified) (2 diffs)
-
icevar.F90 (modified) (5 diffs)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/src/ICE/icecor.F90
r10425 r10994 66 66 WRITE(numout,*) '~~~~~~~' 67 67 ENDIF 68 !69 68 ! !----------------------------------------------------- 70 ! ! ice thickness must exceed himin (for ice diff)!69 ! ! ice thickness must exceed himin (for temp. diff.) ! 71 70 ! !----------------------------------------------------- 72 71 WHERE( a_i(:,:,:) >= epsi20 ) ; h_i(:,:,:) = v_i(:,:,:) / a_i(:,:,:) … … 97 96 END DO 98 97 ENDIF 99 100 98 ! !----------------------------------------------------- 101 99 ! ! Rebin categories with thickness out of bounds ! -
NEMO/trunk/src/ICE/icectl.F90
r10581 r10994 69 69 !! 70 70 REAL(wp) :: zv, zs, zt, zfs, zfv, zft 71 REAL(wp) :: zvmin, zamin, zamax 71 REAL(wp) :: zvmin, zamin, zamax, zeimin, zesmin, zsmin 72 72 REAL(wp) :: zvtrp, zetrp 73 73 REAL(wp) :: zarea, zv_sill, zs_sill, zt_sill … … 141 141 zetrp = glob_sum( 'icectl', ( diag_trp_ei + diag_trp_es ) * e1e2t ) * zconv 142 142 143 zvmin = glob_min( 'icectl', v_i ) 144 zamax = glob_max( 'icectl', SUM( a_i, dim=3 ) ) 145 zamin = glob_min( 'icectl', a_i ) 143 zamax = glob_max( 'icectl', SUM( a_i, dim=3 ) ) 144 zvmin = glob_min( 'icectl', v_i ) 145 zamin = glob_min( 'icectl', a_i ) 146 zsmin = glob_min( 'icectl', sv_i ) 147 zeimin = glob_min( 'icectl', SUM( e_i, dim=3 ) ) 148 zesmin = glob_min( 'icectl', SUM( e_s, dim=3 ) ) 146 149 147 150 ! set threshold values and calculate the ice area (+epsi10 to set a threshold > 0 when there is no ice) … … 152 155 153 156 IF(lwp) THEN 154 IF ( ABS( zv ) > zv_sill ) WRITE(numout,*) 'violation volume [Mt/day] (',cd_routine,') = ',zv 155 IF ( ABS( zs ) > zs_sill ) WRITE(numout,*) 'violation saline [psu*Mt/day] (',cd_routine,') = ',zs 156 IF ( ABS( zt ) > zt_sill ) WRITE(numout,*) 'violation enthalpy [GW] (',cd_routine,') = ',zt 157 IF ( zvmin < -epsi10 ) WRITE(numout,*) 'violation v_i<0 [m] (',cd_routine,') = ',zvmin 158 IF ( zamax > MAX( rn_amax_n, rn_amax_s ) + epsi10 & 159 & .AND. cd_routine /= 'icedyn_adv' .AND. cd_routine /= 'icedyn_rdgrft' .AND. cd_routine /= 'Hbig' ) & 160 & WRITE(numout,*) 'violation a_i>amax (',cd_routine,') = ',zamax 161 IF ( zamin < -epsi10 ) WRITE(numout,*) 'violation a_i<0 (',cd_routine,') = ',zamin 162 !clem: the following check fails when using UMx advection scheme (see comments in icedyn_adv.F90) 157 ! check conservation issues 158 IF ( ABS( zv ) > zv_sill ) WRITE(numout,*) 'violation volume [Mt/day] (',cd_routine,') = ',zv 159 IF ( ABS( zs ) > zs_sill ) WRITE(numout,*) 'violation saline [psu*Mt/day] (',cd_routine,') = ',zs 160 IF ( ABS( zt ) > zt_sill ) WRITE(numout,*) 'violation enthalpy [GW] (',cd_routine,') = ',zt 161 ! check maximum ice concentration 162 IF ( zamax > MAX( rn_amax_n,rn_amax_s)+epsi10 .AND. cd_routine /= 'icedyn_adv' .AND. cd_routine /= 'icedyn_rdgrft' ) & 163 & WRITE(numout,*) 'violation a_i>amax (',cd_routine,') = ',zamax 164 ! check negative values 165 IF ( zvmin < 0. ) WRITE(numout,*) 'violation v_i<0 [m] (',cd_routine,') = ',zvmin 166 IF ( zamin < 0. ) WRITE(numout,*) 'violation a_i<0 (',cd_routine,') = ',zamin 167 IF ( zsmin < 0. ) WRITE(numout,*) 'violation s_i<0 (',cd_routine,') = ',zsmin 168 IF ( zeimin < 0. ) WRITE(numout,*) 'violation e_i<0 (',cd_routine,') = ',zeimin 169 IF ( zesmin < 0. ) WRITE(numout,*) 'violation e_s<0 (',cd_routine,') = ',zesmin 170 !clem: the following check fails (I think...) 163 171 ! IF ( ABS(zvtrp ) > zv_sill .AND. cd_routine == 'icedyn_adv' ) THEN 164 172 ! WRITE(numout,*) 'violation vtrp [Mt/day] (',cd_routine,') = ',zvtrp -
NEMO/trunk/src/ICE/icedyn.F90
r10911 r10994 118 118 CALL ice_dyn_adv ( kt ) ! -- advection of ice 119 119 CALL Hpiling ! -- simple pile-up (replaces ridging/rafting) 120 CALL ice_var_zapsmall ! -- zap small areas 120 121 ! 121 122 CASE ( np_dynADV1D ) !== pure advection ==! (1D) … … 189 190 ! controls 190 191 IF( ln_icediachk ) CALL ice_cons_hsm(0, 'Hpiling', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) ! conservation 191 !192 CALL ice_var_zapsmall !-- zap small areas193 192 ! 194 193 at_i(:,:) = SUM( a_i(:,:,:), dim=3 ) -
NEMO/trunk/src/ICE/icedyn_rdgrft.F90
r10945 r10994 142 142 INTEGER, PARAMETER :: jp_itermax = 20 143 143 !!------------------------------------------------------------------- 144 ! clem: The redistribution of ice between categories can lead to small negative values (as for the remapping in ice_itd_rem)145 ! likely due to truncation error ( i.e. 1. - 1. /= 0 )146 ! I do not think it should be a concern since small areas and volumes are erased (in ice_var_zapsmall.F90)147 148 144 ! controls 149 145 IF( ln_timing ) CALL timing_start('icedyn_rdgrft') ! timing … … 156 152 ENDIF 157 153 158 CALL ice_var_zapsmall ! Zero out categories with very small areas159 160 154 !-------------------------------- 161 155 ! 0) Identify grid cells with ice 162 156 !-------------------------------- 157 at_i(:,:) = SUM( a_i, dim=3 ) 158 ! 163 159 npti = 0 ; nptidx(:) = 0 164 160 ipti = 0 ; iptidx(:) = 0 165 161 DO jj = 1, jpj 166 162 DO ji = 1, jpi 167 IF ( at_i(ji,jj) > 0._wp) THEN163 IF ( at_i(ji,jj) > epsi10 ) THEN 168 164 npti = npti + 1 169 165 nptidx( npti ) = (jj - 1) * jpi + ji … … 178 174 179 175 ! just needed here 180 CALL tab_2d_1d( npti, nptidx(1:npti), zdivu (1:npti), divu_i(:,:))181 CALL tab_2d_1d( npti, nptidx(1:npti), zdelt (1:npti), delta_i(:,:))176 CALL tab_2d_1d( npti, nptidx(1:npti), zdivu (1:npti) , divu_i ) 177 CALL tab_2d_1d( npti, nptidx(1:npti), zdelt (1:npti) , delta_i ) 182 178 ! needed here and in the iteration loop 183 CALL tab_3d_2d( npti, nptidx(1:npti), a_i_2d (1:npti,1:jpl), a_i (:,:,:))184 CALL tab_3d_2d( npti, nptidx(1:npti), v_i_2d (1:npti,1:jpl), v_i (:,:,:))185 CALL tab_2d_1d( npti, nptidx(1:npti), ato_i_1d(1:npti) , ato_i (:,:))179 CALL tab_3d_2d( npti, nptidx(1:npti), a_i_2d (1:npti,1:jpl), a_i ) 180 CALL tab_3d_2d( npti, nptidx(1:npti), v_i_2d (1:npti,1:jpl), v_i ) 181 CALL tab_2d_1d( npti, nptidx(1:npti), ato_i_1d(1:npti) , ato_i ) 186 182 187 183 DO ji = 1, npti … … 310 306 311 307 ! ! Ice thickness needed for rafting 312 WHERE( pa_i(1:npti,:) > 0._wp) ; zhi(1:npti,:) = pv_i(1:npti,:) / pa_i(1:npti,:)313 ELSEWHERE ; zhi(1:npti,:) = 0._wp308 WHERE( pa_i(1:npti,:) > epsi20 ) ; zhi(1:npti,:) = pv_i(1:npti,:) / pa_i(1:npti,:) 309 ELSEWHERE ; zhi(1:npti,:) = 0._wp 314 310 END WHERE 315 311 … … 329 325 zasum(1:npti) = pato_i(1:npti) + SUM( pa_i(1:npti,:), dim=2 ) 330 326 ! 331 WHERE( zasum(1:npti) > 0._wp) ; z1_asum(1:npti) = 1._wp / zasum(1:npti)332 ELSEWHERE ; z1_asum(1:npti) = 0._wp327 WHERE( zasum(1:npti) > epsi20 ) ; z1_asum(1:npti) = 1._wp / zasum(1:npti) 328 ELSEWHERE ; z1_asum(1:npti) = 0._wp 333 329 END WHERE 334 330 ! … … 455 451 ! Based on the ITD of ridging and ridged ice, convert the net closing rate to a gross closing rate. 456 452 ! NOTE: 0 < aksum <= 1 457 WHERE( zaksum(1:npti) > 0._wp) ; closing_gross(1:npti) = pclosing_net(1:npti) / zaksum(1:npti)458 ELSEWHERE ; closing_gross(1:npti) = 0._wp453 WHERE( zaksum(1:npti) > epsi20 ) ; closing_gross(1:npti) = pclosing_net(1:npti) / zaksum(1:npti) 454 ELSEWHERE ; closing_gross(1:npti) = 0._wp 459 455 END WHERE 460 456 … … 466 462 DO ji = 1, npti 467 463 zfac = apartf(ji,jl) * closing_gross(ji) * rdt_ice 468 IF( zfac > pa_i(ji,jl) ) THEN464 IF( zfac > pa_i(ji,jl) .AND. apartf(ji,jl) /= 0._wp ) THEN 469 465 closing_gross(ji) = pa_i(ji,jl) / apartf(ji,jl) * r1_rdtice 470 466 ENDIF … … 510 506 REAL(wp), DIMENSION(jpij) :: zswitch, fvol ! new ridge volume going to jl2 511 507 REAL(wp), DIMENSION(jpij) :: z1_ai ! 1 / a 508 REAL(wp), DIMENSION(jpij) :: zvti ! sum(v_i) 512 509 ! 513 510 REAL(wp), DIMENSION(jpij,nlay_s) :: esrft ! snow energy of rafting ice … … 518 515 INTEGER , DIMENSION(jpij) :: itest_rdg, itest_rft ! test for conservation 519 516 !!------------------------------------------------------------------- 520 517 ! 518 zvti(1:npti) = SUM( v_i_2d(1:npti,:), dim=2 ) ! total ice volume 519 ! 521 520 ! 1) Change in open water area due to closing and opening 522 521 !-------------------------------------------------------- … … 535 534 IF( apartf(ji,jl1) > 0._wp .AND. closing_gross(ji) > 0._wp ) THEN ! only if ice is ridging 536 535 537 z1_ai(ji) = 1._wp / a_i_2d(ji,jl1) 538 536 IF( a_i_2d(ji,jl1) > epsi20 ) THEN ; z1_ai(ji) = 1._wp / a_i_2d(ji,jl1) 537 ELSE ; z1_ai(ji) = 0._wp 538 ENDIF 539 539 540 ! area of ridging / rafting ice (airdg1) and of new ridge (airdg2) 540 541 airdg1 = aridge(ji,jl1) * closing_gross(ji) * rdt_ice … … 549 550 550 551 ! volume and enthalpy (J/m2, >0) of seawater trapped into ridges 551 vsw = v_i_2d(ji,jl1) * afrdg * rn_porordg 552 IF ( zvti(ji) <= 10. ) THEN ; vsw = v_i_2d(ji,jl1) * afrdg * rn_porordg ! v <= 10m then porosity = rn_porordg 553 ELSEIF( zvti(ji) >= 20. ) THEN ; vsw = 0._wp ! v >= 20m then porosity = 0 554 ELSE ; vsw = v_i_2d(ji,jl1) * afrdg * rn_porordg * MAX( 0._wp, 2._wp - 0.1_wp * zvti(ji) ) ! v > 10m and v < 20m then porosity = linear transition to 0 555 ENDIF 552 556 ersw(ji) = -rhoi * vsw * rcp * sst_1d(ji) ! clem: if sst>0, then ersw <0 (is that possible?) 553 557 554 558 ! volume etc of ridging / rafting ice and new ridges (vi, vs, sm, oi, es, ei) 555 559 virdg1 = v_i_2d (ji,jl1) * afrdg 556 virdg2(ji) = v_i_2d (ji,jl1) * afrdg * ( 1. + rn_porordg )560 virdg2(ji) = v_i_2d (ji,jl1) * afrdg + vsw 557 561 vsrdg(ji) = v_s_2d (ji,jl1) * afrdg 558 562 sirdg1 = sv_i_2d(ji,jl1) * afrdg … … 726 730 END DO ! jl1 727 731 ! 732 ! roundoff errors 733 !---------------- 728 734 ! In case ridging/rafting lead to very small negative values (sometimes it happens) 729 WHERE( a_i_2d(1:npti,:) < 0._wp ) a_i_2d(1:npti,:) = 0._wp 730 WHERE( v_i_2d(1:npti,:) < 0._wp ) v_i_2d(1:npti,:) = 0._wp 735 CALL ice_var_roundoff( a_i_2d, v_i_2d, v_s_2d, sv_i_2d, oa_i_2d, a_ip_2d, v_ip_2d, ze_s_2d, ze_i_2d ) 731 736 ! 732 737 END SUBROUTINE rdgrft_shift … … 854 859 CALL tab_2d_1d( npti, nptidx(1:npti), wfx_snw_dyn_1d(1:npti), wfx_snw_dyn(:,:) ) 855 860 CALL tab_2d_1d( npti, nptidx(1:npti), wfx_pnd_1d (1:npti), wfx_pnd (:,:) ) 856 861 ! 857 862 ! !---------------------! 858 863 CASE( 2 ) !== from 1D to 2D ==! -
NEMO/trunk/src/ICE/iceitd.F90
r10069 r10994 21 21 USE ice1D ! sea-ice: thermodynamic variables 22 22 USE ice ! sea-ice: variables 23 USE icevar ! sea-ice: operations 23 24 USE icectl ! sea-ice: conservation tests 24 25 USE icetab ! sea-ice: convert 1D<=>2D … … 91 92 ! 1) Identify grid cells with ice 92 93 !----------------------------------------------------------------------------------------------- 94 at_i(:,:) = SUM( a_i, dim=3 ) 95 ! 93 96 npti = 0 ; nptidx(:) = 0 94 97 DO jj = 1, jpj … … 249 252 ! --- g(h) for each thickness category --- ! 250 253 CALL itd_glinear( zhbnew(1:npti,jl-1), zhbnew(1:npti,jl), h_i_1d(1:npti) , a_i_1d(1:npti) , & ! in 251 & g0 (1:npti,jl ), g1 (1:npti,jl), hL (1:npti,jl), hR(1:npti,jl) ) ! out254 & g0 (1:npti,jl ), g1 (1:npti,jl), hL (1:npti,jl), hR (1:npti,jl) ) ! out 252 255 ! 253 256 END DO … … 389 392 REAL(wp), DIMENSION(:,:), INTENT(in) :: pdvice ! ice volume transferred across boundary 390 393 ! 391 INTEGER :: ji, j j, jl, jk! dummy loop indices392 INTEGER :: ii, ij, jl2, jl1! local integers394 INTEGER :: ji, jl, jk ! dummy loop indices 395 INTEGER :: jl2, jl1 ! local integers 393 396 REAL(wp) :: ztrans ! ice/snow transferred 394 REAL(wp), DIMENSION(jpij) :: zworka, zworkv ! workspace 395 REAL(wp), DIMENSION(jpij,jpl) :: zaTsfn ! - - 397 REAL(wp), DIMENSION(jpij) :: zworka, zworkv ! workspace 398 REAL(wp), DIMENSION(jpij,jpl) :: zaTsfn ! - - 399 REAL(wp), DIMENSION(jpij,nlay_i,jpl) :: ze_i_2d 400 REAL(wp), DIMENSION(jpij,nlay_s,jpl) :: ze_s_2d 396 401 !!------------------------------------------------------------------ 397 402 … … 405 410 CALL tab_3d_2d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip ) 406 411 CALL tab_3d_2d( npti, nptidx(1:npti), t_su_2d(1:npti,1:jpl), t_su ) 412 DO jl = 1, jpl 413 DO jk = 1, nlay_s 414 CALL tab_2d_1d( npti, nptidx(1:npti), ze_s_2d(1:npti,jk,jl), e_s(:,:,jk,jl) ) 415 END DO 416 DO jk = 1, nlay_i 417 CALL tab_2d_1d( npti, nptidx(1:npti), ze_i_2d(1:npti,jk,jl), e_i(:,:,jk,jl) ) 418 END DO 419 END DO 420 ! to correct roundoff errors on a_i 421 CALL tab_2d_1d( npti, nptidx(1:npti), rn_amax_1d(1:npti), rn_amax_2d ) 407 422 408 423 !---------------------------------------------------------------------------------------------- … … 435 450 ELSE ; zworka(ji) = 0._wp 436 451 ENDIF 437 !438 ! clem: The transfer between one category to another can lead to very small negative values (-1.e-20)439 ! because of truncation error ( i.e. 1. - 1. /= 0 )440 ! I do not think it should be a concern since small areas and volumes are erased (in ice_var_zapsmall.F90)441 452 ! 442 453 a_i_2d(ji,jl1) = a_i_2d(ji,jl1) - pdaice(ji,jl) ! Ice areas … … 476 487 ! 477 488 DO jk = 1, nlay_s !--- Snow heat content 478 !479 489 DO ji = 1, npti 480 ii = MOD( nptidx(ji) - 1, jpi ) + 1481 ij = ( nptidx(ji) - 1 ) / jpi + 1482 490 ! 483 491 jl1 = kdonor(ji,jl) … … 487 495 ELSE ; jl2 = jl 488 496 ENDIF 489 ! 490 ztrans = e_s(ii,ij,jk,jl1) * zworkv(ji) 491 e_s(ii,ij,jk,jl1) = e_s(ii,ij,jk,jl1) - ztrans 492 e_s(ii,ij,jk,jl2) = e_s(ii,ij,jk,jl2) + ztrans 497 ztrans = ze_s_2d(ji,jk,jl1) * zworkv(ji) 498 ze_s_2d(ji,jk,jl1) = ze_s_2d(ji,jk,jl1) - ztrans 499 ze_s_2d(ji,jk,jl2) = ze_s_2d(ji,jk,jl2) + ztrans 493 500 ENDIF 494 501 END DO … … 497 504 DO jk = 1, nlay_i !--- Ice heat content 498 505 DO ji = 1, npti 499 ii = MOD( nptidx(ji) - 1, jpi ) + 1500 ij = ( nptidx(ji) - 1 ) / jpi + 1501 506 ! 502 507 jl1 = kdonor(ji,jl) … … 506 511 ELSE ; jl2 = jl 507 512 ENDIF 508 ! 509 ztrans = e_i(ii,ij,jk,jl1) * zworkv(ji) 510 e_i(ii,ij,jk,jl1) = e_i(ii,ij,jk,jl1) - ztrans 511 e_i(ii,ij,jk,jl2) = e_i(ii,ij,jk,jl2) + ztrans 513 ztrans = ze_i_2d(ji,jk,jl1) * zworkv(ji) 514 ze_i_2d(ji,jk,jl1) = ze_i_2d(ji,jk,jl1) - ztrans 515 ze_i_2d(ji,jk,jl2) = ze_i_2d(ji,jk,jl2) + ztrans 512 516 ENDIF 513 517 END DO … … 515 519 ! 516 520 END DO ! boundaries, 1 to jpl-1 521 522 !------------------- 523 ! 3) roundoff errors 524 !------------------- 525 ! clem: The transfer between one category to another can lead to very small negative values (-1.e-20) 526 ! because of truncation error ( i.e. 1. - 1. /= 0 ) 527 CALL ice_var_roundoff( a_i_2d, v_i_2d, v_s_2d, sv_i_2d, oa_i_2d, a_ip_2d, v_ip_2d, ze_s_2d, ze_i_2d ) 528 529 ! at_i must be <= rn_amax 530 zworka(1:npti) = SUM( a_i_2d(1:npti,:), dim=2 ) 531 DO jl = 1, jpl 532 WHERE( zworka(1:npti) > rn_amax_1d(1:npti) ) & 533 & a_i_2d(1:npti,jl) = a_i_2d(1:npti,jl) * rn_amax_1d(1:npti) / zworka(1:npti) 534 END DO 517 535 518 536 !------------------------------------------------------------------------------- 519 ! 3) Update ice thickness and temperature537 ! 4) Update ice thickness and temperature 520 538 !------------------------------------------------------------------------------- 521 539 WHERE( a_i_2d(1:npti,:) >= epsi20 ) … … 536 554 CALL tab_2d_3d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip ) 537 555 CALL tab_2d_3d( npti, nptidx(1:npti), t_su_2d(1:npti,1:jpl), t_su ) 556 DO jl = 1, jpl 557 DO jk = 1, nlay_s 558 CALL tab_1d_2d( npti, nptidx(1:npti), ze_s_2d(1:npti,jk,jl), e_s(:,:,jk,jl) ) 559 END DO 560 DO jk = 1, nlay_i 561 CALL tab_1d_2d( npti, nptidx(1:npti), ze_i_2d(1:npti,jk,jl), e_i(:,:,jk,jl) ) 562 END DO 563 END DO 538 564 ! 539 565 END SUBROUTINE itd_shiftice … … 558 584 ! 559 585 IF( kt == nit000 .AND. lwp ) WRITE(numout,*) '-- ice_itd_reb: rebining ice thickness distribution' 586 ! 587 IF( ln_icediachk ) CALL ice_cons_hsm(0, 'iceitd_reb', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) 560 588 ! 561 589 jdonor(:,:) = 0 … … 635 663 END DO 636 664 ! 665 IF( ln_icediachk ) CALL ice_cons_hsm(1, 'iceitd_reb', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) 666 ! 637 667 END SUBROUTINE ice_itd_reb 638 668 -
NEMO/trunk/src/ICE/icestp.F90
r10932 r10994 323 323 WRITE(numout,*) ' maximum ice concentration for SH = ', rn_amax_s 324 324 ENDIF 325 ! !--- change max ice concentration for roundoff errors 326 rn_amax_n = MIN( rn_amax_n, 1._wp - epsi10 ) 327 rn_amax_s = MIN( rn_amax_s, 1._wp - epsi10 ) 325 328 ! !--- check consistency 326 329 IF ( jpl > 1 .AND. ln_virtual_itd ) THEN -
NEMO/trunk/src/ICE/icethd.F90
r10932 r10994 102 102 ENDIF 103 103 104 CALL ice_var_glo2eqv105 106 104 !---------------------------------------------! 107 105 ! computation of friction velocity at T points … … 162 160 qlead(ji,jj) = MIN( 0._wp , zqld - ( qsb_ice_bot(ji,jj) * at_i(ji,jj) * rdt_ice ) - zqfr ) 163 161 164 ! If there is ice and leads are warming, then transfer energy from the lead budget and use it for bottom melting 165 IF( zqld > 0._wp ) THEN 162 ! If there is ice and leads are warming => transfer energy from the lead budget and use it for bottom melting 163 ! If the grid cell is fully covered by ice (no leads) => transfer energy from the lead budget to the ice bottom budget 164 IF( ( zqld >= 0._wp .AND. at_i(ji,jj) > 0._wp ) .OR. at_i(ji,jj) >= (1._wp - epsi10) ) THEN 166 165 fhld (ji,jj) = rswitch * zqld * r1_rdtice / MAX( at_i(ji,jj), epsi10 ) ! divided by at_i since this is (re)multiplied by a_i in icethd_dh.F90 167 166 qlead(ji,jj) = 0._wp … … 242 241 ! 243 242 END DO 244 243 ! 245 244 IF( ln_icediachk ) CALL ice_cons_hsm(1, 'icethd', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) 246 !247 CALL ice_var_zapsmall ! --- remove very small ice concentration (<1e-10) --- !248 ! ! & make sure at_i=sum(a_i) & ato_i=1 where at_i=0249 245 ! 250 246 IF( jpl > 1 ) CALL ice_itd_rem( kt ) ! --- Transport ice between thickness categories --- ! -
NEMO/trunk/src/ICE/icethd_do.F90
r10425 r10994 114 114 IF( ln_icediachk ) CALL ice_cons_hsm( 0, 'icethd_do', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft ) 115 115 116 CALL ice_var_agg(1) 117 CALL ice_var_glo2eqv 118 116 at_i(:,:) = SUM( a_i, dim=3 ) 119 117 !------------------------------------------------------------------------------! 120 118 ! 1) Collection thickness of ice formed in leads and polynyas … … 319 317 320 318 ! --- lateral ice growth --- ! 321 ! If lateral ice growth gives an ice concentration gt 1, then319 ! If lateral ice growth gives an ice concentration > amax, then 322 320 ! we keep the excessive volume in memory and attribute it later to bottom accretion 323 321 DO ji = 1, npti 324 IF ( za_newice(ji) > ( rn_amax_1d(ji) - at_i_1d(ji) ) ) THEN325 zda_res(ji) = za_newice(ji) - (rn_amax_1d(ji) - at_i_1d(ji) )322 IF ( za_newice(ji) > MAX( 0._wp, rn_amax_1d(ji) - at_i_1d(ji) ) ) THEN ! max is for roundoff error 323 zda_res(ji) = za_newice(ji) - MAX( 0._wp, rn_amax_1d(ji) - at_i_1d(ji) ) 326 324 zdv_res(ji) = zda_res (ji) * zh_newice(ji) 327 za_newice(ji) = za_newice(ji) - zda_res (ji)328 zv_newice(ji) = zv_newice(ji) - zdv_res (ji)325 za_newice(ji) = MAX( 0._wp, za_newice(ji) - zda_res (ji) ) 326 zv_newice(ji) = MAX( 0._wp, zv_newice(ji) - zdv_res (ji) ) 329 327 ELSE 330 328 zda_res(ji) = 0._wp -
NEMO/trunk/src/ICE/icevar.F90
r10945 r10994 44 44 !! ice_var_salprof1d : salinity profile in the ice 1D 45 45 !! ice_var_zapsmall : remove very small area and volume 46 !! ice_var_zapneg : remove negative ice fields (to debug the advection scheme UM3-5) 46 !! ice_var_zapneg : remove negative ice fields 47 !! ice_var_roundoff : remove negative values arising from roundoff erros 47 48 !! ice_var_itd : convert 1-cat to jpl-cat 48 49 !! ice_var_itd2 : convert N-cat to jpl-cat … … 71 72 PUBLIC ice_var_zapsmall 72 73 PUBLIC ice_var_zapneg 74 PUBLIC ice_var_roundoff 73 75 PUBLIC ice_var_itd 74 76 PUBLIC ice_var_itd2 … … 622 624 END SUBROUTINE ice_var_zapneg 623 625 626 627 SUBROUTINE ice_var_roundoff( pa_i, pv_i, pv_s, psv_i, poa_i, pa_ip, pv_ip, pe_s, pe_i ) 628 !!------------------------------------------------------------------- 629 !! *** ROUTINE ice_var_roundoff *** 630 !! 631 !! ** Purpose : Remove negative sea ice values arising from roundoff errors 632 !!------------------------------------------------------------------- 633 REAL(wp), DIMENSION(:,:) , INTENT(inout) :: pa_i ! ice concentration 634 REAL(wp), DIMENSION(:,:) , INTENT(inout) :: pv_i ! ice volume 635 REAL(wp), DIMENSION(:,:) , INTENT(inout) :: pv_s ! snw volume 636 REAL(wp), DIMENSION(:,:) , INTENT(inout) :: psv_i ! salt content 637 REAL(wp), DIMENSION(:,:) , INTENT(inout) :: poa_i ! age content 638 REAL(wp), DIMENSION(:,:) , INTENT(inout) :: pa_ip ! melt pond fraction 639 REAL(wp), DIMENSION(:,:) , INTENT(inout) :: pv_ip ! melt pond volume 640 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: pe_s ! snw heat content 641 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: pe_i ! ice heat content 642 !!------------------------------------------------------------------- 643 ! 644 WHERE( pa_i (1:npti,:) < 0._wp .AND. pa_i (1:npti,:) > -epsi10 ) pa_i (1:npti,:) = 0._wp ! a_i must be >= 0 645 WHERE( pv_i (1:npti,:) < 0._wp .AND. pv_i (1:npti,:) > -epsi10 ) pv_i (1:npti,:) = 0._wp ! v_i must be >= 0 646 WHERE( pv_s (1:npti,:) < 0._wp .AND. pv_s (1:npti,:) > -epsi10 ) pv_s (1:npti,:) = 0._wp ! v_s must be >= 0 647 WHERE( psv_i(1:npti,:) < 0._wp .AND. psv_i(1:npti,:) > -epsi10 ) psv_i(1:npti,:) = 0._wp ! sv_i must be >= 0 648 WHERE( poa_i(1:npti,:) < 0._wp .AND. poa_i(1:npti,:) > -epsi10 ) poa_i(1:npti,:) = 0._wp ! oa_i must be >= 0 649 WHERE( pe_i (1:npti,:,:) < 0._wp .AND. pe_i (1:npti,:,:) > -epsi06 ) pe_i (1:npti,:,:) = 0._wp ! e_i must be >= 0 650 WHERE( pe_s (1:npti,:,:) < 0._wp .AND. pe_s (1:npti,:,:) > -epsi06 ) pe_s (1:npti,:,:) = 0._wp ! e_s must be >= 0 651 IF ( ln_pnd_H12 ) THEN 652 WHERE( pa_ip(1:npti,:) < 0._wp .AND. pa_ip(1:npti,:) > -epsi10 ) pa_ip(1:npti,:) = 0._wp ! a_ip must be >= 0 653 WHERE( pv_ip(1:npti,:) < 0._wp .AND. pv_ip(1:npti,:) > -epsi10 ) pv_ip(1:npti,:) = 0._wp ! v_ip must be >= 0 654 ENDIF 655 ! 656 END SUBROUTINE ice_var_roundoff 657 624 658 625 659 SUBROUTINE ice_var_itd( zhti, zhts, zati, zh_i, zh_s, za_i ) … … 656 690 INTEGER :: idim, i_fill, jl0 657 691 REAL(wp) :: zarg, zV, zconv, zdh, zdv 658 REAL(wp), DIMENSION(:), INTENT(in) :: zhti, zhts, zati ! input ice/snow variables659 REAL(wp), DIMENSION(:,:), INTENT(inout) :: zh_i, zh_s, za_i ! output ice/snow variables692 REAL(wp), DIMENSION(:), INTENT(in) :: zhti, zhts, zati ! input ice/snow variables 693 REAL(wp), DIMENSION(:,:), INTENT(inout) :: zh_i, zh_s, za_i ! output ice/snow variables 660 694 INTEGER , DIMENSION(4) :: itest 661 695 !!------------------------------------------------------------------- … … 786 820 !! 787 821 !! 2) Expand the filling to the cat jlmin-1 and jlmax+1 788 !! by removing 25% ice area from jlmin and jlmax (resp.)822 !! by removing 25% ice area from jlmin and jlmax (resp.) 789 823 !! 790 824 !! 3) Expand the filling to the empty cat between jlmin and jlmax
Note: See TracChangeset
for help on using the changeset viewer.
