New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 13810 for NEMO/branches – NEMO

Changeset 13810 for NEMO/branches


Ignore:
Timestamp:
2020-11-18T11:17:25+01:00 (3 years ago)
Author:
vancop
Message:

Some quick fixes

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2020/SI3-05_MP/src/ICE/icethd_pnd.F90

    r13809 r13810  
    277277             
    278278            ! Calculate the permeability of the ice (Assur 1958, see Flocco 2010) 
    279             ! MV the expression for zsbr has wrong sign!!!! 
    280             ! MV note - the sign-corrected expression is inconsistent  
    281             ! with the rest of the SI3 code which assumes linear liquidus 
    282             ! best expression (most consistent for SI3 should be) 
    283279            ! zsbr = - ( t_i_1d(ji,jk) - rt0 ) / rTmlt 
    284280            DO jk = 1, nlay_i 
     
    287283                  &   - 0.919_wp   * ( t_i_1d(ji,jk) - rt0 )**2 & 
    288284                  &   - 0.0178_wp  * ( t_i_1d(ji,jk) - rt0 )**3 
    289                 
    290                ztmp(jk) = sz_i_1d(ji,jk) / zsbr ! MV -> this is brine fraction and as it reads is always negative 
     285!               zsbr = - ( t_i_1d(ji,jk) - rt0 ) / rTmlt 
     286               ztmp(jk) = sz_i_1d(ji,jk) / zsbr  
    291287            END DO 
    292             zperm = MAX( 0._wp, 3.e-08_wp * MINVAL(ztmp)**3 ) ! MV at present since ztmp is negative, this is always zero 
     288             
     289            zperm = MAX( 0._wp, 3.e-08_wp * MINVAL(ztmp)**3 )  
    293290             
    294291            ! Do the drainage using Darcy's law 
     
    504501 
    505502      !js 03/05/19: we truncate negative values after calculating zvolp, in a 
    506       ! similar manner to the subroutine lim_mp_cesm. Variation dh_i_pnd and 
     503      ! similar manner to the subroutine ice_thd_pnd_cesm. Variation dh_i_pnd and 
    507504      ! dh_s_pnd are negative, indicating a loss of ice or snow. But we can expect them 
    508505      ! to be negative for some reasons. We keep this behaviour as it is, for 
     
    582579          ! OLI 07/2017: Done like for empirical melt pond scheme (CESM). 
    583580          ! Therefore, I chose to put this part of the code before the main 
    584           ! routines lim_mp_area/depth (contrary to the original code), 
     581          ! routines ice_thd_pnd_area/depth (contrary to the original code), 
    585582          ! seeing the freeze-up as a global sink of 
    586583          ! freshwater for melt ponds in the whole grid cell. If this was done 
     
    626623          !-------------------------------------------------------------- 
    627624          zdvn = 0._wp 
    628           CALL lim_mp_area(aice(ji,jj),vice(ji,jj), & 
     625          CALL ice_thd_pnd_area(aice(ji,jj),vice(ji,jj), & 
    629626                    aicen(ji,jj,:), vicen(ji,jj,:), vsnon(ji,jj,:), & 
    630627                    ticen(ji,jj,:,:), salin(ji,jj,:,:), & 
     
    663660   END SUBROUTINE pnd_TOPO 
    664661    
    665        SUBROUTINE lim_mp_area(aice,  vice, & 
     662       SUBROUTINE ice_thd_pnd_area(aice,  vice, & 
    666663                           aicen, vicen, vsnon, ticen, & 
    667664                           salin, zvolpn, zvolp,         & 
     
    669666 
    670667       !!------------------------------------------------------------------- 
    671        !!                ***  ROUTINE lim_mp_area *** 
     668       !!                ***  ROUTINE ice_thd_pnd_area *** 
    672669       !! 
    673670       !! ** Purpose : Given the total volume of meltwater, update 
     
    874871       ! height and area corresponding to the remaining volume 
    875872 
    876       CALL lim_mp_depth(reduced_aicen, asnon, hsnon, rhos, alfan, zvolp, cum_max_vol, hpond, m_index) 
     873      CALL ice_thd_pnd_depth(reduced_aicen, asnon, hsnon, rhos, alfan, zvolp, cum_max_vol, hpond, m_index) 
    877874 
    878875       DO n=1, m_index 
     
    910907             !IF (hicen(n) > z0) THEN           !js: from CICE 5.1.2 
    911908                perm = 0._wp ! MV ugly dummy patch 
    912                 CALL lim_mp_perm(ticen(:,n), salin(:,n), perm) 
     909                CALL ice_thd_pnd_perm(ticen(:,n), salin(:,n), perm) 
    913910                IF (perm > z0) permflag = 1 
    914911 
     
    927924          IF (permflag > 0) THEN 
    928925             ! recompute pond depth 
    929             CALL lim_mp_depth(reduced_aicen, asnon, hsnon, rhos, alfan, zvolp, cum_max_vol, hpond, m_index) 
     926            CALL ice_thd_pnd_depth(reduced_aicen, asnon, hsnon, rhos, alfan, zvolp, cum_max_vol, hpond, m_index) 
    930927             DO n=1, m_index 
    931928                zhpondn(n) = hpond - alfan(n) + alfan(1) 
     
    989986       END DO 
    990987 
    991     END SUBROUTINE lim_mp_area 
    992  
    993  
    994     SUBROUTINE lim_mp_depth(aicen, asnon, hsnon, rhos, alfan, zvolp, cum_max_vol, hpond, m_index) 
     988    END SUBROUTINE ice_thd_pnd_area 
     989 
     990 
     991    SUBROUTINE ice_thd_pnd_depth(aicen, asnon, hsnon, rhos, alfan, zvolp, cum_max_vol, hpond, m_index) 
    995992       !!------------------------------------------------------------------- 
    996        !!                ***  ROUTINE lim_mp_depth  *** 
     993       !!                ***  ROUTINE ice_thd_pnd_depth  *** 
    997994       !! 
    998995       !! ** Purpose :   Compute melt pond depth 
     
    11541151       END IF 
    11551152 
    1156     END SUBROUTINE lim_mp_depth 
    1157  
    1158  
    1159     SUBROUTINE lim_mp_perm(ticen, salin, perm) 
     1153    END SUBROUTINE ice_thd_pnd_depth 
     1154 
     1155 
     1156    SUBROUTINE ice_thd_pnd_perm(ticen, salin, perm) 
    11601157       !!------------------------------------------------------------------- 
    1161        !!                ***  ROUTINE lim_mp_perm *** 
     1158       !!                ***  ROUTINE ice_thd_pnd_perm *** 
    11621159       !! 
    11631160       !! ** Purpose :   Determine the liquid fraction of brine in the ice 
     
    11931190       !----------------------------------------------------------------- 
    11941191 
    1195        IF (maxval(Tin) <= -2.0_wp) THEN 
    1196  
    1197           ! Assur 1958 
    1198           DO k = 1,nlay_i 
    1199              Sbr = - 1.2_wp                 & 
    1200                    - 21.8_wp    * Tin(k)    & 
    1201                    - 0.919_wp   * Tin(k)**2 & 
    1202                    - 0.01878_wp * Tin(k)**3 
    1203              phi(k) = salin(k) / Sbr ! liquid fraction 
    1204  
    1205              !js: Sbr must not be zero 
    1206          END DO ! k 
    1207  
    1208        ELSE 
    1209  
    1210           ! Notz 2005 thesis eq. 3.2 
    1211           !js: "interstitial brine Sbr (in ppt) as a function of temperature T (in degC)". 
    1212           DO k = 1,nlay_i 
    1213              Sbr = - 17.6_wp   * Tin(k)    & 
    1214                    - 0.389_wp  * Tin(k)**2 & 
    1215                    - 0.00362_wp* Tin(k)**3 
    1216              phi(k) = salin(k) / Sbr ! liquid fraction 
    1217  
    1218              !js: Sbr must not be zero 
    1219           END DO 
    1220  
    1221        END IF 
     1192       DO k = 1, nlay_i 
     1193        
     1194          Sbr    = - Tin(k) / rTmlt ! Consistent expression with SI3 (linear liquidus) 
     1195          ! Best expression to date is that one 
     1196          ! Sbr  = - 18.7 * Tin(k) − 0.519 * Tin(k)**2 − 0.00535 * Tin(k) ***3 
     1197          phi(k) = salin(k) / Sbr 
     1198           
     1199       END DO 
    12221200 
    12231201       !----------------------------------------------------------------- 
     
    12271205       perm = 3.0e-08_wp * (minval(phi))**3 ! Golden et al. (2007) 
    12281206 
    1229    END SUBROUTINE lim_mp_perm 
     1207   END SUBROUTINE ice_thd_pnd_perm 
    12301208    
    12311209    
Note: See TracChangeset for help on using the changeset viewer.