Changeset 14302
- Timestamp:
- 2021-01-14T17:21:39+01:00 (3 years ago)
- Location:
- NEMO/branches/UKMO/NEMO_4.0.4_GO8_package
- Files:
-
- 17 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/cfgs/METO_GO/EXPREF/file_def_nemo-ice-BASIC-10D.xml
r12848 r14302 5 5 ============================================================================================================ 6 6 = output files definition = 7 = Define your own files for lim3=7 = Define your own files for sea ice = 8 8 = put the variables you want... = 9 9 ============================================================================================================ … … 28 28 <field field_ref="iceapnd" name="siapnd" /> 29 29 <field field_ref="icevpnd" name="sivpnd" /> 30 <field field_ref="icevlid" name="sivlid" /> 30 31 <!-- sst_m is always the potential temperature even when using teos10 --> 31 32 <field field_ref="sst_m_pot" name="sst_m_pot" /> … … 36 37 <field field_ref="icettop" name="sittop" /> 37 38 <field field_ref="icetbot" name="sitbot" /> 39 <field field_ref="icetsni" name="sitsni" /> 40 41 <!-- ponds --> 42 <field field_ref="dvpn_mlt" name="dvpn_mlt" /> 43 <field field_ref="dvpn_lid" name="dvpn_lid" /> 44 <field field_ref="dvpn_rnf" name="dvpn_rnf" /> 45 <field field_ref="dvpn_drn" name="dvpn_drn" /> 38 46 39 47 <!-- momentum --> … … 85 93 <field field_ref="icesalt_cat" name="sisalcat"/> 86 94 <field field_ref="icetemp_cat" name="sitemcat"/> 95 <field field_ref="iceapnd_cat" name="siapncat"/> 96 <field field_ref="icevpnd_cat" name="sivpncat"/> 87 97 <field field_ref="snwtemp_cat" name="sntemcat"/> 88 98 -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/cfgs/METO_GO/EXPREF/file_def_nemo-ice-BASIC-1M.xml
r12848 r14302 5 5 ============================================================================================================ 6 6 = output files definition = 7 = Define your own files for lim3=7 = Define your own files for sea ice = 8 8 = put the variables you want... = 9 9 ============================================================================================================ … … 28 28 <field field_ref="iceapnd" name="siapnd" /> 29 29 <field field_ref="icevpnd" name="sivpnd" /> 30 <field field_ref="icevlid" name="sivlid" /> 30 31 <!-- sst_m is always the potential temperature even when using teos10 --> 31 32 <field field_ref="sst_m_pot" name="sst_m_pot" /> … … 36 37 <field field_ref="icettop" name="sittop" /> 37 38 <field field_ref="icetbot" name="sitbot" /> 39 <field field_ref="icetsni" name="sitsni" /> 40 41 <!-- ponds --> 42 <field field_ref="dvpn_mlt" name="dvpn_mlt" /> 43 <field field_ref="dvpn_lid" name="dvpn_lid" /> 44 <field field_ref="dvpn_rnf" name="dvpn_rnf" /> 45 <field field_ref="dvpn_drn" name="dvpn_drn" /> 38 46 39 47 <!-- momentum --> … … 85 93 <field field_ref="icesalt_cat" name="sisalcat"/> 86 94 <field field_ref="icetemp_cat" name="sitemcat"/> 95 <field field_ref="iceapnd_cat" name="siapncat"/> 96 <field field_ref="icevpnd_cat" name="sivpncat"/> 87 97 <field field_ref="snwtemp_cat" name="sntemcat"/> 88 98 -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/cfgs/METO_GO/EXPREF/file_def_nemo-ice-CMIP6-10D.xml
r12848 r14302 5 5 ============================================================================================================ 6 6 = output files definition = 7 = Define your own files for lim3=7 = Define your own files for sea ice = 8 8 = put the variables you want... = 9 9 ============================================================================================================ … … 28 28 <field field_ref="iceapnd" name="siapnd" /> 29 29 <field field_ref="icevpnd" name="sivpnd" /> 30 <field field_ref="icevlid" name="sivlid" /> 30 31 <!-- sst_m is always the potential temperature even when using teos10 --> 31 32 <field field_ref="sst_m_pot" name="sst_m_pot" /> … … 36 37 <field field_ref="icettop" name="sittop" /> 37 38 <field field_ref="icetbot" name="sitbot" /> 39 <field field_ref="icetsni" name="sitsni" /> 40 41 <!-- ponds --> 42 <field field_ref="dvpn_mlt" name="dvpn_mlt" /> 43 <field field_ref="dvpn_lid" name="dvpn_lid" /> 44 <field field_ref="dvpn_rnf" name="dvpn_rnf" /> 45 <field field_ref="dvpn_drn" name="dvpn_drn" /> 38 46 39 47 <!-- momentum --> … … 85 93 <field field_ref="icesalt_cat" name="sisalcat"/> 86 94 <field field_ref="icetemp_cat" name="sitemcat"/> 95 <field field_ref="iceapnd_cat" name="siapncat"/> 96 <field field_ref="icevpnd_cat" name="sivpncat"/> 87 97 <field field_ref="snwtemp_cat" name="sntemcat"/> 88 98 89 99 <!-- mass balance --> 90 100 <field field_ref="dmithd" name="sidmassth" /> -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/cfgs/METO_GO/EXPREF/file_def_nemo-ice-CMIP6-1M.xml
r12848 r14302 5 5 ============================================================================================================ 6 6 = output files definition = 7 = Define your own files for lim3=7 = Define your own files for sea ice = 8 8 = put the variables you want... = 9 9 ============================================================================================================ … … 28 28 <field field_ref="iceapnd" name="siapnd" /> 29 29 <field field_ref="icevpnd" name="sivpnd" /> 30 <field field_ref="icevlid" name="sivlid" /> 30 31 <!-- sst_m is always the potential temperature even when using teos10 --> 31 32 <field field_ref="sst_m_pot" name="sst_m_pot" /> … … 36 37 <field field_ref="icettop" name="sittop" /> 37 38 <field field_ref="icetbot" name="sitbot" /> 39 <field field_ref="icetsni" name="sitsni" /> 40 41 <!-- ponds --> 42 <field field_ref="dvpn_mlt" name="dvpn_mlt" /> 43 <field field_ref="dvpn_lid" name="dvpn_lid" /> 44 <field field_ref="dvpn_rnf" name="dvpn_rnf" /> 45 <field field_ref="dvpn_drn" name="dvpn_drn" /> 38 46 39 47 <!-- momentum --> … … 85 93 <field field_ref="icesalt_cat" name="sisalcat"/> 86 94 <field field_ref="icetemp_cat" name="sitemcat"/> 95 <field field_ref="iceapnd_cat" name="siapncat"/> 96 <field field_ref="icevpnd_cat" name="sivpncat"/> 87 97 <field field_ref="snwtemp_cat" name="sntemcat"/> 88 98 89 99 <!-- mass balance --> 90 100 <field field_ref="dmithd" name="sidmassth" /> -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/cfgs/ORCA2_SAS_ICE/EXPREF/file_def_nemo-ice.xml
r9896 r14302 28 28 <field field_ref="iceapnd" name="siapnd" /> 29 29 <field field_ref="icevpnd" name="sivpnd" /> 30 <field field_ref="icevlid" name="sivlid" /> 30 31 <field field_ref="sst_m" name="sst_m" /> 31 32 <field field_ref="sss_m" name="sss_m" /> … … 36 37 <field field_ref="icetbot" name="sitbot" /> 37 38 <field field_ref="icetsni" name="sitsni" /> 39 40 <!-- ponds --> 41 <field field_ref="dvpn_mlt" name="dvpn_mlt" /> 42 <field field_ref="dvpn_lid" name="dvpn_lid" /> 43 <field field_ref="dvpn_rnf" name="dvpn_rnf" /> 44 <field field_ref="dvpn_drn" name="dvpn_drn" /> 38 45 39 46 <!-- momentum --> … … 85 92 <field field_ref="icesalt_cat" name="sisalcat"/> 86 93 <field field_ref="icetemp_cat" name="sitemcat"/> 94 <field field_ref="iceapnd_cat" name="siapncat"/> 95 <field field_ref="icevpnd_cat" name="sivpncat"/> 87 96 88 97 </file> -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/cfgs/SHARED/field_def_nemo-ice.xml
r14078 r14302 51 51 <field id="icehlid" long_name="melt pond lid depth" standard_name="sea_ice_meltpondlid_depth" unit="m" /> 52 52 <field id="icevlid" long_name="melt pond lid volume" standard_name="sea_ice_meltpondlid_volume" unit="m" /> 53 54 <field id="dvpn_mlt" long_name="pond volume tendency due to surface melt" standard_name="sea_ice_pondvolume_tendency_melt" unit="cm/d" /> 55 <field id="dvpn_lid" long_name="pond volume tendency due to exchanges with lid" standard_name="sea_ice_pondvolume_tendency_lids" unit="cm/d" /> 56 <field id="dvpn_rnf" long_name="pond volume tendency due to runoff" standard_name="sea_ice_pondvolume_tendency_runoff" unit="cm/d" /> 57 <field id="dvpn_drn" long_name="pond volume tendency due to drainage" standard_name="sea_ice_pondvolume_tendency_drainage" unit="cm/d" /> 53 58 54 59 <!-- heat --> … … 284 289 <field id="snwtemp_cat" long_name="Snow temperature per category" unit="degC" detect_missing_value="true" /> 285 290 <field id="icettop_cat" long_name="Ice/snow surface temperature per category" unit="degC" detect_missing_value="true" /> 286 <field id="iceapnd_cat" long_name="Ice melt pond concentration per category" unit="" /> 291 <field id="icevpnd_cat" long_name="Ice melt pond volume per grid area per category" unit="m" /> 292 <field id="iceapnd_cat" long_name="Ice melt pond grid fraction" unit="" /> 287 293 <field id="icehpnd_cat" long_name="Ice melt pond thickness per category" unit="m" detect_missing_value="true" /> 288 294 <field id="icehlid_cat" long_name="Ice melt pond lid thickness per category" unit="m" detect_missing_value="true" /> 289 <field id="iceafpnd_cat" long_name="Ice melt pond fraction per category"unit="" />295 <field id="iceafpnd_cat" long_name="Ice melt pond ice fraction per category" unit="" /> 290 296 <field id="iceaepnd_cat" long_name="Ice melt pond effective fraction per category" unit="" /> 291 297 <field id="icemask_cat" long_name="Fraction of time step with sea ice (per category)" unit="" /> -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/cfgs/SHARED/namelist_ice_ref
r14158 r14302 195 195 !------------------------------------------------------------------------------ 196 196 ln_pnd = .true. ! activate melt ponds or not 197 ln_pnd_LEV = .true. ! level ice melt ponds (from Flocco et al 2007,2010 & Holland et al 2012) 198 rn_apnd_min = 0.15 ! minimum ice fraction that contributes to melt pond. range: 0.0 -- 0.15 ?? 199 rn_apnd_max = 0.85 ! maximum ice fraction that contributes to melt pond. range: 0.7 -- 0.85 ?? 197 ln_pnd_TOPO = .true. ! topographic melt ponds 198 ln_pnd_LEV = .false. ! level ice melt ponds 199 rn_apnd_min = 0.15 ! minimum meltwater fraction contributing to pond growth (TOPO and LEV) 200 rn_apnd_max = 0.85 ! maximum meltwater fraction contributing to pond growth (TOPO and LEV) 200 201 rn_pnd_flush= 0.01 ! pond flushing efficiency (tuning parameter) (LEV) 201 202 ln_pnd_CST = .false. ! constant melt ponds 202 203 rn_apnd = 0.2 ! prescribed pond fraction, at Tsu=0 degC 203 204 rn_hpnd = 0.05 ! prescribed pond depth, at Tsu=0 degC 204 ln_pnd_lids = . true.! frozen lids on top of the ponds (only for ln_pnd_LEV)205 ln_pnd_lids = .false. ! frozen lids on top of the ponds (only for ln_pnd_LEV) 205 206 ln_pnd_alb = .true. ! effect of melt ponds on ice albedo 206 207 / … … 227 228 rn_tms_ini_n = 270. ! initial snw temperature (K), North 228 229 rn_tms_ini_s = 270. ! " " South 229 rn_apd_ini_n = 0. 2! initial pond fraction (-), North230 rn_apd_ini_n = 0.0 ! initial pond fraction (-), North 230 231 rn_apd_ini_s = 0.2 ! " " South 231 rn_hpd_ini_n = 0.0 5! initial pond depth (m), North232 rn_hpd_ini_n = 0.00 ! initial pond depth (m), North 232 233 rn_hpd_ini_s = 0.05 ! " " South 233 234 rn_hld_ini_n = 0.0 ! initial pond lid depth (m), North -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/src/ICE/ice.F90
r14158 r14302 208 208 ! !!** ice-ponds namelist (namthd_pnd) 209 209 LOGICAL , PUBLIC :: ln_pnd !: Melt ponds (T) or not (F) 210 LOGICAL , PUBLIC :: ln_pnd_LEV !: Melt ponds scheme from Holland et al (2012), Flocco et al (2007, 2010) 211 REAL(wp), PUBLIC :: rn_apnd_min !: Minimum ice fraction that contributes to melt ponds 212 REAL(wp), PUBLIC :: rn_apnd_max !: Maximum ice fraction that contributes to melt ponds 210 LOGICAL , PUBLIC :: ln_pnd_TOPO !: Topographic Melt ponds scheme (Flocco et al 2007, 2010) 211 LOGICAL , PUBLIC :: ln_pnd_LEV !: Simple melt pond scheme 212 REAL(wp), PUBLIC :: rn_apnd_min !: Minimum fraction of melt water contributing to ponds 213 REAL(wp), PUBLIC :: rn_apnd_max !: Maximum fraction of melt water contributing to ponds 214 REAL(wp), PUBLIC :: rn_pnd_flush !: Pond flushing efficiency (tuning parameter) 213 215 REAL(wp), PUBLIC :: rn_pnd_flush !: Pond flushing efficiency (tuning parameter) 214 216 LOGICAL , PUBLIC :: ln_pnd_CST !: Melt ponds scheme with constant fraction and depth … … 309 311 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: t1_ice !: temperature of the first layer (ln_cndflx=T) [K] 310 312 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: cnd_ice !: effective conductivity of the 1st layer (ln_cndflx=T) [W.m-2.K-1] 313 314 ! meltwater arrays to save for melt ponds 315 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: dh_i_sum_2d !: surface melt (2d arrays for ponds) [m] 316 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: dh_s_mlt_2d !: snow surf melt (2d arrays for ponds) [m] 311 317 312 318 !!---------------------------------------------------------------------- … … 459 465 ii = ii + 1 460 466 ALLOCATE( qtr_ice_bot(jpi,jpj,jpl) , cnd_ice(jpi,jpj,jpl) , t1_ice(jpi,jpj,jpl) , & 467 & dh_i_sum_2d(jpi,jpj,jpl) , dh_s_mlt_2d(jpi,jpj,jpl) , & 461 468 & h_i (jpi,jpj,jpl) , a_i (jpi,jpj,jpl) , v_i (jpi,jpj,jpl) , & 462 469 & v_s (jpi,jpj,jpl) , h_s (jpi,jpj,jpl) , t_su (jpi,jpj,jpl) , & -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/src/ICE/icedyn_adv_pra.F90
r14075 r14302 178 178 z0ei(:,:,jk,jl) = pe_i(:,:,jk,jl) * e1e2t(:,:) ! Ice heat content 179 179 END DO 180 IF ( ln_pnd_LEV ) THEN180 IF ( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN 181 181 z0ap(:,:,jl) = pa_ip(:,:,jl) * e1e2t(:,:) ! Melt pond fraction 182 182 z0vp(:,:,jl) = pv_ip(:,:,jl) * e1e2t(:,:) ! Melt pond volume … … 214 214 END DO 215 215 ! 216 IF ( ln_pnd_LEV ) THEN216 IF ( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN 217 217 CALL adv_x( zdt , zudy , 1._wp , zarea , z0ap , sxap , sxxap , syap , syyap , sxyap ) !--- melt pond fraction 218 218 CALL adv_y( zdt , zvdx , 0._wp , zarea , z0ap , sxap , sxxap , syap , syyap , sxyap ) … … 249 249 & sxxe(:,:,jk,:), sye(:,:,jk,:), syye(:,:,jk,:), sxye(:,:,jk,:) ) 250 250 END DO 251 IF ( ln_pnd_LEV ) THEN251 IF ( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN 252 252 CALL adv_y( zdt , zvdx , 1._wp , zarea , z0ap , sxap , sxxap , syap , syyap , sxyap ) !--- melt pond fraction 253 253 CALL adv_x( zdt , zudy , 0._wp , zarea , z0ap , sxap , sxxap , syap , syyap , sxyap ) … … 278 278 CALL lbc_lnk_multi( 'icedyn_adv_pra', z0ei , 'T', 1._wp, sxe , 'T', -1._wp, sye , 'T', -1._wp & ! ice enthalpy 279 279 & , sxxe , 'T', 1._wp, syye , 'T', 1._wp, sxye , 'T', 1._wp ) 280 IF ( ln_pnd_LEV ) THEN280 IF ( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN 281 281 CALL lbc_lnk_multi( 'icedyn_adv_pra', z0ap , 'T', 1._wp, sxap , 'T', -1._wp, syap , 'T', -1._wp & ! melt pond fraction 282 282 & , sxxap, 'T', 1._wp, syyap, 'T', 1._wp, sxyap, 'T', 1._wp & … … 302 302 pe_i(:,:,jk,jl) = z0ei(:,:,jk,jl) * r1_e1e2t(:,:) * tmask(:,:,1) 303 303 END DO 304 IF ( ln_pnd_LEV ) THEN304 IF ( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN 305 305 pa_ip(:,:,jl) = z0ap(:,:,jl) * r1_e1e2t(:,:) * tmask(:,:,1) 306 306 pv_ip(:,:,jl) = z0vp(:,:,jl) * r1_e1e2t(:,:) * tmask(:,:,1) … … 780 780 ! ! -- check h_ip -- ! 781 781 ! if h_ip is larger than the surrounding 9 pts => reduce h_ip and increase a_ip 782 IF( ln_pnd_LEV.AND. pv_ip(ji,jj,jl) > 0._wp ) THEN782 IF( ( ln_pnd_LEV .OR. ln_pnd_TOPO ) .AND. pv_ip(ji,jj,jl) > 0._wp ) THEN 783 783 zhip = pv_ip(ji,jj,jl) / MAX( epsi20, pa_ip(ji,jj,jl) ) 784 784 IF( zhip > phip_max(ji,jj,jl) .AND. pa_ip(ji,jj,jl) < 0.15 ) THEN … … 1027 1027 END DO 1028 1028 ! 1029 IF( ln_pnd_LEV ) THEN ! melt pond fraction1029 IF( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN ! melt pond fraction 1030 1030 IF( iom_varid( numrir, 'sxap', ldstop = .FALSE. ) > 0 ) THEN 1031 1031 CALL iom_get( numrir, jpdom_autoglo, 'sxap' , sxap ) … … 1069 1069 sxc0 = 0._wp ; syc0 = 0._wp ; sxxc0 = 0._wp ; syyc0 = 0._wp ; sxyc0 = 0._wp ! snow layers heat content 1070 1070 sxe = 0._wp ; sye = 0._wp ; sxxe = 0._wp ; syye = 0._wp ; sxye = 0._wp ! ice layers heat content 1071 IF( ln_pnd_LEV ) THEN1071 IF( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN 1072 1072 sxap = 0._wp ; syap = 0._wp ; sxxap = 0._wp ; syyap = 0._wp ; sxyap = 0._wp ! melt pond fraction 1073 1073 sxvp = 0._wp ; syvp = 0._wp ; sxxvp = 0._wp ; syyvp = 0._wp ; sxyvp = 0._wp ! melt pond volume … … 1137 1137 END DO 1138 1138 ! 1139 IF( ln_pnd_LEV ) THEN ! melt pond fraction1139 IF( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN ! melt pond fraction 1140 1140 CALL iom_rstput( iter, nitrst, numriw, 'sxap' , sxap ) 1141 1141 CALL iom_rstput( iter, nitrst, numriw, 'syap' , syap ) -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/src/ICE/icedyn_adv_umx.F90
r14075 r14302 341 341 ! 342 342 !== melt ponds ==! 343 IF ( ln_pnd_LEV ) THEN343 IF ( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN 344 344 ! concentration 345 345 zamsk = 1._wp … … 361 361 ! 362 362 ! --- Lateral boundary conditions --- ! 363 IF ( ln_pnd_LEV.AND. ln_pnd_lids ) THEN363 IF ( ( ln_pnd_LEV .OR. ln_pnd_TOPO ) .AND. ln_pnd_lids ) THEN 364 364 CALL lbc_lnk_multi( 'icedyn_adv_umx', pa_i,'T',1._wp, pv_i,'T',1._wp, pv_s,'T',1._wp, psv_i,'T',1._wp, poa_i,'T',1._wp & 365 365 & , pa_ip,'T',1._wp, pv_ip,'T',1._wp, pv_il,'T',1._wp ) 366 ELSEIF( ln_pnd_LEV.AND. .NOT.ln_pnd_lids ) THEN366 ELSEIF( ( ln_pnd_LEV .OR. ln_pnd_TOPO ) .AND. .NOT.ln_pnd_lids ) THEN 367 367 CALL lbc_lnk_multi( 'icedyn_adv_umx', pa_i,'T',1._wp, pv_i,'T',1._wp, pv_s,'T',1._wp, psv_i,'T',1._wp, poa_i,'T',1._wp & 368 368 & , pa_ip,'T',1._wp, pv_ip,'T',1._wp ) … … 1611 1611 ! ! -- check h_ip -- ! 1612 1612 ! if h_ip is larger than the surrounding 9 pts => reduce h_ip and increase a_ip 1613 IF( ln_pnd_LEV.AND. pv_ip(ji,jj,jl) > 0._wp ) THEN1613 IF( ( ln_pnd_LEV .OR. ln_pnd_TOPO ) .AND. pv_ip(ji,jj,jl) > 0._wp ) THEN 1614 1614 zhip = pv_ip(ji,jj,jl) / MAX( epsi20, pa_ip(ji,jj,jl) ) 1615 1615 IF( zhip > phip_max(ji,jj,jl) .AND. pa_ip(ji,jj,jl) < 0.15 ) THEN -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/src/ICE/icedyn_rdgrft.F90
r14075 r14302 567 567 oirft2(ji) = oa_i_2d(ji,jl1) * afrft * hi_hrft 568 568 569 IF ( ln_pnd_LEV ) THEN569 IF ( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN 570 570 aprdg1 = a_ip_2d(ji,jl1) * afrdg 571 571 aprdg2(ji) = a_ip_2d(ji,jl1) * afrdg * hi_hrdg(ji,jl1) … … 604 604 sv_i_2d(ji,jl1) = sv_i_2d(ji,jl1) - sirdg1 - sirft(ji) 605 605 oa_i_2d(ji,jl1) = oa_i_2d(ji,jl1) - oirdg1 - oirft1 606 IF ( ln_pnd_LEV ) THEN606 IF ( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN 607 607 a_ip_2d(ji,jl1) = a_ip_2d(ji,jl1) - aprdg1 - aprft1 608 608 v_ip_2d(ji,jl1) = v_ip_2d(ji,jl1) - vprdg(ji) - vprft(ji) … … 701 701 v_s_2d (ji,jl2) = v_s_2d (ji,jl2) + ( vsrdg (ji) * rn_fsnwrdg * fvol(ji) + & 702 702 & vsrft (ji) * rn_fsnwrft * zswitch(ji) ) 703 IF ( ln_pnd_LEV ) THEN703 IF ( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN 704 704 v_ip_2d (ji,jl2) = v_ip_2d(ji,jl2) + ( vprdg (ji) * rn_fpndrdg * fvol (ji) & 705 705 & + vprft (ji) * rn_fpndrft * zswitch(ji) ) -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/src/ICE/iceitd.F90
r14075 r14302 305 305 IF ( a_i_1d(ji) > epsi10 .AND. h_i_1d(ji) < rn_himin ) THEN 306 306 a_i_1d(ji) = a_i_1d(ji) * h_i_1d(ji) / rn_himin 307 IF( ln_pnd_LEV ) a_ip_1d(ji) = a_ip_1d(ji) * h_i_1d(ji) / rn_himin307 IF( ln_pnd_LEV .OR. ln_pnd_TOPO ) a_ip_1d(ji) = a_ip_1d(ji) * h_i_1d(ji) / rn_himin 308 308 h_i_1d(ji) = rn_himin 309 309 ENDIF … … 476 476 zaTsfn(ji,jl2) = zaTsfn(ji,jl2) + ztrans 477 477 ! 478 IF ( ln_pnd_LEV ) THEN478 IF ( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN 479 479 ztrans = a_ip_2d(ji,jl1) * zworka(ji) ! Pond fraction 480 480 a_ip_2d(ji,jl1) = a_ip_2d(ji,jl1) - ztrans -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/src/ICE/icestp.F90
r14075 r14302 414 414 wfx_pnd(ji,jj) = 0._wp 415 415 416 dh_i_sum_2d(:,:,:) = 0._wp 417 dh_s_mlt_2d(:,:,:) = 0._wp 418 416 419 hfx_thd(ji,jj) = 0._wp ; 417 420 hfx_snw(ji,jj) = 0._wp ; hfx_opw(ji,jj) = 0._wp -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/src/ICE/icethd.F90
r14075 r14302 247 247 IF( ln_icedH ) THEN ! --- Growing/Melting --- ! 248 248 CALL ice_thd_dh ! Ice-Snow thickness 249 CALL ice_thd_pnd ! Melt ponds formation250 249 CALL ice_thd_ent( e_i_1d(1:npti,:) ) ! Ice enthalpy remapping 251 250 ENDIF … … 268 267 IF( ln_icediachk ) CALL ice_cons2D (1, 'icethd', diag_v, diag_s, diag_t, diag_fv, diag_fs, diag_ft) 269 268 ! 269 IF ( ln_pnd ) CALL ice_thd_pnd ! --- Melt ponds 270 270 271 IF( jpl > 1 ) CALL ice_itd_rem( kt ) ! --- Transport ice between thickness categories --- ! 271 272 ! … … 536 537 CALL tab_1d_2d( npti, nptidx(1:npti), cnd_ice_1d(1:npti), cnd_ice(:,:,kl) ) 537 538 CALL tab_1d_2d( npti, nptidx(1:npti), t1_ice_1d (1:npti), t1_ice (:,:,kl) ) 539 ! ponds 540 CALL tab_1d_2d( npti, nptidx(1:npti), dh_i_sum (1:npti) , dh_i_sum_2d(:,:,kl) ) 541 CALL tab_1d_2d( npti, nptidx(1:npti), dh_s_mlt (1:npti) , dh_s_mlt_2d(:,:,kl) ) 538 542 ! SIMIP diagnostics 539 543 CALL tab_1d_2d( npti, nptidx(1:npti), t_si_1d (1:npti), t_si (:,:,kl) ) -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/src/ICE/icethd_pnd.F90
r14158 r14302 20 20 USE ice1D ! sea-ice: thermodynamics variables 21 21 USE icetab ! sea-ice: 1D <==> 2D transformation 22 USE sbc_ice ! surface energy budget 22 23 ! 23 24 USE in_out_manager ! I/O manager 25 USE iom ! I/O manager library 24 26 USE lib_mpp ! MPP library 25 27 USE lib_fortran ! fortran utilities (glob_sum + no signed zero) … … 34 36 INTEGER :: nice_pnd ! choice of the type of pond scheme 35 37 ! ! associated indices: 36 INTEGER, PARAMETER :: np_pndNO = 0 ! No pond scheme 37 INTEGER, PARAMETER :: np_pndCST = 1 ! Constant ice pond scheme 38 INTEGER, PARAMETER :: np_pndLEV = 2 ! Level ice pond scheme 38 INTEGER, PARAMETER :: np_pndNO = 0 ! No pond scheme 39 INTEGER, PARAMETER :: np_pndCST = 1 ! Constant ice pond scheme 40 INTEGER, PARAMETER :: np_pndLEV = 2 ! Level ice pond scheme 41 INTEGER, PARAMETER :: np_pndTOPO = 3 ! Level ice pond scheme 42 43 REAL(wp), PARAMETER :: & ! shared parameters for topographic melt ponds 44 zhi_min = 0.1_wp , & ! minimum ice thickness with ponds (m) 45 zTd = 273.0_wp , & ! temperature difference for freeze-up (C) 46 zvp_min = 1.e-4_wp ! minimum pond volume (m) 47 48 !-------------------------------------------------------------------------- 49 ! 50 ! Pond volume per area budget diags 51 ! 52 ! The idea of diags is the volume of ponds per grid cell area is 53 ! 54 ! dV/dt = mlt + drn + lid + rnf 55 ! mlt = input from surface melting 56 ! drn = drainage through brine network 57 ! lid = lid growth & melt 58 ! rnf = runoff (water directly removed out of surface melting + overflow) 59 ! 60 ! In topo mode, the pond water lost because it is in the snow is not included in the budget 61 ! 62 ! In level mode, all terms are incorporated 63 64 REAL(wp), ALLOCATABLE, DIMENSION(:,:) :: & ! pond volume budget diagnostics 65 diag_dvpn_mlt, & !! meltwater pond volume input [m/day] 66 diag_dvpn_drn, & !! pond volume lost by drainage [m/day] 67 diag_dvpn_lid, & !! exchange with lid / refreezing [m/day] 68 diag_dvpn_rnf !! meltwater pond lost to runoff [m/day] 69 70 REAL(wp), ALLOCATABLE, DIMENSION(:) :: & ! pond volume budget diagnostics (1d) 71 diag_dvpn_mlt_1d, & !! meltwater pond volume input [m/day] 72 diag_dvpn_drn_1d, & !! pond volume lost by drainage [m/day] 73 diag_dvpn_lid_1d, & !! exchange with lid / refreezing [m/day] 74 diag_dvpn_rnf_1d !! meltwater pond lost to runoff [m/day] 75 ! 76 !-------------------------------------------------------------------------- 39 77 40 78 !! * Substitutions … … 52 90 !! 53 91 !! ** Purpose : change melt pond fraction and thickness 54 !! 92 !! 93 !! Note: Melt ponds affect only radiative transfer for now 94 !! 95 !! No heat, no salt. 96 !! The melt water they carry is collected but 97 !! not removed from fw budget or released to the ocean 98 !! 99 !! A wfx_pnd has been coded for diagnostic purposes 100 !! It is not fully consistent yet. 101 !! 102 !! The current diagnostic lacks a contribution from drainage 103 !! 55 104 !!------------------------------------------------------------------- 56 ! 105 !! 106 107 ALLOCATE( diag_dvpn_mlt(jpi,jpj), diag_dvpn_lid(jpi,jpj), diag_dvpn_drn(jpi,jpj), diag_dvpn_rnf(jpi,jpj) ) 108 ALLOCATE( diag_dvpn_mlt_1d(jpij), diag_dvpn_lid_1d(jpij), diag_dvpn_drn_1d(jpij), diag_dvpn_rnf_1d(jpij) ) 109 110 diag_dvpn_mlt(:,:) = 0._wp ; diag_dvpn_drn(:,:) = 0._wp 111 diag_dvpn_lid(:,:) = 0._wp ; diag_dvpn_rnf(:,:) = 0._wp 112 diag_dvpn_mlt_1d(:) = 0._wp ; diag_dvpn_drn_1d(:) = 0._wp 113 diag_dvpn_lid_1d(:) = 0._wp ; diag_dvpn_rnf_1d(:) = 0._wp 114 57 115 SELECT CASE ( nice_pnd ) 58 116 ! 59 CASE (np_pndCST) ; CALL pnd_CST !== Constant melt ponds ==!117 CASE (np_pndCST) ; CALL pnd_CST !== Constant melt ponds ==! 60 118 ! 61 CASE (np_pndLEV) ; CALL pnd_LEV !== Level ice melt ponds ==! 119 CASE (np_pndLEV) ; CALL pnd_LEV !== Level ice melt ponds ==! 120 ! 121 CASE (np_pndTOPO) ; CALL pnd_TOPO !== Topographic melt ponds ==! 62 122 ! 63 123 END SELECT 64 124 ! 125 126 IF( iom_use('dvpn_mlt' ) ) CALL iom_put( 'dvpn_mlt', diag_dvpn_mlt * 100._wp * 86400._wp ) ! input from melting 127 IF( iom_use('dvpn_lid' ) ) CALL iom_put( 'dvpn_lid', diag_dvpn_lid * 100._wp * 86400._wp ) ! exchanges with lid 128 IF( iom_use('dvpn_drn' ) ) CALL iom_put( 'dvpn_drn', diag_dvpn_drn * 100._wp * 86400._wp ) ! vertical drainage 129 IF( iom_use('dvpn_rnf' ) ) CALL iom_put( 'dvpn_rnf', diag_dvpn_rnf * 100._wp * 86400._wp ) ! runoff + overflow 130 131 DEALLOCATE( diag_dvpn_mlt, diag_dvpn_lid, diag_dvpn_drn, diag_dvpn_rnf ) 132 DEALLOCATE( diag_dvpn_mlt_1d, diag_dvpn_lid_1d, diag_dvpn_drn_1d, diag_dvpn_rnf_1d ) 133 65 134 END SUBROUTINE ice_thd_pnd 66 135 … … 75 144 !! to non-zero values when t_su = 0C 76 145 !! 77 !! ** Tunable parameters : pond fraction (rn_apnd), ponddepth (rn_hpnd)146 !! ** Tunable parameters : Pond fraction (rn_apnd) & depth (rn_hpnd) 78 147 !! 79 148 !! ** Note : Coupling with such melt ponds is only radiative … … 147 216 !! ** Tunable parameters : rn_apnd_max, rn_apnd_min, rn_pnd_flush 148 217 !! 149 !! ** Note : mostly stolen from CICE 150 !! 151 !! ** References : Flocco and Feltham (JGR, 2007) 152 !! Flocco et al (JGR, 2010) 153 !! Holland et al (J. Clim, 2012) 218 !! ** Note : Mostly stolen from CICE but not only 219 !! 220 !! ** References : Holland et al (J. Clim, 2012), Hunke et al (OM 2012) 221 !! 154 222 !!------------------------------------------------------------------- 155 223 REAL(wp), DIMENSION(nlay_i) :: ztmp ! temporary array … … 168 236 REAL(wp) :: zfac, zdum ! temporary arrays 169 237 REAL(wp) :: z1_rhow, z1_aspect, z1_Tp ! inverse 170 !! 171 INTEGER :: ji, jk ! loop indices 238 REAL(wp) :: zvold ! 239 !! 240 INTEGER :: ji, jj, jk, jl ! loop indices 241 172 242 !!------------------------------------------------------------------- 173 243 z1_rhow = 1._wp / rhow 174 244 z1_aspect = 1._wp / zaspect 175 245 z1_Tp = 1._wp / zTp 176 177 DO ji = 1, npti 178 ! !----------------------------------------------------! 179 IF( h_i_1d(ji) < rn_himin .OR. a_i_1d(ji) < epsi10 ) THEN ! Case ice thickness < rn_himin or tiny ice fraction ! 180 ! !----------------------------------------------------! 181 !--- Remove ponds on thin ice or tiny ice fractions 182 a_ip_1d(ji) = 0._wp 183 h_ip_1d(ji) = 0._wp 184 h_il_1d(ji) = 0._wp 185 ! !--------------------------------! 186 ELSE ! Case ice thickness >= rn_himin ! 187 ! !--------------------------------! 188 v_ip_1d(ji) = h_ip_1d(ji) * a_ip_1d(ji) ! retrieve volume from thickness 189 v_il_1d(ji) = h_il_1d(ji) * a_ip_1d(ji) 190 ! 191 !------------------! 192 ! case ice melting ! 193 !------------------! 194 ! 195 !--- available meltwater for melt ponding ---! 196 zdum = -( dh_i_sum(ji)*rhoi + dh_s_mlt(ji)*rhos ) * z1_rhow * a_i_1d(ji) 197 zfr_mlt = rn_apnd_min + ( rn_apnd_max - rn_apnd_min ) * at_i_1d(ji) ! = ( 1 - r ) = fraction of melt water that is not flushed 198 zdv_mlt = MAX( 0._wp, zfr_mlt * zdum ) ! max for roundoff errors? 199 ! 200 !--- overflow ---! 201 ! If pond area exceeds zfr_mlt * a_i_1d(ji) then reduce the pond volume 202 ! a_ip_max = zfr_mlt * a_i 203 ! => from zaspect = h_ip / (a_ip / a_i), set v_ip_max as: 204 zv_ip_max = zfr_mlt**2 * a_i_1d(ji) * zaspect 205 zdv_mlt = MAX( 0._wp, MIN( zdv_mlt, zv_ip_max - v_ip_1d(ji) ) ) 206 207 ! If pond depth exceeds half the ice thickness then reduce the pond volume 208 ! h_ip_max = 0.5 * h_i 209 ! => from zaspect = h_ip / (a_ip / a_i), set v_ip_max as: 210 zv_ip_max = z1_aspect * a_i_1d(ji) * 0.25 * h_i_1d(ji) * h_i_1d(ji) 211 zdv_mlt = MAX( 0._wp, MIN( zdv_mlt, zv_ip_max - v_ip_1d(ji) ) ) 246 247 !----------------------------------------------------------------------------------------------- 248 ! Identify grid cells with ice 249 !----------------------------------------------------------------------------------------------- 250 at_i(:,:) = SUM( a_i, dim=3 ) 251 ! 252 npti = 0 ; nptidx(:) = 0 253 DO jj = 1, jpj 254 DO ji = 1, jpi 255 IF ( at_i(ji,jj) > epsi10 ) THEN 256 npti = npti + 1 257 nptidx( npti ) = (jj - 1) * jpi + ji 258 ENDIF 259 END DO 260 END DO 261 262 !----------------------------------------------------------------------------------------------- 263 ! Prepare 1D arrays 264 !----------------------------------------------------------------------------------------------- 265 266 IF( npti > 0 ) THEN 267 268 CALL tab_2d_1d( npti, nptidx(1:npti), wfx_snw_sum_1d(1:npti), wfx_snw_sum ) 269 CALL tab_2d_1d( npti, nptidx(1:npti), wfx_sum_1d (1:npti) , wfx_sum ) 270 CALL tab_2d_1d( npti, nptidx(1:npti), wfx_pnd_1d (1:npti) , wfx_pnd ) 271 272 CALL tab_2d_1d( npti, nptidx(1:npti), diag_dvpn_mlt_1d (1:npti), diag_dvpn_mlt ) 273 CALL tab_2d_1d( npti, nptidx(1:npti), diag_dvpn_drn_1d (1:npti), diag_dvpn_drn ) 274 CALL tab_2d_1d( npti, nptidx(1:npti), diag_dvpn_lid_1d (1:npti), diag_dvpn_lid ) 275 CALL tab_2d_1d( npti, nptidx(1:npti), diag_dvpn_rnf_1d (1:npti), diag_dvpn_rnf ) 276 277 DO jl = 1, jpl 278 279 CALL tab_2d_1d( npti, nptidx(1:npti), a_i_1d (1:npti), a_i (:,:,jl) ) 280 CALL tab_2d_1d( npti, nptidx(1:npti), h_i_1d (1:npti), h_i (:,:,jl) ) 281 CALL tab_2d_1d( npti, nptidx(1:npti), t_su_1d (1:npti), t_su (:,:,jl) ) 282 CALL tab_2d_1d( npti, nptidx(1:npti), a_ip_1d (1:npti), a_ip (:,:,jl) ) 283 CALL tab_2d_1d( npti, nptidx(1:npti), h_ip_1d (1:npti), h_ip (:,:,jl) ) 284 CALL tab_2d_1d( npti, nptidx(1:npti), h_il_1d (1:npti), h_il (:,:,jl) ) 285 286 CALL tab_2d_1d( npti, nptidx(1:npti), dh_i_sum (1:npti), dh_i_sum_2d (:,:,jl) ) 287 CALL tab_2d_1d( npti, nptidx(1:npti), dh_s_mlt (1:npti), dh_s_mlt_2d (:,:,jl) ) 212 288 213 !--- Pond growing ---! 214 v_ip_1d(ji) = v_ip_1d(ji) + zdv_mlt 215 ! 216 !--- Lid melting ---! 217 IF( ln_pnd_lids ) v_il_1d(ji) = MAX( 0._wp, v_il_1d(ji) - zdv_mlt ) ! must be bounded by 0 218 ! 219 !--- mass flux ---! 220 IF( zdv_mlt > 0._wp ) THEN 221 zfac = zdv_mlt * rhow * r1_rdtice ! melt pond mass flux < 0 [kg.m-2.s-1] 222 wfx_pnd_1d(ji) = wfx_pnd_1d(ji) - zfac 223 ! 224 zdum = zfac / ( wfx_snw_sum_1d(ji) + wfx_sum_1d(ji) ) ! adjust ice/snow melting flux > 0 to balance melt pond flux 225 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) * (1._wp + zdum) 226 wfx_sum_1d(ji) = wfx_sum_1d(ji) * (1._wp + zdum) 227 ENDIF 228 229 !-------------------! 230 ! case ice freezing ! i.e. t_su_1d(ji) < (zTp+rt0) 231 !-------------------! 232 ! 233 zdT = MAX( zTp+rt0 - t_su_1d(ji), 0._wp ) 234 ! 235 !--- Pond contraction (due to refreezing) ---! 236 IF( ln_pnd_lids ) THEN 237 ! 238 !--- Lid growing and subsequent pond shrinking ---! 239 zdv_frz = 0.5_wp * MAX( 0._wp, -v_il_1d(ji) + & ! Flocco 2010 (eq. 5) solved implicitly as aH**2 + bH + c = 0 240 & SQRT( v_il_1d(ji)**2 + a_ip_1d(ji)**2 * 4._wp * rcnd_i * zdT * rdt_ice / (rLfus * rhow) ) ) ! max for roundoff errors 241 242 ! Lid growing 243 v_il_1d(ji) = MAX( 0._wp, v_il_1d(ji) + zdv_frz ) 244 245 ! Pond shrinking 246 v_ip_1d(ji) = MAX( 0._wp, v_ip_1d(ji) - zdv_frz ) 247 248 ELSE 249 ! Pond shrinking 250 v_ip_1d(ji) = v_ip_1d(ji) * EXP( 0.01_wp * zdT * z1_Tp ) ! Holland 2012 (eq. 6) 251 ENDIF 252 ! 253 !--- Set new pond area and depth ---! assuming linear relation between h_ip and a_ip_frac 254 ! v_ip = h_ip * a_ip 255 ! a_ip/a_i = a_ip_frac = h_ip / zaspect (cf Holland 2012, fitting SHEBA so that knowing v_ip we can distribute it to a_ip and h_ip) 256 a_ip_1d(ji) = MIN( a_i_1d(ji), SQRT( v_ip_1d(ji) * z1_aspect * a_i_1d(ji) ) ) ! make sure a_ip < a_i 257 h_ip_1d(ji) = zaspect * a_ip_1d(ji) / a_i_1d(ji) 258 259 !---------------! 260 ! Pond flushing ! 261 !---------------! 262 ! height of top of the pond above sea-level 263 zhp = ( h_i_1d(ji) * ( rau0 - rhoi ) + h_ip_1d(ji) * ( rau0 - rhow * a_ip_1d(ji) / a_i_1d(ji) ) ) * r1_rau0 289 DO jk = 1, nlay_i 290 CALL tab_2d_1d( npti, nptidx(1:npti), sz_i_1d(1:npti,jk), sz_i(:,:,jk,jl) ) 291 END DO 264 292 265 ! Calculate the permeability of the ice (Assur 1958, see Flocco 2010) 266 DO jk = 1, nlay_i 267 zsbr = - 1.2_wp & 268 & - 21.8_wp * ( t_i_1d(ji,jk) - rt0 ) & 269 & - 0.919_wp * ( t_i_1d(ji,jk) - rt0 )**2 & 270 & - 0.0178_wp * ( t_i_1d(ji,jk) - rt0 )**3 271 ztmp(jk) = sz_i_1d(ji,jk) / zsbr 272 END DO 273 zperm = MAX( 0._wp, 3.e-08_wp * MINVAL(ztmp)**3 ) 274 275 ! Do the drainage using Darcy's law 276 zdv_flush = -zperm * rau0 * grav * zhp * rdt_ice / (zvisc * h_i_1d(ji)) * a_ip_1d(ji) * rn_pnd_flush ! tunable rn_pnd_flush from hunke et al. (2013) 277 zdv_flush = MAX( zdv_flush, -v_ip_1d(ji) ) 278 v_ip_1d(ji) = v_ip_1d(ji) + zdv_flush 279 280 !--- Set new pond area and depth ---! assuming linear relation between h_ip and a_ip_frac 281 a_ip_1d(ji) = MIN( a_i_1d(ji), SQRT( v_ip_1d(ji) * z1_aspect * a_i_1d(ji) ) ) ! make sure a_ip < a_i 282 h_ip_1d(ji) = zaspect * a_ip_1d(ji) / a_i_1d(ji) 283 284 !--- Corrections and lid thickness ---! 285 IF( ln_pnd_lids ) THEN 286 !--- retrieve lid thickness from volume ---! 287 IF( a_ip_1d(ji) > epsi10 ) THEN ; h_il_1d(ji) = v_il_1d(ji) / a_ip_1d(ji) 288 ELSE ; h_il_1d(ji) = 0._wp 289 ENDIF 290 !--- remove ponds if lids are much larger than ponds ---! 291 IF ( h_il_1d(ji) > h_ip_1d(ji) * 10._wp ) THEN 293 !----------------------------------------------------------------------------------------------- 294 ! Go for ponds 295 !----------------------------------------------------------------------------------------------- 296 297 DO ji = 1, npti 298 ! !----------------------------------------------------! 299 IF( h_i_1d(ji) < rn_himin .OR. a_i_1d(ji) < epsi10 ) THEN ! Case ice thickness < rn_himin or tiny ice fraction ! 300 ! !----------------------------------------------------! 301 !--- Remove ponds on thin ice or tiny ice fractions 292 302 a_ip_1d(ji) = 0._wp 293 303 h_ip_1d(ji) = 0._wp 294 304 h_il_1d(ji) = 0._wp 305 ! !--------------------------------! 306 ELSE ! Case ice thickness >= rn_himin ! 307 ! !--------------------------------! 308 v_ip_1d(ji) = h_ip_1d(ji) * a_ip_1d(ji) ! retrieve volume from thickness 309 v_il_1d(ji) = h_il_1d(ji) * a_ip_1d(ji) 310 ! 311 !------------------! 312 ! case ice melting ! 313 !------------------! 314 ! 315 !--- available meltwater for melt ponding ---! 316 zdum = -( dh_i_sum(ji)*rhoi + dh_s_mlt(ji)*rhos ) * z1_rhow * a_i_1d(ji) 317 zfr_mlt = rn_apnd_min + ( rn_apnd_max - rn_apnd_min ) * at_i_1d(ji) ! = ( 1 - r ) = fraction of melt water that is not flushed 318 zdv_mlt = MAX( 0._wp, zfr_mlt * zdum ) ! max for roundoff errors? 319 320 diag_dvpn_mlt_1d(ji) = diag_dvpn_mlt_1d(ji) + zdum * r1_rdtice ! surface melt input diag 321 diag_dvpn_rnf_1d(ji) = diag_dvpn_rnf_1d(ji) + ( 1. - zfr_mlt ) * zdum * r1_rdtice ! runoff diag 322 ! 323 !--- overflow ---! 324 ! 325 ! 1) area driven overflow 326 ! 327 ! If pond area exceeds zfr_mlt * a_i_1d(ji) then reduce the pond water volume 328 ! a_ip_max = zfr_mlt * a_i 329 ! => from zaspect = h_ip / (a_ip / a_i), set v_ip_max as: 330 zv_ip_max = zfr_mlt**2 * a_i_1d(ji) * zaspect 331 zvold = zdv_mlt 332 zdv_mlt = MAX( 0._wp, MIN( zdv_mlt, zv_ip_max - v_ip_1d(ji) ) ) 333 334 ! 335 ! 2) depth driven overflow 336 ! 337 ! If pond depth exceeds half the ice thickness then reduce the pond volume 338 ! h_ip_max = 0.5 * h_i 339 ! => from zaspect = h_ip / (a_ip / a_i), set v_ip_max as: 340 zv_ip_max = z1_aspect * a_i_1d(ji) * 0.25 * h_i_1d(ji) * h_i_1d(ji) ! MV dimensions are wrong here or comment is unclear 341 zdv_mlt = MAX( 0._wp, MIN( zdv_mlt, zv_ip_max - v_ip_1d(ji) ) ) 342 diag_dvpn_rnf_1d(ji) = diag_dvpn_rnf_1d(ji) + ( zdv_mlt - zvold ) * r1_rdtice ! runoff diag - overflow contribution 343 344 !--- Pond growing ---! 345 v_ip_1d(ji) = v_ip_1d(ji) + zdv_mlt 346 ! 347 !--- Lid melting ---! 348 IF( ln_pnd_lids ) v_il_1d(ji) = MAX( 0._wp, v_il_1d(ji) - zdv_mlt ) ! must be bounded by 0 349 ! 350 !--- mass flux ---! 351 ! MV I would recommend to remove that 352 ! Since melt ponds carry no freshwater there is no point in modifying water fluxes 353 354 IF( zdv_mlt > 0._wp ) THEN 355 zfac = zdv_mlt * rhow * r1_rdtice ! melt pond mass flux < 0 [kg.m-2.s-1] 356 wfx_pnd_1d(ji) = wfx_pnd_1d(ji) - zfac 357 ! 358 zdum = zfac / ( wfx_snw_sum_1d(ji) + wfx_sum_1d(ji) ) ! adjust ice/snow melting flux > 0 to balance melt pond flux 359 wfx_snw_sum_1d(ji) = wfx_snw_sum_1d(ji) * (1._wp + zdum) 360 wfx_sum_1d(ji) = wfx_sum_1d(ji) * (1._wp + zdum) 361 ENDIF 362 363 !-------------------! 364 ! case ice freezing ! i.e. t_su_1d(ji) < (zTp+rt0) 365 !-------------------! 366 ! 367 zdT = MAX( zTp+rt0 - t_su_1d(ji), 0._wp ) 368 ! 369 !--- Pond contraction (due to refreezing) ---! 370 zvold = v_ip_1d(ji) ! for diag 371 372 IF( ln_pnd_lids ) THEN 373 ! 374 !--- Lid growing and subsequent pond shrinking ---! 375 zdv_frz = 0.5_wp * MAX( 0._wp, -v_il_1d(ji) + & ! Flocco 2010 (eq. 5) solved implicitly as aH**2 + bH + c = 0 376 & SQRT( v_il_1d(ji)**2 + a_ip_1d(ji)**2 * 4._wp * rcnd_i * zdT * rdt_ice / (rLfus * rhow) ) ) ! max for roundoff errors 377 378 ! Lid growing 379 v_il_1d(ji) = MAX( 0._wp, v_il_1d(ji) + zdv_frz ) 380 381 ! Pond shrinking 382 v_ip_1d(ji) = MAX( 0._wp, v_ip_1d(ji) - zdv_frz ) 383 384 ELSE 385 386 ! Pond shrinking 387 v_ip_1d(ji) = v_ip_1d(ji) * EXP( 0.01_wp * zdT * z1_Tp ) ! Holland 2012 (eq. 6) 388 389 ENDIF 390 391 diag_dvpn_lid_1d(ji) = diag_dvpn_lid_1d(ji) + ( v_ip_1d(ji) - zvold ) * r1_rdtice ! shrinking counted as lid diagnostic 392 393 ! 394 !--- Set new pond area and depth ---! assuming linear relation between h_ip and a_ip_frac 395 ! v_ip = h_ip * a_ip 396 ! a_ip/a_i = a_ip_frac = h_ip / zaspect (cf Holland 2012, fitting SHEBA so that knowing v_ip we can distribute it to a_ip and h_ip) 397 a_ip_1d(ji) = MIN( a_i_1d(ji), SQRT( v_ip_1d(ji) * z1_aspect * a_i_1d(ji) ) ) ! make sure a_ip < a_i 398 h_ip_1d(ji) = zaspect * a_ip_1d(ji) / a_i_1d(ji) 399 400 !------------------------------------------------! 401 ! Pond drainage through brine network (flushing) ! 402 !------------------------------------------------! 403 ! height of top of the pond above sea-level 404 zhp = ( h_i_1d(ji) * ( rau0 - rhoi ) + h_ip_1d(ji) * ( rau0 - rhow * a_ip_1d(ji) / a_i_1d(ji) ) ) * r1_rau0 405 406 ! Calculate the permeability of the ice (Assur 1958, see Flocco 2010) 407 DO jk = 1, nlay_i 408 ! MV Assur is inconsistent with SI3 409 zsbr = - 1.2_wp & 410 & - 21.8_wp * ( t_i_1d(ji,jk) - rt0 ) & 411 & - 0.919_wp * ( t_i_1d(ji,jk) - rt0 )**2 & 412 & - 0.0178_wp * ( t_i_1d(ji,jk) - rt0 )**3 413 ! MV linear expression more consistent & simpler zsbr = - ( t_i_1d(ji,jk) - rt0 ) / rTmlt 414 ztmp(jk) = sz_i_1d(ji,jk) / zsbr 415 END DO 416 zperm = MAX( 0._wp, 3.e-08_wp * MINVAL(ztmp)**3 ) 417 418 ! Do the drainage using Darcy's law 419 zdv_flush = - zperm * rau0 * grav * zhp * rdt_ice / (zvisc * h_i_1d(ji)) * a_ip_1d(ji) * rn_pnd_flush 420 zdv_flush = MAX( zdv_flush, -v_ip_1d(ji) ) 421 ! zdv_flush = 0._wp ! MV remove pond drainage for now 422 v_ip_1d(ji) = v_ip_1d(ji) + zdv_flush 423 424 diag_dvpn_drn_1d(ji) = diag_dvpn_drn_1d(ji) + zdv_flush * r1_rdtice ! shrinking counted as lid diagnostic 425 426 ! MV --- why pond drainage does not give back water into freshwater flux ? 427 !--- Set new pond area and depth ---! assuming linear relation between h_ip and a_ip_frac 428 a_ip_1d(ji) = MIN( a_i_1d(ji), SQRT( v_ip_1d(ji) * z1_aspect * a_i_1d(ji) ) ) ! make sure a_ip < a_i 429 h_ip_1d(ji) = zaspect * a_ip_1d(ji) / a_i_1d(ji) 430 431 !--- Corrections and lid thickness ---! 432 IF( ln_pnd_lids ) THEN 433 !--- retrieve lid thickness from volume ---! 434 IF( a_ip_1d(ji) > epsi10 ) THEN ; h_il_1d(ji) = v_il_1d(ji) / a_ip_1d(ji) 435 ELSE ; h_il_1d(ji) = 0._wp 436 ENDIF 437 !--- remove ponds if lids are much larger than ponds ---! 438 IF ( h_il_1d(ji) > h_ip_1d(ji) * 10._wp ) THEN 439 a_ip_1d(ji) = 0._wp 440 h_ip_1d(ji) = 0._wp 441 h_il_1d(ji) = 0._wp 442 ENDIF 443 ENDIF 444 ! 295 445 ENDIF 446 447 END DO ! ji 448 449 !----------------------------------------------------------------------------------------------- 450 ! Retrieve 2D arrays 451 !----------------------------------------------------------------------------------------------- 452 453 v_ip_1d(1:npti) = h_ip_1d(1:npti) * a_ip_1d(1:npti) 454 v_il_1d(1:npti) = h_il_1d(1:npti) * a_ip_1d(1:npti) 455 456 CALL tab_1d_2d( npti, nptidx(1:npti), a_ip_1d (1:npti), a_ip (:,:,jl) ) 457 CALL tab_1d_2d( npti, nptidx(1:npti), h_ip_1d (1:npti), h_ip (:,:,jl) ) 458 CALL tab_1d_2d( npti, nptidx(1:npti), h_il_1d (1:npti), h_il (:,:,jl) ) 459 CALL tab_1d_2d( npti, nptidx(1:npti), v_ip_1d (1:npti), v_ip (:,:,jl) ) 460 CALL tab_1d_2d( npti, nptidx(1:npti), v_il_1d (1:npti), v_il (:,:,jl) ) 461 DO jk = 1, nlay_i 462 CALL tab_1d_2d( npti, nptidx(1:npti), sz_i_1d(1:npti,jk), sz_i(:,:,jk,jl) ) 463 END DO 464 465 END DO ! ji 466 467 CALL tab_1d_2d( npti, nptidx(1:npti), wfx_snw_sum_1d(1:npti), wfx_snw_sum ) 468 CALL tab_1d_2d( npti, nptidx(1:npti), wfx_sum_1d (1:npti), wfx_sum ) 469 CALL tab_1d_2d( npti, nptidx(1:npti), wfx_pnd_1d (1:npti), wfx_pnd ) 470 471 CALL tab_1d_2d( npti, nptidx(1:npti), diag_dvpn_mlt_1d (1:npti), diag_dvpn_mlt ) 472 CALL tab_1d_2d( npti, nptidx(1:npti), diag_dvpn_drn_1d (1:npti), diag_dvpn_drn ) 473 CALL tab_1d_2d( npti, nptidx(1:npti), diag_dvpn_lid_1d (1:npti), diag_dvpn_lid ) 474 CALL tab_1d_2d( npti, nptidx(1:npti), diag_dvpn_rnf_1d (1:npti), diag_dvpn_rnf ) 475 476 ! 477 ENDIF 478 479 END SUBROUTINE pnd_LEV 480 481 SUBROUTINE pnd_TOPO 482 483 !!------------------------------------------------------------------- 484 !! *** ROUTINE pnd_TOPO *** 485 !! 486 !! ** Purpose : Compute melt pond evolution based on the ice 487 !! topography inferred from the ice thickness distribution 488 !! 489 !! ** Method : This code is initially based on Flocco and Feltham 490 !! (2007) and Flocco et al. (2010). 491 !! 492 !! - Calculate available pond water base on surface meltwater 493 !! - Redistribute water as a function of topography, drain water 494 !! - Exchange water with the lid 495 !! 496 !! ** Tunable parameters : 497 !! 498 !! ** Note : 499 !! 500 !! ** References 501 !! 502 !! Flocco, D. and D. L. Feltham, 2007. A continuum model of melt pond 503 !! evolution on Arctic sea ice. J. Geophys. Res. 112, C08016, doi: 504 !! 10.1029/2006JC003836. 505 !! 506 !! Flocco, D., D. L. Feltham and A. K. Turner, 2010. Incorporation of 507 !! a physically based melt pond scheme into the sea ice component of a 508 !! climate model. J. Geophys. Res. 115, C08012, 509 !! doi: 10.1029/2009JC005568. 510 !! 511 !!------------------------------------------------------------------- 512 513 ! local variables 514 REAL(wp) :: & 515 zdHui, & ! change in thickness of ice lid (m) 516 zomega, & ! conduction 517 zdTice, & ! temperature difference across ice lid (C) 518 zdvice, & ! change in ice volume (m) 519 zTavg, & ! mean surface temperature across categories (C) 520 zfsurf, & ! net heat flux, excluding conduction and transmitted radiation (W/m2) 521 zTp, & ! pond freezing temperature (C) 522 zrhoi_L, & ! volumetric latent heat of sea ice (J/m^3) 523 zfr_mlt, & ! fraction and volume of available meltwater retained for melt ponding 524 z1_rhow, & ! inverse water density 525 zv_pnd , & ! volume of meltwater contributing to ponds 526 zv_mlt ! total amount of meltwater produced 527 528 REAL(wp), DIMENSION(jpi,jpj) :: zvolp_ini, & !! total melt pond water available before redistribution and drainage 529 zvolp , & !! total melt pond water volume 530 zvolp_res !! remaining melt pond water available after drainage 531 532 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z1_a_i 533 534 INTEGER :: ji, jj, jk, jl ! loop indices 535 536 INTEGER :: i_test 537 538 ! Note 539 ! equivalent for CICE translation 540 ! a_ip -> apond 541 ! a_ip_frac -> apnd 542 543 !--------------------------------------------------------------- 544 ! Initialise 545 !--------------------------------------------------------------- 546 547 ! Parameters & constants (move to parameters) 548 zrhoi_L = rhoi * rLfus ! volumetric latent heat (J/m^3) 549 zTp = rt0 - 0.15_wp ! pond freezing point, slightly below 0C (ponds are bid saline) 550 z1_rhow = 1._wp / rhow 551 552 ! Set required ice variables (hard-coded here for now) 553 ! zfpond(:,:) = 0._wp ! contributing freshwater flux (?) 554 555 at_i (:,:) = SUM( a_i (:,:,:), dim=3 ) ! ice fraction 556 vt_i (:,:) = SUM( v_i (:,:,:), dim=3 ) ! volume per grid area 557 vt_ip(:,:) = SUM( v_ip(:,:,:), dim=3 ) ! pond volume per grid area 558 vt_il(:,:) = SUM( v_il(:,:,:), dim=3 ) ! lid volume per grid area 559 560 ! thickness 561 WHERE( a_i(:,:,:) > epsi20 ) ; z1_a_i(:,:,:) = 1._wp / a_i(:,:,:) 562 ELSEWHERE ; z1_a_i(:,:,:) = 0._wp 563 END WHERE 564 h_i(:,:,:) = v_i (:,:,:) * z1_a_i(:,:,:) 565 566 !--------------------------------------------------------------- 567 ! Change 2D to 1D 568 !--------------------------------------------------------------- 569 ! MV 570 ! a less computing-intensive version would have 2D-1D passage here 571 ! use what we have in iceitd.F90 (incremental remapping) 572 573 !-------------------------------------------------------------- 574 ! Collect total available pond water volume 575 !-------------------------------------------------------------- 576 ! Assuming that meltwater (+rain in principle) runsoff the surface 577 ! Holland et al (2012) suggest that the fraction of runoff decreases with total ice fraction 578 ! I cite her words, they are very talkative 579 ! "grid cells with very little ice cover (and hence more open water area) 580 ! have a higher runoff fraction to rep- resent the greater proximity of ice to open water." 581 ! "This results in the same runoff fraction r for each ice category within a grid cell" 582 583 zvolp(:,:) = 0._wp 584 585 DO jl = 1, jpl 586 DO jj = 1, jpj 587 DO ji = 1, jpi 588 589 IF ( a_i(ji,jj,jl) > epsi10 ) THEN 590 591 !--- Available and contributing meltwater for melt ponding ---! 592 zv_mlt = - ( dh_i_sum_2d(ji,jj,jl) * rhoi + dh_s_mlt_2d(ji,jj,jl) * rhos ) & ! available volume of surface melt water per grid area 593 & * z1_rhow * a_i(ji,jj,jl) 594 ! MV -> could move this directly in ice_thd_dh and get an array (ji,jj,jl) for surface melt water volume per grid area 595 zfr_mlt = rn_apnd_min + ( rn_apnd_max - rn_apnd_min ) * at_i(ji,jj) ! fraction of surface meltwater going to ponds 596 zv_pnd = zfr_mlt * zv_mlt ! contributing meltwater volume for category jl 597 598 diag_dvpn_mlt(ji,jj) = diag_dvpn_mlt(ji,jj) + zv_mlt * r1_rdtice ! diags 599 diag_dvpn_rnf(ji,jj) = diag_dvpn_rnf(ji,jj) + ( 1. - zfr_mlt ) * zv_mlt * r1_rdtice 600 601 !--- Create possible new ponds 602 ! if pond does not exist, create new pond over full ice area 603 !IF ( a_ip_frac(ji,jj,jl) < epsi10 ) THEN 604 IF ( a_ip(ji,jj,jl) < epsi10 ) THEN 605 a_ip(ji,jj,jl) = a_i(ji,jj,jl) 606 a_ip_frac(ji,jj,jl) = 1.0_wp ! pond fraction of sea ice (apnd for CICE) 607 ENDIF 608 609 !--- Deepen existing ponds with no change in pond fraction, before redistribution and drainage 610 v_ip(ji,jj,jl) = v_ip(ji,jj,jl) + zv_pnd ! use pond water to increase thickness 611 h_ip(ji,jj,jl) = v_ip(ji,jj,jl) / a_ip(ji,jj,jl) 612 613 !--- Total available pond water volume (pre-existing + newly produced)j 614 zvolp(ji,jj) = zvolp(ji,jj) + v_ip(ji,jj,jl) 615 ! zfpond(ji,jj) = zfpond(ji,jj) + zpond * a_ip_frac(ji,jj,jl) ! useless for now 616 617 ENDIF ! a_i 618 619 END DO! jl 620 END DO ! jj 621 END DO ! ji 622 623 zvolp_ini(:,:) = zvolp(:,:) 624 625 !-------------------------------------------------------------- 626 ! Redistribute and drain water from ponds 627 !-------------------------------------------------------------- 628 CALL ice_thd_pnd_area( zvolp, zvolp_res ) 629 630 !-------------------------------------------------------------- 631 ! Melt pond lid growth and melt 632 !-------------------------------------------------------------- 633 634 IF( ln_pnd_lids ) THEN 635 636 DO jj = 1, jpj 637 DO ji = 1, jpi 638 639 IF ( at_i(ji,jj) > 0.01 .AND. hm_i(ji,jj) > zhi_min .AND. zvolp_ini(ji,jj) > zvp_min * at_i(ji,jj) ) THEN 640 641 !-------------------------- 642 ! Pond lid growth and melt 643 !-------------------------- 644 ! Mean surface temperature 645 zTavg = 0._wp 646 DO jl = 1, jpl 647 zTavg = zTavg + t_su(ji,jj,jl)*a_i(ji,jj,jl) 648 END DO 649 zTavg = zTavg / a_i(ji,jj,jl) !!! could get a division by zero here 650 651 DO jl = 1, jpl-1 652 653 IF ( v_il(ji,jj,jl) > epsi10 ) THEN 654 655 !---------------------------------------------------------------- 656 ! Lid melting: floating upper ice layer melts in whole or part 657 !---------------------------------------------------------------- 658 ! Use Tsfc for each category 659 IF ( t_su(ji,jj,jl) > zTp ) THEN 660 661 zdvice = MIN( dh_i_sum_2d(ji,jj,jl)*a_ip(ji,jj,jl), v_il(ji,jj,jl) ) 662 663 IF ( zdvice > epsi10 ) THEN 664 665 v_il (ji,jj,jl) = v_il (ji,jj,jl) - zdvice 666 v_ip(ji,jj,jl) = v_ip(ji,jj,jl) + zdvice ! MV: not sure i understand dh_i_sum seems counted twice - 667 ! as it is already counted in surface melt 668 ! zvolp(ji,jj) = zvolp(ji,jj) + zdvice ! pointless to calculate total volume (done in icevar) 669 ! zfpond(ji,jj) = fpond(ji,jj) + zdvice ! pointless to follow fw budget (ponds have no fw) 670 671 IF ( v_il(ji,jj,jl) < epsi10 .AND. v_ip(ji,jj,jl) > epsi10) THEN 672 ! ice lid melted and category is pond covered 673 v_ip(ji,jj,jl) = v_ip(ji,jj,jl) + v_il(ji,jj,jl) 674 ! zfpond(ji,jj) = zfpond (ji,jj) + v_il(ji,jj,jl) 675 v_il(ji,jj,jl) = 0._wp 676 ENDIF 677 h_ip(ji,jj,jl) = v_ip(ji,jj,jl) / a_ip(ji,jj,jl) !!! could get a division by zero here 678 679 diag_dvpn_lid(ji,jj) = diag_dvpn_lid(ji,jj) + zdvice ! diag 680 681 ENDIF 682 683 !---------------------------------------------------------------- 684 ! Freeze pre-existing lid 685 !---------------------------------------------------------------- 686 687 ELSE IF ( v_ip(ji,jj,jl) > epsi10 ) THEN ! Tsfcn(i,j,n) <= Tp 688 689 ! differential growth of base of surface floating ice layer 690 ! zdTice = MAX( - t_su(ji,jj,jl) - zTd , 0._wp ) ! > 0 ! line valid for temp in C 691 zdTice = MAX( - ( t_su(ji,jj,jl) - zTd ) , 0._wp ) ! > 0 692 zomega = rcnd_i * zdTice / zrhoi_L 693 zdHui = SQRT( 2._wp * zomega * rdt_ice + ( v_il(ji,jj,jl) / a_i(ji,jj,jl) )**2 ) & 694 - v_il(ji,jj,jl) / a_i(ji,jj,jl) 695 zdvice = min( zdHui*a_ip(ji,jj,jl) , v_ip(ji,jj,jl) ) 696 697 IF ( zdvice > epsi10 ) THEN 698 v_il (ji,jj,jl) = v_il(ji,jj,jl) + zdvice 699 v_ip(ji,jj,jl) = v_ip(ji,jj,jl) - zdvice 700 ! zvolp(ji,jj) = zvolp(ji,jj) - zdvice 701 ! zfpond(ji,jj) = zfpond(ji,jj) - zdvice 702 h_ip(ji,jj,jl) = v_ip(ji,jj,jl) / a_ip(ji,jj,jl) 703 704 diag_dvpn_lid(ji,jj) = diag_dvpn_lid(ji,jj) - zdvice ! diag 705 706 ENDIF 707 708 ENDIF ! Tsfcn(i,j,n) 709 710 !---------------------------------------------------------------- 711 ! Freeze new lids 712 !---------------------------------------------------------------- 713 ! upper ice layer begins to form 714 ! note: albedo does not change 715 716 ELSE ! v_il < epsi10 717 718 ! thickness of newly formed ice 719 ! the surface temperature of a meltpond is the same as that 720 ! of the ice underneath (0C), and the thermodynamic surface 721 ! flux is the same 722 723 !!! we need net surface energy flux, excluding conduction 724 !!! fsurf is summed over categories in CICE 725 !!! we have the category-dependent flux, let us use it ? 726 zfsurf = qns_ice(ji,jj,jl) + qsr_ice(ji,jj,jl) 727 zdHui = MAX ( -zfsurf * rdt_ice / zrhoi_L , 0._wp ) 728 zdvice = MIN ( zdHui * a_ip(ji,jj,jl) , v_ip(ji,jj,jl) ) 729 730 IF ( zdvice > epsi10 ) THEN 731 732 v_il (ji,jj,jl) = v_il(ji,jj,jl) + zdvice 733 v_ip(ji,jj,jl) = v_ip(ji,jj,jl) - zdvice 734 diag_dvpn_lid(ji,jj) = diag_dvpn_lid(ji,jj) - zdvice ! diag 735 736 ! zvolp(ji,jj) = zvolp(ji,jj) - zdvice 737 ! zfpond(ji,jj) = zfpond(ji,jj) - zdvice 738 739 h_ip(ji,jj,jl) = v_ip(ji,jj,jl) / a_ip(ji,jj,jl) ! MV - in principle, this is useless as h_ip is computed in icevar 740 ENDIF 741 742 ENDIF ! v_il 743 744 END DO ! jl 745 746 ELSE ! remove ponds on thin ice 747 748 v_ip(ji,jj,:) = 0._wp 749 v_il(ji,jj,:) = 0._wp 750 ! zfpond(ji,jj) = zfpond(ji,jj)- zvolp(ji,jj) 751 ! zvolp(ji,jj) = 0._wp 752 296 753 ENDIF 297 ! 298 ENDIF 299 300 END DO 301 ! 302 END SUBROUTINE pnd_LEV 303 754 755 END DO ! jj 756 END DO ! ji 757 758 ENDIF ! ln_pnd_lids 759 760 !--------------------------------------------------------------- 761 ! Clean-up variables (probably duplicates what icevar would do) 762 !--------------------------------------------------------------- 763 ! MV comment 764 ! In the ideal world, the lines above should update only v_ip, a_ip, v_il 765 ! icevar should recompute all other variables (if needed at all) 766 767 DO jl = 1, jpl 768 769 DO jj = 1, jpj 770 DO ji = 1, jpi 771 772 ! ! zap lids on small ponds 773 ! IF ( a_i(ji,jj,jl) > epsi10 .AND. v_ip(ji,jj,jl) < epsi10 & 774 ! .AND. v_il(ji,jj,jl) > epsi10) THEN 775 ! v_il(ji,jj,jl) = 0._wp ! probably uselesss now since we get zap_small 776 ! ENDIF 777 778 ! recalculate equivalent pond variables 779 IF ( a_ip(ji,jj,jl) > epsi10) THEN 780 h_ip(ji,jj,jl) = v_ip(ji,jj,jl) / a_i(ji,jj,jl) 781 a_ip_frac(ji,jj,jl) = a_ip(ji,jj,jl) / a_i(ji,jj,jl) ! MV in principle, useless as computed in icevar 782 h_il(ji,jj,jl) = v_il(ji,jj,jl) / a_ip(ji,jj,jl) ! MV in principle, useless as computed in icevar 783 ENDIF 784 ! h_ip(ji,jj,jl) = 0._wp ! MV in principle, useless as computed in icevar 785 ! h_il(ji,jj,jl) = 0._wp ! MV in principle, useless as omputed in icevar 786 ! ENDIF 787 788 END DO ! ji 789 END DO ! jj 790 791 END DO ! jl 792 793 END SUBROUTINE pnd_TOPO 794 795 796 SUBROUTINE ice_thd_pnd_area( zvolp , zdvolp ) 797 798 !!------------------------------------------------------------------- 799 !! *** ROUTINE ice_thd_pnd_area *** 800 !! 801 !! ** Purpose : Given the total volume of available pond water, 802 !! redistribute and drain water 803 !! 804 !! ** Method 805 !! 806 !-----------| 807 ! | 808 ! |-----------| 809 !___________|___________|______________________________________sea-level 810 ! | | 811 ! | |---^--------| 812 ! | | | | 813 ! | | | |-----------| |------- 814 ! | | | alfan | | | 815 ! | | | | |--------------| 816 ! | | | | | | 817 !---------------------------v------------------------------------------- 818 ! | | ^ | | | 819 ! | | | | |--------------| 820 ! | | | betan | | | 821 ! | | | |-----------| |------- 822 ! | | | | 823 ! | |---v------- | 824 ! | | 825 ! |-----------| 826 ! | 827 !-----------| 828 ! 829 !! 830 !!------------------------------------------------------------------ 831 832 REAL (wp), DIMENSION(jpi,jpj), INTENT(INOUT) :: & 833 zvolp, & ! total available pond water 834 zdvolp ! remaining meltwater after redistribution 835 836 INTEGER :: & 837 ns, & 838 m_index, & 839 permflag 840 841 REAL (wp), DIMENSION(jpl) :: & 842 hicen, & 843 hsnon, & 844 asnon, & 845 alfan, & 846 betan, & 847 cum_max_vol, & 848 reduced_aicen 849 850 REAL (wp), DIMENSION(0:jpl) :: & 851 cum_max_vol_tmp 852 853 REAL (wp) :: & 854 hpond, & 855 drain, & 856 floe_weight, & 857 pressure_head, & 858 hsl_rel, & 859 deltah, & 860 perm, & 861 msno 862 863 REAL (wp), parameter :: & 864 viscosity = 1.79e-3_wp ! kinematic water viscosity in kg/m/s 865 866 INTEGER :: ji, jj, jk, jl ! loop indices 867 868 a_ip(:,:,:) = 0._wp 869 h_ip(:,:,:) = 0._wp 870 871 DO jj = 1, jpj 872 DO ji = 1, jpi 873 874 IF ( at_i(ji,jj) > 0.01 .AND. hm_i(ji,jj) > zhi_min .AND. zvolp(ji,jj) > zvp_min * at_i(ji,jj) ) THEN 875 876 !------------------------------------------------------------------- 877 ! initialize 878 !------------------------------------------------------------------- 879 880 DO jl = 1, jpl 881 882 !---------------------------------------- 883 ! compute the effective snow fraction 884 !---------------------------------------- 885 886 IF (a_i(ji,jj,jl) < epsi10) THEN 887 hicen(jl) = 0._wp 888 hsnon(jl) = 0._wp 889 reduced_aicen(jl) = 0._wp 890 asnon(jl) = 0._wp !js: in CICE 5.1.2: make sense as the compiler may not initiate the variables 891 ELSE 892 hicen(jl) = v_i(ji,jj,jl) / a_i(ji,jj,jl) 893 hsnon(jl) = v_s(ji,jj,jl) / a_i(ji,jj,jl) 894 reduced_aicen(jl) = 1._wp ! n=jpl 895 896 !js: initial code in NEMO_DEV 897 !IF (n < jpl) reduced_aicen(jl) = aicen(jl) & 898 ! * (-0.024_wp*hicen(jl) + 0.832_wp) 899 900 !js: from CICE 5.1.2: this limit reduced_aicen to 0.2 when hicen is too large 901 IF (jl < jpl) reduced_aicen(jl) = a_i(ji,jj,jl) & 902 * max(0.2_wp,(-0.024_wp*hicen(jl) + 0.832_wp)) 903 904 asnon(jl) = reduced_aicen(jl) ! effective snow fraction (empirical) 905 ! MV should check whether this makes sense to have the same effective snow fraction in here 906 ! OLI: it probably doesn't 907 END IF 908 909 ! This choice for alfa and beta ignores hydrostatic equilibium of categories. 910 ! Hydrostatic equilibium of the entire ITD is accounted for below, assuming 911 ! a surface topography implied by alfa=0.6 and beta=0.4, and rigidity across all 912 ! categories. alfa and beta partition the ITD - they are areas not thicknesses! 913 ! Multiplying by hicen, alfan and betan (below) are thus volumes per unit area. 914 ! Here, alfa = 60% of the ice area (and since hice is constant in a category, 915 ! alfan = 60% of the ice volume) in each category lies above the reference line, 916 ! and 40% below. Note: p6 is an arbitrary choice, but alfa+beta=1 is required. 917 918 ! MV: 919 ! Note that this choice is not in the original FF07 paper and has been adopted in CICE 920 ! No reason why is explained in the doc, but I guess there is a reason. I'll try to investigate, maybe 921 922 ! Where does that choice come from ? => OLI : Coz' Chuck Norris said so... 923 924 alfan(jl) = 0.6 * hicen(jl) 925 betan(jl) = 0.4 * hicen(jl) 926 927 cum_max_vol(jl) = 0._wp 928 cum_max_vol_tmp(jl) = 0._wp 929 930 END DO ! jpl 931 932 cum_max_vol_tmp(0) = 0._wp 933 drain = 0._wp 934 zdvolp(ji,jj) = 0._wp 935 936 !---------------------------------------------------------- 937 ! Drain overflow water, update pond fraction and volume 938 !---------------------------------------------------------- 939 940 !-------------------------------------------------------------------------- 941 ! the maximum amount of water that can be contained up to each ice category 942 !-------------------------------------------------------------------------- 943 ! If melt ponds are too deep to be sustainable given the ITD (OVERFLOW) 944 ! Then the excess volume cum_max_vol(jl) drains out of the system 945 ! It should be added to wfx_pnd_out 946 947 DO jl = 1, jpl-1 ! last category can not hold any volume 948 949 IF (alfan(jl+1) >= alfan(jl) .AND. alfan(jl+1) > 0._wp ) THEN 950 951 ! total volume in level including snow 952 cum_max_vol_tmp(jl) = cum_max_vol_tmp(jl-1) + & 953 (alfan(jl+1) - alfan(jl)) * sum(reduced_aicen(1:jl)) 954 955 ! subtract snow solid volumes from lower categories in current level 956 DO ns = 1, jl 957 cum_max_vol_tmp(jl) = cum_max_vol_tmp(jl) & 958 - rhos/rhow * & ! free air fraction that can be filled by water 959 asnon(ns) * & ! effective areal fraction of snow in that category 960 max(min(hsnon(ns)+alfan(ns)-alfan(jl), alfan(jl+1)-alfan(jl)), 0._wp) 961 END DO 962 963 ELSE ! assume higher categories unoccupied 964 cum_max_vol_tmp(jl) = cum_max_vol_tmp(jl-1) 965 END IF 966 !IF (cum_max_vol_tmp(jl) < z0) THEN 967 ! CALL abort_ice('negative melt pond volume') 968 !END IF 969 END DO 970 cum_max_vol_tmp(jpl) = cum_max_vol_tmp(jpl-1) ! last category holds no volume 971 cum_max_vol (1:jpl) = cum_max_vol_tmp(1:jpl) 972 973 !---------------------------------------------------------------- 974 ! is there more meltwater than can be held in the floe? 975 !---------------------------------------------------------------- 976 IF (zvolp(ji,jj) >= cum_max_vol(jpl)) THEN 977 drain = zvolp(ji,jj) - cum_max_vol(jpl) + epsi10 978 zvolp(ji,jj) = zvolp(ji,jj) - drain ! update meltwater volume available 979 980 diag_dvpn_rnf(ji,jj) = - drain ! diag - overflow counted in the runoff part (arbitrary choice) 981 982 zdvolp(ji,jj) = drain ! this is the drained water 983 IF (zvolp(ji,jj) < epsi10) THEN 984 zdvolp(ji,jj) = zdvolp(ji,jj) + zvolp(ji,jj) 985 zvolp(ji,jj) = 0._wp 986 END IF 987 END IF 988 989 ! height and area corresponding to the remaining volume 990 ! routine leaves zvolp unchanged 991 CALL ice_thd_pnd_depth(reduced_aicen, asnon, hsnon, alfan, zvolp(ji,jj), cum_max_vol, hpond, m_index) 992 993 DO jl = 1, m_index 994 !h_ip(jl) = hpond - alfan(jl) + alfan(1) ! here oui choulde update 995 ! ! volume instead, no ? 996 h_ip(ji,jj,jl) = max((hpond - alfan(jl) + alfan(1)), 0._wp) !js: from CICE 5.1.2 997 a_ip(ji,jj,jl) = reduced_aicen(jl) 998 ! in practise, pond fraction depends on the empirical snow fraction 999 ! so in turn on ice thickness 1000 END DO 1001 !zapond = sum(a_ip(1:m_index)) !js: from CICE 5.1.2; not in Icepack1.1.0-6-gac6195d 1002 1003 !------------------------------------------------------------------------ 1004 ! Drainage through brine network (permeability) 1005 !------------------------------------------------------------------------ 1006 !!! drainage due to ice permeability - Darcy's law 1007 1008 ! sea water level 1009 msno = 0._wp 1010 DO jl = 1 , jpl 1011 msno = msno + v_s(ji,jj,jl) * rhos 1012 END DO 1013 floe_weight = ( msno + rhoi*vt_i(ji,jj) + rau0*zvolp(ji,jj) ) / at_i(ji,jj) 1014 hsl_rel = floe_weight / rau0 & 1015 - ( ( sum(betan(:)*a_i(ji,jj,:)) / at_i(ji,jj) ) + alfan(1) ) 1016 1017 deltah = hpond - hsl_rel 1018 pressure_head = grav * rau0 * max(deltah, 0._wp) 1019 1020 ! drain if ice is permeable 1021 permflag = 0 1022 1023 IF (pressure_head > 0._wp) THEN 1024 DO jl = 1, jpl-1 1025 IF ( hicen(jl) /= 0._wp ) THEN 1026 1027 !IF (hicen(jl) > 0._wp) THEN !js: from CICE 5.1.2 1028 1029 perm = 0._wp ! MV ugly dummy patch 1030 CALL ice_thd_pnd_perm(t_i(ji,jj,:,jl), sz_i(ji,jj,:,jl), perm) ! bof 1031 IF (perm > 0._wp) permflag = 1 1032 1033 drain = perm*a_ip(ji,jj,jl)*pressure_head*rdt_ice / & 1034 (viscosity*hicen(jl)) 1035 zdvolp(ji,jj) = zdvolp(ji,jj) + min(drain, zvolp(ji,jj)) 1036 zvolp(ji,jj) = max(zvolp(ji,jj) - drain, 0._wp) 1037 1038 diag_dvpn_drn(ji,jj) = - drain ! diag (could be better coded) 1039 1040 IF (zvolp(ji,jj) < epsi10) THEN 1041 zdvolp(ji,jj) = zdvolp(ji,jj) + zvolp(ji,jj) 1042 zvolp(ji,jj) = 0._wp 1043 END IF 1044 END IF 1045 END DO 1046 1047 ! adjust melt pond dimensions 1048 IF (permflag > 0) THEN 1049 ! recompute pond depth 1050 CALL ice_thd_pnd_depth(reduced_aicen, asnon, hsnon, alfan, zvolp(ji,jj), cum_max_vol, hpond, m_index) 1051 DO jl = 1, m_index 1052 h_ip(ji,jj,jl) = hpond - alfan(jl) + alfan(1) 1053 a_ip(ji,jj,jl) = reduced_aicen(jl) 1054 END DO 1055 !zapond = sum(a_ip(1:m_index)) !js: from CICE 5.1.2; not in Icepack1.1.0-6-gac6195d 1056 END IF 1057 END IF ! pressure_head 1058 1059 !------------------------------- 1060 ! remove water from the snow 1061 !------------------------------- 1062 !------------------------------------------------------------------------ 1063 ! total melt pond volume in category does not include snow volume 1064 ! snow in melt ponds is not melted 1065 !------------------------------------------------------------------------ 1066 1067 ! MV here, it seems that we remove some meltwater from the ponds, but I can't really tell 1068 ! how much, so I did not diagnose it 1069 ! so if there is a problem here, nobody is going to see it... 1070 1071 1072 ! Calculate pond volume for lower categories 1073 DO jl = 1,m_index-1 1074 v_ip(ji,jj,jl) = a_ip(ji,jj,jl) * h_ip(ji,jj,jl) & ! what is not in the snow 1075 - (rhos/rhow) * asnon(jl) * min(hsnon(jl), h_ip(ji,jj,jl)) 1076 END DO 1077 1078 ! Calculate pond volume for highest category = remaining pond volume 1079 1080 ! The following is completely unclear to Martin at least 1081 ! Could we redefine properly and recode in a more readable way ? 1082 1083 ! m_index = last category with melt pond 1084 1085 IF (m_index == 1) v_ip(ji,jj,m_index) = zvolp(ji,jj) ! volume of mw in 1st category is the total volume of melt water 1086 1087 IF (m_index > 1) THEN 1088 IF (zvolp(ji,jj) > sum( v_ip(ji,jj,1:m_index-1))) THEN 1089 v_ip(ji,jj,m_index) = zvolp(ji,jj) - sum(v_ip(ji,jj,1:m_index-1)) 1090 ELSE 1091 v_ip(ji,jj,m_index) = 0._wp 1092 h_ip(ji,jj,m_index) = 0._wp 1093 a_ip(ji,jj,m_index) = 0._wp 1094 ! If remaining pond volume is negative reduce pond volume of 1095 ! lower category 1096 IF ( zvolp(ji,jj) + epsi10 < SUM(v_ip(ji,jj,1:m_index-1))) & 1097 v_ip(ji,jj,m_index-1) = v_ip(ji,jj,m_index-1) - sum(v_ip(ji,jj,1:m_index-1)) + zvolp(ji,jj) 1098 END IF 1099 END IF 1100 1101 DO jl = 1,m_index 1102 IF (a_ip(ji,jj,jl) > epsi10) THEN 1103 h_ip(ji,jj,jl) = v_ip(ji,jj,jl) / a_ip(ji,jj,jl) 1104 ELSE 1105 zdvolp(ji,jj) = zdvolp(ji,jj) + v_ip(ji,jj,jl) 1106 h_ip(ji,jj,jl) = 0._wp 1107 v_ip(ji,jj,jl) = 0._wp 1108 a_ip(ji,jj,jl) = 0._wp 1109 END IF 1110 END DO 1111 DO jl = m_index+1, jpl 1112 h_ip(ji,jj,jl) = 0._wp 1113 a_ip(ji,jj,jl) = 0._wp 1114 v_ip(ji,jj,jl) = 0._wp 1115 END DO 1116 1117 ENDIF 1118 END DO ! ji 1119 END DO ! jj 1120 1121 END SUBROUTINE ice_thd_pnd_area 1122 1123 1124 SUBROUTINE ice_thd_pnd_depth(aicen, asnon, hsnon, alfan, zvolp, cum_max_vol, hpond, m_index) 1125 !!------------------------------------------------------------------- 1126 !! *** ROUTINE ice_thd_pnd_depth *** 1127 !! 1128 !! ** Purpose : Compute melt pond depth 1129 !!------------------------------------------------------------------- 1130 1131 REAL (wp), DIMENSION(jpl), INTENT(IN) :: & 1132 aicen, & 1133 asnon, & 1134 hsnon, & 1135 alfan, & 1136 cum_max_vol 1137 1138 REAL (wp), INTENT(IN) :: & 1139 zvolp 1140 1141 REAL (wp), INTENT(OUT) :: & 1142 hpond 1143 1144 INTEGER, INTENT(OUT) :: & 1145 m_index 1146 1147 INTEGER :: n, ns 1148 1149 REAL (wp), DIMENSION(0:jpl+1) :: & 1150 hitl, & 1151 aicetl 1152 1153 REAL (wp) :: & 1154 rem_vol, & 1155 area, & 1156 vol, & 1157 tmp, & 1158 z0 = 0.0_wp 1159 1160 !---------------------------------------------------------------- 1161 ! hpond is zero if zvolp is zero - have we fully drained? 1162 !---------------------------------------------------------------- 1163 1164 IF (zvolp < epsi10) THEN 1165 hpond = z0 1166 m_index = 0 1167 ELSE 1168 1169 !---------------------------------------------------------------- 1170 ! Calculate the category where water fills up to 1171 !---------------------------------------------------------------- 1172 1173 !----------| 1174 ! | 1175 ! | 1176 ! |----------| -- -- 1177 !__________|__________|_________________________________________ ^ 1178 ! | | rem_vol ^ | Semi-filled 1179 ! | |----------|-- -- -- - ---|-- ---- -- -- --v layer 1180 ! | | | | 1181 ! | | | |hpond 1182 ! | | |----------| | |------- 1183 ! | | | | | | 1184 ! | | | |---v-----| 1185 ! | | m_index | | | 1186 !------------------------------------------------------------- 1187 1188 m_index = 0 ! 1:m_index categories have water in them 1189 DO n = 1, jpl 1190 IF (zvolp <= cum_max_vol(n)) THEN 1191 m_index = n 1192 IF (n == 1) THEN 1193 rem_vol = zvolp 1194 ELSE 1195 rem_vol = zvolp - cum_max_vol(n-1) 1196 END IF 1197 exit ! to break out of the loop 1198 END IF 1199 END DO 1200 m_index = min(jpl-1, m_index) 1201 1202 !---------------------------------------------------------------- 1203 ! semi-filled layer may have m_index different snow in it 1204 !---------------------------------------------------------------- 1205 1206 !----------------------------------------------------------- ^ 1207 ! | alfan(m_index+1) 1208 ! | 1209 !hitl(3)--> |----------| | 1210 !hitl(2)--> |------------| * * * * *| | 1211 !hitl(1)--> |----------|* * * * * * |* * * * * | | 1212 !hitl(0)-->------------------------------------------------- | ^ 1213 ! various snow from lower categories | |alfa(m_index) 1214 1215 ! hitl - heights of the snow layers from thinner and current categories 1216 ! aicetl - area of each snow depth in this layer 1217 1218 hitl(:) = z0 1219 aicetl(:) = z0 1220 DO n = 1, m_index 1221 hitl(n) = max(min(hsnon(n) + alfan(n) - alfan(m_index), & 1222 alfan(m_index+1) - alfan(m_index)), z0) 1223 aicetl(n) = asnon(n) 1224 1225 aicetl(0) = aicetl(0) + (aicen(n) - asnon(n)) 1226 END DO 1227 1228 hitl(m_index+1) = alfan(m_index+1) - alfan(m_index) 1229 aicetl(m_index+1) = z0 1230 1231 !---------------------------------------------------------------- 1232 ! reorder array according to hitl 1233 ! snow heights not necessarily in height order 1234 !---------------------------------------------------------------- 1235 1236 DO ns = 1, m_index+1 1237 DO n = 0, m_index - ns + 1 1238 IF (hitl(n) > hitl(n+1)) THEN ! swap order 1239 tmp = hitl(n) 1240 hitl(n) = hitl(n+1) 1241 hitl(n+1) = tmp 1242 tmp = aicetl(n) 1243 aicetl(n) = aicetl(n+1) 1244 aicetl(n+1) = tmp 1245 END IF 1246 END DO 1247 END DO 1248 1249 !---------------------------------------------------------------- 1250 ! divide semi-filled layer into set of sublayers each vertically homogenous 1251 !---------------------------------------------------------------- 1252 1253 !hitl(3)---------------------------------------------------------------- 1254 ! | * * * * * * * * 1255 ! |* * * * * * * * * 1256 !hitl(2)---------------------------------------------------------------- 1257 ! | * * * * * * * * | * * * * * * * * 1258 ! |* * * * * * * * * |* * * * * * * * * 1259 !hitl(1)---------------------------------------------------------------- 1260 ! | * * * * * * * * | * * * * * * * * | * * * * * * * * 1261 ! |* * * * * * * * * |* * * * * * * * * |* * * * * * * * * 1262 !hitl(0)---------------------------------------------------------------- 1263 ! aicetl(0) aicetl(1) aicetl(2) aicetl(3) 1264 1265 ! move up over layers incrementing volume 1266 DO n = 1, m_index+1 1267 1268 area = sum(aicetl(:)) - & ! total area of sub-layer 1269 (rhos/rau0) * sum(aicetl(n:jpl+1)) ! area of sub-layer occupied by snow 1270 1271 vol = (hitl(n) - hitl(n-1)) * area ! thickness of sub-layer times area 1272 1273 IF (vol >= rem_vol) THEN ! have reached the sub-layer with the depth within 1274 hpond = rem_vol / area + hitl(n-1) + alfan(m_index) - alfan(1) 1275 1276 exit 1277 ELSE ! still in sub-layer below the sub-layer with the depth 1278 rem_vol = rem_vol - vol 1279 END IF 1280 1281 END DO 1282 1283 END IF 1284 1285 END SUBROUTINE ice_thd_pnd_depth 1286 1287 1288 SUBROUTINE ice_thd_pnd_perm(ticen, salin, perm) 1289 !!------------------------------------------------------------------- 1290 !! *** ROUTINE ice_thd_pnd_perm *** 1291 !! 1292 !! ** Purpose : Determine the liquid fraction of brine in the ice 1293 !! and its permeability 1294 !!------------------------------------------------------------------- 1295 1296 REAL (wp), DIMENSION(nlay_i), INTENT(IN) :: & 1297 ticen, & ! internal ice temperature (K) 1298 salin ! salinity (ppt) !js: ppt according to cice 1299 1300 REAL (wp), INTENT(OUT) :: & 1301 perm ! permeability 1302 1303 REAL (wp) :: & 1304 Sbr ! brine salinity 1305 1306 REAL (wp), DIMENSION(nlay_i) :: & 1307 Tin, & ! ice temperature 1308 phi ! liquid fraction 1309 1310 INTEGER :: k 1311 1312 !----------------------------------------------------------------- 1313 ! Compute ice temperatures from enthalpies using quadratic formula 1314 !----------------------------------------------------------------- 1315 1316 DO k = 1,nlay_i 1317 Tin(k) = ticen(k) - rt0 !js: from K to degC 1318 END DO 1319 1320 !----------------------------------------------------------------- 1321 ! brine salinity and liquid fraction 1322 !----------------------------------------------------------------- 1323 1324 DO k = 1, nlay_i 1325 1326 ! Sbr = - Tin(k) / rTmlt ! Consistent expression with SI3 (linear liquidus) 1327 ! Best expression to date is that one 1328 Sbr = - 18.7 * Tin(k) - 0.519 * Tin(k)**2 - 0.00535 * Tin(k) **3 1329 phi(k) = salin(k) / Sbr 1330 1331 END DO 1332 1333 !----------------------------------------------------------------- 1334 ! permeability 1335 !----------------------------------------------------------------- 1336 1337 perm = 3.0e-08_wp * (minval(phi))**3 ! Golden et al. (2007) 1338 1339 END SUBROUTINE ice_thd_pnd_perm 1340 1341 1342 !---------------------------------------------------------------------------------------------------------------------- 304 1343 305 1344 SUBROUTINE ice_thd_pnd_init … … 318 1357 !! 319 1358 NAMELIST/namthd_pnd/ ln_pnd, ln_pnd_LEV , rn_apnd_min, rn_apnd_max, rn_pnd_flush, & 320 & ln_pnd_CST , rn_apnd, rn_hpnd, & 1359 & ln_pnd_CST , rn_apnd, rn_hpnd, & 1360 & ln_pnd_TOPO , & 321 1361 & ln_pnd_lids, ln_pnd_alb 322 1362 !!------------------------------------------------------------------- … … 336 1376 WRITE(numout,*) ' Namelist namicethd_pnd:' 337 1377 WRITE(numout,*) ' Melt ponds activated or not ln_pnd = ', ln_pnd 1378 WRITE(numout,*) ' Topographic melt pond scheme ln_pnd_TOPO = ', ln_pnd_TOPO 338 1379 WRITE(numout,*) ' Level ice melt pond scheme ln_pnd_LEV = ', ln_pnd_LEV 339 1380 WRITE(numout,*) ' Minimum ice fraction that contributes to melt ponds rn_apnd_min = ', rn_apnd_min … … 352 1393 IF( ln_pnd_CST ) THEN ; ioptio = ioptio + 1 ; nice_pnd = np_pndCST ; ENDIF 353 1394 IF( ln_pnd_LEV ) THEN ; ioptio = ioptio + 1 ; nice_pnd = np_pndLEV ; ENDIF 1395 IF( ln_pnd_TOPO ) THEN ; ioptio = ioptio + 1 ; nice_pnd = np_pndTOPO ; ENDIF 354 1396 IF( ioptio /= 1 ) & 355 & CALL ctl_stop( 'ice_thd_pnd_init: choose either none (ln_pnd=F) or only one pond scheme (ln_pnd_LEV or ln_pnd_CST)' )1397 & CALL ctl_stop( 'ice_thd_pnd_init: choose either none (ln_pnd=F) or only one pond scheme (ln_pnd_LEV, ln_pnd_CST or ln_pnd_TOPO)' ) 356 1398 ! 357 1399 SELECT CASE( nice_pnd ) -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/src/ICE/icevar.F90
r14075 r14302 565 565 END DO 566 566 END DO 567 568 !----------------------------------------------------------------- 569 ! zap small ponds 570 !----------------------------------------------------------------- 571 DO jj = 1, jpj 572 DO ji = 1, jpi 573 IF ( v_ip(ji,jj,jl) <= epsi10 ) THEN 574 a_ip(ji,jj,jl) = 0._wp 575 a_ip_frac(ji,jj,jl) = 0._wp 576 v_ip(ji,jj,jl) = 0._wp 577 h_ip(ji,jj,jl) = 0._wp 578 v_il(ji,jj,jl) = 0._wp 579 h_il(ji,jj,jl) = 0._wp 580 ENDIF 581 END DO 582 END DO 567 583 ! 568 584 END DO … … 696 712 WHERE( pe_i (1:npti,:,:) < 0._wp .AND. pe_i (1:npti,:,:) > -epsi06 ) pe_i (1:npti,:,:) = 0._wp ! e_i must be >= 0 697 713 WHERE( pe_s (1:npti,:,:) < 0._wp .AND. pe_s (1:npti,:,:) > -epsi06 ) pe_s (1:npti,:,:) = 0._wp ! e_s must be >= 0 698 IF( ln_pnd_LEV ) THEN714 IF( ln_pnd_LEV .OR. ln_pnd_TOPO ) THEN 699 715 WHERE( pa_ip(1:npti,:) < 0._wp .AND. pa_ip(1:npti,:) > -epsi10 ) pa_ip(1:npti,:) = 0._wp ! a_ip must be >= 0 700 716 WHERE( pv_ip(1:npti,:) < 0._wp .AND. pv_ip(1:npti,:) > -epsi10 ) pv_ip(1:npti,:) = 0._wp ! v_ip must be >= 0 -
NEMO/branches/UKMO/NEMO_4.0.4_GO8_package/src/ICE/icewri.F90
r14078 r14302 164 164 IF( iom_use('icebrv_cat' ) ) CALL iom_put( 'icebrv_cat' , bv_i * 100. * zmsk00l ) ! brine volume 165 165 IF( iom_use('iceapnd_cat' ) ) CALL iom_put( 'iceapnd_cat' , a_ip * zmsk00l ) ! melt pond frac for categories 166 IF( iom_use('icevpnd_cat' ) ) CALL iom_put( 'icevpnd_cat' , v_ip * zmsk00l ) ! melt pond volume for categories 166 167 IF( iom_use('icehpnd_cat' ) ) CALL iom_put( 'icehpnd_cat' , h_ip * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! melt pond thickness for categories 167 168 IF( iom_use('icehlid_cat' ) ) CALL iom_put( 'icehlid_cat' , h_il * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! melt pond lid thickness for categories 168 IF( iom_use('iceafpnd_cat') ) CALL iom_put( 'iceafpnd_cat', a_ip_frac * zmsk00l ) ! melt pond frac for categories169 IF( iom_use('iceafpnd_cat') ) CALL iom_put( 'iceafpnd_cat', a_ip_frac * zmsk00l ) ! melt pond frac per ice area for categories 169 170 IF( iom_use('iceaepnd_cat') ) CALL iom_put( 'iceaepnd_cat', a_ip_eff * zmsk00l ) ! melt pond effective frac for categories 170 171 IF( iom_use('icealb_cat' ) ) CALL iom_put( 'icealb_cat' , alb_ice * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! ice albedo for categories
Note: See TracChangeset
for help on using the changeset viewer.