- Timestamp:
- 2017-09-26T15:24:17+02:00 (7 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icevar.F90
r8559 r8563 15 15 !! - oa_i (jpi,jpj,jpl) 16 16 !! VEQV : equivalent variables sometimes used in the model 17 !! - h t_i(jpi,jpj,jpl)18 !! - h t_s(jpi,jpj,jpl)19 !! - t_i 17 !! - h_i(jpi,jpj,jpl) 18 !! - h_s(jpi,jpj,jpl) 19 !! - t_i(jpi,jpj,nlay_i,jpl) 20 20 !! ... 21 21 !! VAGG : aggregate variables, averaged/summed over all … … 174 174 END WHERE 175 175 ! 176 h t_i(:,:,:) = v_i (:,:,:) * z1_a_i(:,:,:) !--- ice thickness176 h_i(:,:,:) = v_i (:,:,:) * z1_a_i(:,:,:) !--- ice thickness 177 177 178 178 zhmax = hi_max(jpl) 179 179 z1_zhmax = 1._wp / hi_max(jpl) 180 WHERE( h t_i(:,:,jpl) > zhmax ) !--- bound ht_i by hi_max (i.e. 99 m) with associated update of ice area181 h t_i (:,:,jpl) = zhmax180 WHERE( h_i(:,:,jpl) > zhmax ) !--- bound h_i by hi_max (i.e. 99 m) with associated update of ice area 181 h_i (:,:,jpl) = zhmax 182 182 a_i (:,:,jpl) = v_i(:,:,jpl) * z1_zhmax 183 z1_a_i(:,:,jpl) = zhmax * z1_v_i(:,:,jpl) ! NB: v_i always /=0 as h t_i > hi_max183 z1_a_i(:,:,jpl) = zhmax * z1_v_i(:,:,jpl) ! NB: v_i always /=0 as h_i > hi_max 184 184 END WHERE 185 185 186 h t_s(:,:,:) = v_s (:,:,:) * z1_a_i(:,:,:) !--- snow thickness186 h_s(:,:,:) = v_s (:,:,:) * z1_a_i(:,:,:) !--- snow thickness 187 187 188 188 o_i(:,:,:) = oa_i(:,:,:) * z1_a_i(:,:,:) !--- ice age … … 253 253 !!------------------------------------------------------------------- 254 254 ! 255 v_i (:,:,:) = h t_i(:,:,:) * a_i(:,:,:)256 v_s (:,:,:) = h t_s(:,:,:) * a_i(:,:,:)255 v_i (:,:,:) = h_i(:,:,:) * a_i(:,:,:) 256 v_s (:,:,:) = h_s(:,:,:) * a_i(:,:,:) 257 257 smv_i(:,:,:) = sm_i(:,:,:) * v_i(:,:,:) 258 258 ! … … 305 305 END DO 306 306 ! ! Slope of the linear profile 307 WHERE( h t_i(:,:,:) > epsi20 ) ; z_slope_s(:,:,:) = 2._wp * sm_i(:,:,:) / ht_i(:,:,:)307 WHERE( h_i(:,:,:) > epsi20 ) ; z_slope_s(:,:,:) = 2._wp * sm_i(:,:,:) / h_i(:,:,:) 308 308 ELSEWHERE ; z_slope_s(:,:,:) = 0._wp 309 309 END WHERE … … 326 326 DO ji = 1, jpi 327 327 ! ! linear profile with 0 surface value 328 zs0 = z_slope_s(ji,jj,jl) * ( REAL(jk,wp) - 0.5_wp ) * h t_i(ji,jj,jl) * r1_nlay_i328 zs0 = z_slope_s(ji,jj,jl) * ( REAL(jk,wp) - 0.5_wp ) * h_i(ji,jj,jl) * r1_nlay_i 329 329 zs = zalpha(ji,jj,jl) * zs0 + ( 1._wp - zalpha(ji,jj,jl) ) * sm_i(ji,jj,jl) ! weighting the profile 330 330 s_i(ji,jj,jk,jl) = MIN( rn_simax, MAX( zs, rn_simin ) ) … … 389 389 ! 390 390 ! ! Slope of the linear profile 391 WHERE( h t_i_1d(1:nidx) > epsi20 ) ; z_slope_s(1:nidx) = 2._wp * sm_i_1d(1:nidx) / ht_i_1d(1:nidx)391 WHERE( h_i_1d(1:nidx) > epsi20 ) ; z_slope_s(1:nidx) = 2._wp * sm_i_1d(1:nidx) / h_i_1d(1:nidx) 392 392 ELSEWHERE ; z_slope_s(1:nidx) = 0._wp 393 393 END WHERE … … 404 404 DO ji = 1, nidx 405 405 ! ! linear profile with 0 surface value 406 zs0 = z_slope_s(ji) * ( REAL(jk,wp) - 0.5_wp ) * h t_i_1d(ji) * r1_nlay_i406 zs0 = z_slope_s(ji) * ( REAL(jk,wp) - 0.5_wp ) * h_i_1d(ji) * r1_nlay_i 407 407 zs = zalpha(ji) * zs0 + ( 1._wp - zalpha(ji) ) * sm_i_1d(ji) 408 408 s_i_1d(ji,jk) = MIN( rn_simax , MAX( zs , rn_simin ) ) … … 447 447 ! Zap ice energy and use ocean heat to melt ice 448 448 !----------------------------------------------------------------- 449 WHERE( a_i(:,:,jl) > epsi20 ) ; h t_i(:,:,jl) = v_i(:,:,jl) / a_i(:,:,jl)450 ELSEWHERE ; h t_i(:,:,jl) = 0._wp449 WHERE( a_i(:,:,jl) > epsi20 ) ; h_i(:,:,jl) = v_i(:,:,jl) / a_i(:,:,jl) 450 ELSEWHERE ; h_i(:,:,jl) = 0._wp 451 451 END WHERE 452 452 ! 453 WHERE( a_i(:,:,jl) < epsi10 .OR. v_i(:,:,jl) < epsi10 .OR. h t_i(:,:,jl) < epsi10 ) ; zswitch(:,:) = 0._wp453 WHERE( a_i(:,:,jl) < epsi10 .OR. v_i(:,:,jl) < epsi10 .OR. h_i(:,:,jl) < epsi10 ) ; zswitch(:,:) = 0._wp 454 454 ELSEWHERE ; zswitch(:,:) = 1._wp 455 455 END WHERE … … 490 490 smv_i(ji,jj,jl) = smv_i(ji,jj,jl) * zswitch(ji,jj) 491 491 492 h t_i (ji,jj,jl) = ht_i (ji,jj,jl) * zswitch(ji,jj)493 h t_s (ji,jj,jl) = ht_s (ji,jj,jl) * zswitch(ji,jj)492 h_i (ji,jj,jl) = h_i (ji,jj,jl) * zswitch(ji,jj) 493 h_s (ji,jj,jl) = h_s (ji,jj,jl) * zswitch(ji,jj) 494 494 495 495 END DO … … 519 519 520 520 521 SUBROUTINE ice_var_itd( zhti, zhts, za i, zht_i, zht_s, za_i )521 SUBROUTINE ice_var_itd( zhti, zhts, zati, zh_i, zh_s, za_i ) 522 522 !!------------------------------------------------------------------- 523 523 !! *** ROUTINE ice_var_itd *** … … 543 543 !! ** Arguments : zhti: 1-cat ice thickness 544 544 !! zhts: 1-cat snow depth 545 !! za i: 1-cat ice concentration545 !! zati: 1-cat ice concentration 546 546 !! 547 547 !! ** Output : jpl-cat … … 552 552 INTEGER :: ijpij, i_fill, jl0 553 553 REAL(wp) :: zarg, zV, zconv, zdh, zdv 554 REAL(wp), DIMENSION(:), INTENT(in) :: zhti, zhts, za i ! input ice/snow variables555 REAL(wp), DIMENSION(:,:), INTENT(inout) :: zh t_i, zht_s, za_i ! output ice/snow variables554 REAL(wp), DIMENSION(:), INTENT(in) :: zhti, zhts, zati ! input ice/snow variables 555 REAL(wp), DIMENSION(:,:), INTENT(inout) :: zh_i, zh_s, za_i ! output ice/snow variables 556 556 INTEGER , DIMENSION(4) :: itest 557 557 !!------------------------------------------------------------------- … … 564 564 ! volume and area conservation, positivity and ice categories bounds 565 565 ijpij = SIZE( zhti , 1 ) 566 zh t_i(1:ijpij,1:jpl) = 0._wp567 zh t_s(1:ijpij,1:jpl) = 0._wp566 zh_i(1:ijpij,1:jpl) = 0._wp 567 zh_s(1:ijpij,1:jpl) = 0._wp 568 568 za_i (1:ijpij,1:jpl) = 0._wp 569 569 … … 587 587 i_fill = i_fill - 1 588 588 ! 589 zh t_i(ji,1:jpl) = 0._wp589 zh_i(ji,1:jpl) = 0._wp 590 590 za_i (ji,1:jpl) = 0._wp 591 591 itest(:) = 0 592 592 593 593 IF ( i_fill == 1 ) THEN !-- case very thin ice: fill only category 1 594 zh t_i(ji,1) = zhti(ji)595 za_i (ji,1) = za i (ji)594 zh_i(ji,1) = zhti(ji) 595 za_i (ji,1) = zati (ji) 596 596 ELSE !-- case ice is thicker: fill categories >1 597 597 ! thickness 598 598 DO jl = 1, i_fill - 1 599 zh t_i(ji,jl) = hi_mean(jl)599 zh_i(ji,jl) = hi_mean(jl) 600 600 END DO 601 601 602 602 ! concentration 603 za_i(ji,jl0) = za i(ji) / SQRT(REAL(jpl))603 za_i(ji,jl0) = zati(ji) / SQRT(REAL(jpl)) 604 604 DO jl = 1, i_fill - 1 605 605 IF ( jl /= jl0 ) THEN 606 zarg = ( zh t_i(ji,jl) - zhti(ji) ) / ( zhti(ji) * 0.5_wp )606 zarg = ( zh_i(ji,jl) - zhti(ji) ) / ( zhti(ji) * 0.5_wp ) 607 607 za_i(ji,jl) = za_i (ji,jl0) * EXP(-zarg**2) 608 608 ENDIF … … 610 610 611 611 ! last category 612 za_i(ji,i_fill) = za i(ji) - SUM( za_i(ji,1:i_fill-1) )613 zV = SUM( za_i(ji,1:i_fill-1) * zh t_i(ji,1:i_fill-1) )614 zh t_i(ji,i_fill) = ( zhti(ji) * zai(ji) - zV ) / MAX( za_i(ji,i_fill), epsi10 )612 za_i(ji,i_fill) = zati(ji) - SUM( za_i(ji,1:i_fill-1) ) 613 zV = SUM( za_i(ji,1:i_fill-1) * zh_i(ji,1:i_fill-1) ) 614 zh_i(ji,i_fill) = ( zhti(ji) * zati(ji) - zV ) / MAX( za_i(ji,i_fill), epsi10 ) 615 615 616 616 ! clem: correction if concentration of upper cat is greater than lower cat … … 619 619 DO jl = jpl, jl0+1, -1 620 620 IF ( za_i(ji,jl) > za_i(ji,jl-1) ) THEN 621 zdv = zh t_i(ji,jl) * za_i(ji,jl)622 zh t_i(ji,jl ) = 0._wp621 zdv = zh_i(ji,jl) * za_i(ji,jl) 622 zh_i(ji,jl ) = 0._wp 623 623 za_i (ji,jl ) = 0._wp 624 624 za_i (ji,1:jl-1) = za_i(ji,1:jl-1) + zdv / MAX( REAL(jl-1) * zhti(ji), epsi10 ) … … 630 630 631 631 ! Compatibility tests 632 zconv = ABS( za i(ji) - SUM( za_i(ji,1:jpl) ) )632 zconv = ABS( zati(ji) - SUM( za_i(ji,1:jpl) ) ) 633 633 IF ( zconv < epsi06 ) itest(1) = 1 ! Test 1: area conservation 634 634 635 zconv = ABS( zhti(ji)*za i(ji) - SUM( za_i(ji,1:jpl)*zht_i(ji,1:jpl) ) )635 zconv = ABS( zhti(ji)*zati(ji) - SUM( za_i(ji,1:jpl)*zh_i(ji,1:jpl) ) ) 636 636 IF ( zconv < epsi06 ) itest(2) = 1 ! Test 2: volume conservation 637 637 638 IF ( zh t_i(ji,i_fill) >= hi_max(i_fill-1) ) itest(3) = 1 ! Test 3: thickness of the last category is in-bounds ?638 IF ( zh_i(ji,i_fill) >= hi_max(i_fill-1) ) itest(3) = 1 ! Test 3: thickness of the last category is in-bounds ? 639 639 640 640 itest(4) = 1 … … 652 652 DO ji = 1, ijpij 653 653 IF( za_i(ji,jl) > 0._wp ) THEN 654 zh t_s(ji,jl) = zht_i(ji,jl) * ( zhts(ji) / zhti(ji) )654 zh_s(ji,jl) = zh_i(ji,jl) * ( zhts(ji) / zhti(ji) ) 655 655 ! In case snow load is in excess that would lead to transformation from snow to ice 656 656 ! Then, transfer the snow excess into the ice (different from icethd_dh) 657 zdh = MAX( 0._wp, ( rhosn * zh t_s(ji,jl) + ( rhoic - rau0 ) * zht_i(ji,jl) ) * r1_rau0 )658 ! recompute h t_i, ht_s avoiding out of bounds values659 zh t_i(ji,jl) = MIN( hi_max(jl), zht_i(ji,jl) + zdh )660 zh t_s(ji,jl) = MAX( 0._wp, zht_s(ji,jl) - zdh * rhoic * r1_rhosn )657 zdh = MAX( 0._wp, ( rhosn * zh_s(ji,jl) + ( rhoic - rau0 ) * zh_i(ji,jl) ) * r1_rau0 ) 658 ! recompute h_i, h_s avoiding out of bounds values 659 zh_i(ji,jl) = MIN( hi_max(jl), zh_i(ji,jl) + zdh ) 660 zh_s(ji,jl) = MAX( 0._wp, zh_s(ji,jl) - zdh * rhoic * r1_rhosn ) 661 661 ENDIF 662 662 END DO
Note: See TracChangeset
for help on using the changeset viewer.