Changeset 13080
- Timestamp:
- 2020-06-09T18:45:11+02:00 (5 years ago)
- Location:
- NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt
- Files:
-
- 25 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/cfgs/SHARED/field_def_nemo-ice.xml
r11575 r13080 158 158 <field id="hfxcndtop" long_name="Net conductive heat flux at the ice surface (neg = ice cooling)" standard_name="conductive_heat_flux_at_sea_ice_surface" unit="W/m2" /> 159 159 <field id="hfxcndbot" long_name="Net conductive heat flux at the ice bottom (neg = ice cooling)" standard_name="conductive_heat_flux_at_sea_ice_bottom" unit="W/m2" /> 160 <field id="hfxcndcpl" long_name="Conductive heat flux coming through the coupler" standard_name="conductive_heat_flux_from_coupler" unit="W/m2" /> 160 161 161 162 <!-- diags --> … … 245 246 <field id="iceconc_cat" long_name="Sea-ice concentration per category" unit="" /> 246 247 <field id="icethic_cat" long_name="Sea-ice thickness per category" unit="m" detect_missing_value="true" /> 248 <field id="icevol_cat" long_name="Sea-ice volume per category" unit="m" detect_missing_value="true" /> 247 249 <field id="snwthic_cat" long_name="Snow thickness per category" unit="m" detect_missing_value="true" /> 248 250 <field id="icesalt_cat" long_name="Sea-Ice Bulk salinity per category" unit="g/kg" detect_missing_value="true" /> … … 253 255 <field id="icehpnd_cat" long_name="Ice melt pond thickness per category" unit="m" detect_missing_value="true" /> 254 256 <field id="iceafpnd_cat" long_name="Ice melt pond fraction per category" unit="" /> 257 <field id="iceaepnd_cat" long_name="Ice melt pond effective fraction per category" unit="" /> 258 <field id="icelhpnd_cat" long_name="Ice melt pond lid thickness per category" unit="m" detect_missing_value="true" /> 255 259 <field id="icemask_cat" long_name="Fraction of time step with sea ice (per category)" unit="" /> 256 260 <field id="iceage_cat" long_name="Ice age per category" unit="days" detect_missing_value="true" /> -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/cfgs/SHARED/namelist_ice_ref
r11649 r13080 177 177 ln_pnd = .false. ! activate melt ponds or not 178 178 ln_pnd_H12 = .false. ! activate evolutive melt ponds (from Holland et al 2012) 179 rn_pnd_min = 0.15 ! Minimum ice fraction that contributes to melt ponds 180 rn_pnd_max = 0.7 ! Maximum ice fraction that contributes to melt ponds 181 ln_pnd_totfrac = .false. ! Use total ice fraction instead of category ice fraction in melt pond contribution fracton 182 ln_pnd_overflow = .false. ! Allow ponds to overflow and have vertical flushing 183 ln_pnd_lids = .false. ! Melt ponds can have frozen lids 184 ln_use_pndmass = .true. ! Use melt pond mass flux diagnostic, passing it to the ocean for emp 179 185 ln_pnd_CST = .false. ! activate constant melt ponds 180 186 rn_apnd = 0.2 ! prescribed pond fraction, at Tsu=0 degC … … 186 192 !------------------------------------------------------------------------------ 187 193 ln_iceini = .true. ! activate ice initialization (T) or not (F) 188 ln_iceini_file = .false. ! netcdf file provided for initialization (T) or not (F) 194 nn_iceini_file = 0 ! 0 = Initialise sea ice based on SSTs 195 ! 1 = Initialise sea ice from single category netcdf file 196 ! 2 = Initialise sea ice from multi category restart file 189 197 rn_thres_sst = 2.0 ! max temp. above Tfreeze with initial ice = (sst - tfreeze) 190 198 rn_hti_ini_n = 3.0 ! initial ice thickness (m), North -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/doc/namelists/namthd_pnd
r11536 r13080 4 4 ln_pnd = .false. ! activate melt ponds or not 5 5 ln_pnd_H12 = .false. ! activate evolutive melt ponds (from Holland et al 2012) 6 rn_pnd_min = 0.15 ! Minimum ice fraction that contributes to melt ponds 7 rn_pnd_max = 0.7 ! Maximum ice fraction that contributes to melt ponds 8 ln_pnd_totfrac = .false. ! Use total ice fraction instead of category ice fraction in melt pond contribution fracton 9 ln_pnd_overflow = .false. ! Allow ponds to overflow and have vertical flushing 10 ln_pnd_lids = .false. ! Melt ponds can have frozen lids 11 ln_use_pndmass = .true. ! Use melt pond mass flux diagnostic, passing it to the ocean for emp 6 12 ln_pnd_CST = .false. ! activate constant melt ponds 7 13 rn_apnd = 0.2 ! prescribed pond fraction, at Tsu=0 degC -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/ice.F90
r11715 r13080 70 70 !! a_ip | - | Ice pond concentration | | 71 71 !! v_ip | - | Ice pond volume per unit area| m | 72 !! lh_ip ! lh_ip_1d ! Ice pond lid thickness ! m ! 72 73 !! | 73 74 !!-------------|-------------|---------------------------------|-------| … … 195 196 REAL(wp), PUBLIC :: rn_hpnd !: prescribed pond depth (0<rn_hpnd<1) 196 197 LOGICAL , PUBLIC :: ln_pnd_alb !: melt ponds affect albedo 198 REAL(wp), PUBLIC :: rn_pnd_min !: Minimum ice fraction that contributes to melt ponds 199 REAL(wp), PUBLIC :: rn_pnd_max !: Maximum ice fraction that contributes to melt ponds 200 LOGICAL, PUBLIC :: ln_pnd_totfrac !: Use total ice fraction instead of category ice fraction 201 !: when determining ice fraction that contributes to melt ponds 202 LOGICAL, PUBLIC :: ln_pnd_overflow !: Allow ponds to overflow and have vertical flushing 203 LOGICAL, PUBLIC :: ln_pnd_lids !: Melt ponds can have frozen lids 204 LOGICAL, PUBLIC :: ln_use_pndmass !: Use melt pond mass flux diagnostic, passing it to the ocean for emp 197 205 198 206 ! !!** ice-diagnostics namelist (namdia) ** … … 295 303 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: h_i !: Ice thickness (m) 296 304 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: a_i !: Ice fractional areas (concentration) 305 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: a_i_last_couple !: Ice fractional area at last coupling time 297 306 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: v_i !: Ice volume per unit area (m) 298 307 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: v_s !: Snow volume per unit area (m) … … 331 340 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: v_ip !: melt pond volume per grid cell area [m] 332 341 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: a_ip_frac !: melt pond fraction (a_ip/a_i) 342 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: a_ip_eff !: melt pond effective fraction (not covered up by lid) (a_ip/a_i) 333 343 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: h_ip !: melt pond depth [m] 344 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: lh_ip !: melt pond lid thickness [m] 334 345 335 346 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: at_ip !: total melt pond concentration … … 435 446 436 447 ii = ii + 1 448 ALLOCATE( a_i_last_couple(jpi,jpj,jpl) , STAT=ierr(ii) ) 449 450 ii = ii + 1 437 451 ALLOCATE( u_ice(jpi,jpj) , v_ice(jpi,jpj) , & 438 452 & vt_i (jpi,jpj) , vt_s (jpi,jpj) , st_i(jpi,jpj) , at_i(jpi,jpj) , ato_i(jpi,jpj) , & … … 448 462 449 463 ii = ii + 1 450 ALLOCATE( a_ip(jpi,jpj,jpl) , v_ip(jpi,jpj,jpl) , a_ip_frac(jpi,jpj,jpl) , h_ip(jpi,jpj,jpl) , STAT = ierr(ii) ) 464 ALLOCATE( a_ip(jpi,jpj,jpl) , v_ip(jpi,jpj,jpl) , a_ip_frac(jpi,jpj,jpl) , h_ip(jpi,jpj,jpl), & 465 & lh_ip(jpi,jpj,jpl), a_ip_eff(jpi,jpj,jpl) , STAT = ierr(ii) ) 451 466 452 467 ii = ii + 1 -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/ice1d.F90
r11715 r13080 128 128 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: h_ip_1d !: 129 129 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: a_ip_frac_1d !: 130 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: a_ip_eff_1d !: 131 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:) :: lh_ip_1d !: Ice pond lid thickness [m] 130 132 131 133 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: t_s_1d !: corresponding to the 2D var t_s … … 157 159 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: a_ip_2d 158 160 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: v_ip_2d 161 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: lh_ip_2d 159 162 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: t_su_2d 160 163 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: h_i_2d … … 208 211 & dh_s_tot(jpij) , dh_i_sum(jpij) , dh_i_itm (jpij) , dh_i_bom(jpij) , dh_i_bog(jpij) , & 209 212 & dh_i_sub(jpij) , dh_s_mlt(jpij) , dh_snowice(jpij) , s_i_1d (jpij) , s_i_new (jpij) , & 210 & a_ip_1d (jpij) , v_ip_1d (jpij) , v_i_1d (jpij) , v_s_1d (jpij) , 211 & h_ip_1d (jpij) , a_ip_frac_1d(jpij) , 213 & a_ip_1d (jpij) , v_ip_1d (jpij) , v_i_1d (jpij) , v_s_1d (jpij) , lh_ip_1d(jpij) , & 214 & h_ip_1d (jpij) , a_ip_frac_1d(jpij) , a_ip_eff_1d(jpij) , & 212 215 & sv_i_1d (jpij) , oa_i_1d (jpij) , o_i_1d (jpij) , STAT=ierr(ii) ) 213 216 ! … … 226 229 ALLOCATE( a_i_2d (jpij,jpl) , a_ib_2d(jpij,jpl) , h_i_2d (jpij,jpl) , h_ib_2d(jpij,jpl) , & 227 230 & v_i_2d (jpij,jpl) , v_s_2d (jpij,jpl) , oa_i_2d(jpij,jpl) , sv_i_2d(jpij,jpl) , & 228 & a_ip_2d(jpij,jpl) , v_ip_2d(jpij,jpl) , t_su_2d(jpij,jpl) , 231 & a_ip_2d(jpij,jpl) , v_ip_2d(jpij,jpl) , t_su_2d(jpij,jpl) , lh_ip_2d(jpij,jpl), & 229 232 & STAT=ierr(ii) ) 230 233 -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icealb.F90
r11715 r13080 96 96 LOGICAL , INTENT(in ) :: ld_pnd_alb ! effect of melt ponds on albedo 97 97 REAL(wp), INTENT(in ), DIMENSION(:,:,:) :: pafrac_pnd ! melt pond relative fraction (per unit ice area) 98 ! This is the effective fraction not covered up by a pond lid 98 99 REAL(wp), INTENT(in ), DIMENSION(:,:,:) :: ph_pnd ! melt pond depth 99 100 REAL(wp), INTENT( out), DIMENSION(:,:,:) :: palb_cs ! albedo of ice under clear sky -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icedyn.F90
r11715 r13080 101 101 ELSEWHERE 102 102 h_ip(:,:,:) = 0._wp 103 lh_ip(:,:,:) = 0._wp 103 104 END WHERE 104 105 ! -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icedyn_adv.F90
r11715 r13080 84 84 ! !-----------------------! 85 85 CALL ice_dyn_adv_umx( nn_UMx, kt, u_ice, v_ice, h_i, h_s, h_ip, & 86 & ato_i, v_i, v_s, sv_i, oa_i, a_i, a_ip, v_ip, e_s, e_i )86 & ato_i, v_i, v_s, sv_i, oa_i, a_i, a_ip, v_ip, lh_ip, e_s, e_i ) 87 87 ! !-----------------------! 88 88 CASE( np_advPRA ) ! PRATHER scheme ! 89 89 ! !-----------------------! 90 90 CALL ice_dyn_adv_pra( kt, u_ice, v_ice, & 91 & ato_i, v_i, v_s, sv_i, oa_i, a_i, a_ip, v_ip, e_s, e_i )91 & ato_i, v_i, v_s, sv_i, oa_i, a_i, a_ip, v_ip, lh_ip, e_s, e_i ) 92 92 END SELECT 93 93 -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icedyn_adv_pra.F90
r11715 r13080 44 44 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: sxap , syap , sxxap , syyap , sxyap ! melt pond fraction 45 45 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: sxvp , syvp , sxxvp , syyvp , sxyvp ! melt pond volume 46 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: sxlhp, sylhp, sxxlhp, syylhp, sxylhp ! melt pond lid thickness 46 47 47 48 !! * Substitutions … … 55 56 56 57 SUBROUTINE ice_dyn_adv_pra( kt, pu_ice, pv_ice, & 57 & pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, p e_s, pe_i )58 & pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, plh_ip, pe_s, pe_i ) 58 59 !!---------------------------------------------------------------------- 59 60 !! ** routine ice_dyn_adv_pra ** … … 78 79 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: pa_ip ! melt pond fraction 79 80 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: pv_ip ! melt pond volume 81 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: plh_ip ! melt pond lid thickness 80 82 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pe_s ! snw heat content 81 83 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pe_i ! ice heat content … … 88 90 REAL(wp), DIMENSION(jpi,jpj,1) :: z0opw 89 91 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z0ice, z0snw, z0ai, z0smi, z0oi 90 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z0ap , z0vp 92 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z0ap , z0vp, z0lhp 91 93 REAL(wp), DIMENSION(jpi,jpj,nlay_s,jpl) :: z0es 92 94 REAL(wp), DIMENSION(jpi,jpj,nlay_i,jpl) :: z0ei … … 129 131 z0ap(:,:,jl) = pa_ip(:,:,jl) * e1e2t(:,:) ! Melt pond fraction 130 132 z0vp(:,:,jl) = pv_ip(:,:,jl) * e1e2t(:,:) ! Melt pond volume 133 z0lhp(:,:,jl) = plh_ip(:,:,jl) * e1e2t(:,:) ! Melt pond lid thickness 131 134 ENDIF 132 135 END DO … … 167 170 CALL adv_x( zdt , pu_ice , 1._wp , zarea , z0vp , sxvp , sxxvp , syvp , syyvp , sxyvp ) !--- melt pond volume 168 171 CALL adv_y( zdt , pv_ice , 0._wp , zarea , z0vp , sxvp , sxxvp , syvp , syyvp , sxyvp ) 172 CALL adv_x( zdt , pu_ice , 1._wp , zarea , z0lhp, sxlhp, sxxlhp, sylhp, syylhp, sxylhp) !--- melt pond lid thickness 173 CALL adv_y( zdt , pv_ice , 0._wp , zarea , z0lhp, sxlhp, sxxlhp, sylhp, syylhp, sxylhp) 169 174 ENDIF 170 175 END DO … … 202 207 CALL adv_y( zdt , pv_ice , 1._wp , zarea , z0vp , sxvp , sxxvp , syvp , syyvp , sxyvp ) !--- melt pond volume 203 208 CALL adv_x( zdt , pu_ice , 0._wp , zarea , z0vp , sxvp , sxxvp , syvp , syyvp , sxyvp ) 209 CALL adv_y( zdt , pv_ice , 1._wp , zarea , z0lhp, sxlhp, sxxlhp, sylhp, syylhp, sxylhp) !--- melt pond lid thickness 210 CALL adv_x( zdt , pu_ice , 0._wp , zarea , z0lhp, sxlhp, sxxlhp, sylhp, syylhp, sxylhp) 204 211 ENDIF 205 212 END DO … … 225 232 pa_ip(:,:,jl) = z0ap(:,:,jl) * r1_e1e2t(:,:) * tmask(:,:,1) 226 233 pv_ip(:,:,jl) = z0vp(:,:,jl) * r1_e1e2t(:,:) * tmask(:,:,1) 234 plh_ip(:,:,jl) = z0lhp (:,:,jl) * r1_e1e2t(:,:) * tmask(:,:,1) 227 235 ENDIF 228 236 END DO … … 231 239 ! Remove negative values (conservation is ensured) 232 240 ! (because advected fields are not perfectly bounded and tiny negative values can occur, e.g. -1.e-20) 233 CALL ice_var_zapneg( zdt, pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, p e_s, pe_i )241 CALL ice_var_zapneg( zdt, pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, plh_ip, pe_s, pe_i ) 234 242 ! 235 243 ! --- Ensure snow load is not too big --- ! … … 653 661 & sxap(jpi,jpj,jpl) , syap (jpi,jpj,jpl) , sxxap (jpi,jpj,jpl) , syyap (jpi,jpj,jpl) , sxyap (jpi,jpj,jpl) , & 654 662 & sxvp(jpi,jpj,jpl) , syvp (jpi,jpj,jpl) , sxxvp (jpi,jpj,jpl) , syyvp (jpi,jpj,jpl) , sxyvp (jpi,jpj,jpl) , & 663 & sxlhp(jpi,jpj,jpl) , sylhp (jpi,jpj,jpl), sxxlhp (jpi,jpj,jpl), syylhp (jpi,jpj,jpl), sxylhp (jpi,jpj,jpl), & 655 664 ! 656 665 & sxc0 (jpi,jpj,nlay_s,jpl) , syc0 (jpi,jpj,nlay_s,jpl) , sxxc0(jpi,jpj,nlay_s,jpl) , & … … 765 774 CALL iom_get( numrir, jpdom_autoglo, 'syyvp', syyvp ) 766 775 CALL iom_get( numrir, jpdom_autoglo, 'sxyvp', sxyvp ) 776 ! ! melt pond lid thickness 777 CALL iom_get( numrir, jpdom_autoglo, 'sxlhp' , sxlhp ) 778 CALL iom_get( numrir, jpdom_autoglo, 'sylhp' , sylhp ) 779 CALL iom_get( numrir, jpdom_autoglo, 'sxxlhp', sxxlhp ) 780 CALL iom_get( numrir, jpdom_autoglo, 'syylhp', syylhp ) 781 CALL iom_get( numrir, jpdom_autoglo, 'sxylhp', sxylhp ) 767 782 ENDIF 768 783 ! … … 782 797 sxap = 0._wp ; syap = 0._wp ; sxxap = 0._wp ; syyap = 0._wp ; sxyap = 0._wp ! melt pond fraction 783 798 sxvp = 0._wp ; syvp = 0._wp ; sxxvp = 0._wp ; syyvp = 0._wp ; sxyvp = 0._wp ! melt pond volume 799 sxlhp = 0._wp ; sylhp = 0._wp ; sxxlhp = 0._wp ; syylhp = 0._wp ; sxylhp = 0._wp ! melt pond lid thickness 784 800 ENDIF 785 801 ENDIF … … 862 878 CALL iom_rstput( iter, nitrst, numriw, 'syyvp', syyvp ) 863 879 CALL iom_rstput( iter, nitrst, numriw, 'sxyvp', sxyvp ) 880 ! ! melt pond lid thickness 881 CALL iom_rstput( iter, nitrst, numriw, 'sxlhp' , sxlhp ) 882 CALL iom_rstput( iter, nitrst, numriw, 'sylhp' , sylhp ) 883 CALL iom_rstput( iter, nitrst, numriw, 'sxxlhp', sxxlhp ) 884 CALL iom_rstput( iter, nitrst, numriw, 'syylhp', syylhp ) 885 CALL iom_rstput( iter, nitrst, numriw, 'sxylhp', sxylhp ) 864 886 ENDIF 865 887 ! -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icedyn_adv_umx.F90
r11715 r13080 60 60 61 61 SUBROUTINE ice_dyn_adv_umx( kn_umx, kt, pu_ice, pv_ice, ph_i, ph_s, ph_ip, & 62 & pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, p e_s, pe_i )62 & pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, plh_ip, pe_s, pe_i ) 63 63 !!---------------------------------------------------------------------- 64 64 !! *** ROUTINE ice_dyn_adv_umx *** … … 85 85 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: pa_ip ! melt pond concentration 86 86 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: pv_ip ! melt pond volume 87 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: plh_ip ! melt pond lid thickness 87 88 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pe_s ! snw heat content 88 89 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pe_i ! ice heat content … … 334 335 CALL adv_umx( zamsk, kn_umx, jt, kt, zdt, zudy , zvdx , zua_ho , zva_ho , zcu_box, zcv_box, & 335 336 & zhvar, pv_ip, zua_ups, zva_ups ) 337 ! lid thickness 338 zamsk = 0._wp 339 zhvar(:,:,:) = plh_ip(:,:,:) * z1_aip(:,:,:) 340 CALL adv_umx( zamsk, kn_umx, jt, kt, zdt, zudy , zvdx , zua_ho , zva_ho , zcu_box, zcv_box, & 341 & zhvar, plh_ip, zua_ups, zva_ups ) 342 336 343 ENDIF 337 344 ! … … 350 357 ! Remove negative values (conservation is ensured) 351 358 ! (because advected fields are not perfectly bounded and tiny negative values can occur, e.g. -1.e-20) 352 CALL ice_var_zapneg( zdt, pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, p e_s, pe_i )359 CALL ice_var_zapneg( zdt, pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, plh_ip, pe_s, pe_i ) 353 360 ! 354 361 ! Make sure ice thickness is not too big … … 1519 1526 !! 1520 1527 !! ** Purpose : Thickness correction in case advection scheme creates 1521 !! abnormally t ick ice or snow1528 !! abnormally thick ice or snow 1522 1529 !! 1523 1530 !! ** Method : 1- check whether ice thickness is larger than the surrounding 9-points -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icedyn_rdgrft.F90
r11715 r13080 503 503 REAL(wp) :: airdg1, oirdg1, aprdg1, virdg1, sirdg1 504 504 REAL(wp) :: airft1, oirft1, aprft1 505 REAL(wp), DIMENSION(jpij) :: airdg2, oirdg2, aprdg2, virdg2, sirdg2, vsrdg, vprdg ! area etc of new ridges506 REAL(wp), DIMENSION(jpij) :: airft2, oirft2, aprft2, virft , sirft , vsrft, vprft ! area etc of rafted ice505 REAL(wp), DIMENSION(jpij) :: airdg2, oirdg2, aprdg2, virdg2, sirdg2, vsrdg, vprdg, lhprdg ! area etc of new ridges 506 REAL(wp), DIMENSION(jpij) :: airft2, oirft2, aprft2, virft , sirft , vsrft, vprft, lhprft ! area etc of rafted ice 507 507 ! 508 508 REAL(wp), DIMENSION(jpij) :: ersw ! enth of water trapped into ridges … … 578 578 aprdg2(ji) = a_ip_2d(ji,jl1) * afrdg * hi_hrdg(ji,jl1) 579 579 vprdg (ji) = v_ip_2d(ji,jl1) * afrdg 580 lhprdg(ji) = lh_ip_2d(ji,jl1) * afrdg 580 581 aprft1 = a_ip_2d(ji,jl1) * afrft 581 582 aprft2(ji) = a_ip_2d(ji,jl1) * afrft * hi_hrft 582 583 vprft (ji) = v_ip_2d(ji,jl1) * afrft 584 lhprft(ji) = lh_ip_2d(ji,jl1) * afrft 583 585 ENDIF 584 586 … … 610 612 a_ip_2d(ji,jl1) = a_ip_2d(ji,jl1) - aprdg1 - aprft1 611 613 v_ip_2d(ji,jl1) = v_ip_2d(ji,jl1) - vprdg(ji) - vprft(ji) 614 lh_ip_2d(ji,jl1) = lh_ip_2d(ji,jl1) - lhprdg(ji) - lhprft(ji) 612 615 ENDIF 613 616 ENDIF … … 706 709 a_ip_2d (ji,jl2) = a_ip_2d(ji,jl2) + ( aprdg2(ji) * rn_fpndrdg * farea & 707 710 & + aprft2(ji) * rn_fpndrft * zswitch(ji) ) 711 lh_ip_2d (ji,jl2) = lh_ip_2d(ji,jl2) + ( lhprdg (ji) * rn_fpndrdg * fvol (ji) & 712 & + lhprft (ji) * rn_fpndrft * zswitch(ji) ) 708 713 ENDIF 709 714 … … 736 741 !---------------- 737 742 ! In case ridging/rafting lead to very small negative values (sometimes it happens) 738 CALL ice_var_roundoff( a_i_2d, v_i_2d, v_s_2d, sv_i_2d, oa_i_2d, a_ip_2d, v_ip_2d, ze_s_2d, ze_i_2d )743 CALL ice_var_roundoff( a_i_2d, v_i_2d, v_s_2d, sv_i_2d, oa_i_2d, a_ip_2d, v_ip_2d, lh_ip_2d, ze_s_2d, ze_i_2d ) 739 744 ! 740 745 END SUBROUTINE rdgrft_shift … … 848 853 CALL tab_3d_2d( npti, nptidx(1:npti), a_ip_2d(1:npti,1:jpl), a_ip(:,:,:) ) 849 854 CALL tab_3d_2d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip(:,:,:) ) 855 CALL tab_3d_2d( npti, nptidx(1:npti), lh_ip_2d(1:npti,1:jpl), lh_ip(:,:,:) ) 850 856 DO jl = 1, jpl 851 857 DO jk = 1, nlay_s … … 874 880 CALL tab_2d_3d( npti, nptidx(1:npti), a_ip_2d(1:npti,1:jpl), a_ip(:,:,:) ) 875 881 CALL tab_2d_3d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip(:,:,:) ) 882 CALL tab_2d_3d( npti, nptidx(1:npti), lh_ip_2d(1:npti,1:jpl), lh_ip(:,:,:) ) 876 883 DO jl = 1, jpl 877 884 DO jk = 1, nlay_s -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/iceistate.F90
r11715 r13080 41 41 ! !! ** namelist (namini) ** 42 42 LOGICAL, PUBLIC :: ln_iceini !: Ice initialization or not 43 LOGICAL, PUBLIC :: ln_iceini_file !: Ice initialization from 2D netcdf file 43 INTEGER, PUBLIC :: nn_iceini_file ! Ice initialization: 44 ! 0 = Initialise sea ice based on SSTs 45 ! 1 = Initialise sea ice from single category netcdf file 46 ! 2 = Initialise sea ice from multi category restart file 44 47 REAL(wp) :: rn_thres_sst 45 48 REAL(wp) :: rn_hti_ini_n, rn_hts_ini_n, rn_ati_ini_n, rn_smi_ini_n, rn_tmi_ini_n, rn_tsu_ini_n, rn_tms_ini_n … … 48 51 REAL(wp) :: rn_apd_ini_s, rn_hpd_ini_s 49 52 ! 50 ! ! if ln_iceini_file = T53 ! ! if nn_iceini_file = 1 51 54 INTEGER , PARAMETER :: jpfldi = 9 ! maximum number of files to read 52 55 INTEGER , PARAMETER :: jp_hti = 1 ! index of ice thickness (m) … … 155 158 a_ip (:,:,:) = 0._wp 156 159 v_ip (:,:,:) = 0._wp 160 lh_ip (:,:,:) = 0._wp 157 161 a_ip_frac(:,:,:) = 0._wp 162 a_ip_eff (:,:,:) = 0._wp 158 163 h_ip (:,:,:) = 0._wp 159 164 ! … … 167 172 IF( ln_iceini ) THEN 168 173 ! !---------------! 169 IF( ln_iceini_file)THEN ! Read a file !174 IF( nn_iceini_file == 1 )THEN ! Read a file ! 170 175 ! !---------------! 171 176 WHERE( ff_t(:,:) >= 0._wp ) ; zswitch(:,:) = 1._wp … … 362 367 a_ip_frac(:,:,:) = 0._wp 363 368 END WHERE 369 a_ip_eff(:,:,:) = a_ip_frac(:,:,:) 364 370 v_ip(:,:,:) = h_ip(:,:,:) * a_ip(:,:,:) 365 371 … … 466 472 TYPE(FLD_N), DIMENSION(jpfldi) :: slf_i ! array of namelist informations on the fields to read 467 473 ! 468 NAMELIST/namini/ ln_iceini, ln_iceini_file, rn_thres_sst, &474 NAMELIST/namini/ ln_iceini, nn_iceini_file, rn_thres_sst, & 469 475 & rn_hti_ini_n, rn_hti_ini_s, rn_hts_ini_n, rn_hts_ini_s, & 470 476 & rn_ati_ini_n, rn_ati_ini_s, rn_smi_ini_n, rn_smi_ini_s, & … … 493 499 WRITE(numout,*) ' Namelist namini:' 494 500 WRITE(numout,*) ' ice initialization (T) or not (F) ln_iceini = ', ln_iceini 495 WRITE(numout,*) ' ice initialization from a netcdf file ln_iceini_file = ', ln_iceini_file501 WRITE(numout,*) ' ice initialization from a netcdf file nn_iceini_file = ', nn_iceini_file 496 502 WRITE(numout,*) ' max ocean temp. above Tfreeze with initial ice rn_thres_sst = ', rn_thres_sst 497 IF( ln_iceini .AND. .NOT.ln_iceini_file) THEN503 IF( ln_iceini .AND. nn_iceini_file == 0 ) THEN 498 504 WRITE(numout,*) ' initial snw thickness in the north-south rn_hts_ini = ', rn_hts_ini_n,rn_hts_ini_s 499 505 WRITE(numout,*) ' initial ice thickness in the north-south rn_hti_ini = ', rn_hti_ini_n,rn_hti_ini_s … … 508 514 ENDIF 509 515 ! 510 IF( ln_iceini_file) THEN ! Ice initialization using input file516 IF( nn_iceini_file == 1 ) THEN ! Ice initialization using input file 511 517 ! 512 518 ! set si structure -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/iceitd.F90
r11715 r13080 411 411 CALL tab_3d_2d( npti, nptidx(1:npti), a_ip_2d(1:npti,1:jpl), a_ip ) 412 412 CALL tab_3d_2d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip ) 413 CALL tab_3d_2d( npti, nptidx(1:npti), lh_ip_2d(1:npti,1:jpl), lh_ip ) 413 414 CALL tab_3d_2d( npti, nptidx(1:npti), t_su_2d(1:npti,1:jpl), t_su ) 414 415 DO jl = 1, jpl … … 483 484 v_ip_2d(ji,jl1) = v_ip_2d(ji,jl1) - ztrans 484 485 v_ip_2d(ji,jl2) = v_ip_2d(ji,jl2) + ztrans 486 ! 487 ztrans = lh_ip_2d(ji,jl1) * zworka(ji) ! Pond lid thickness 488 lh_ip_2d(ji,jl1) = lh_ip_2d(ji,jl1) - ztrans 489 lh_ip_2d(ji,jl2) = lh_ip_2d(ji,jl2) + ztrans 485 490 ENDIF 486 491 ! … … 527 532 ! clem: The transfer between one category to another can lead to very small negative values (-1.e-20) 528 533 ! because of truncation error ( i.e. 1. - 1. /= 0 ) 529 CALL ice_var_roundoff( a_i_2d, v_i_2d, v_s_2d, sv_i_2d, oa_i_2d, a_ip_2d, v_ip_2d, ze_s_2d, ze_i_2d )534 CALL ice_var_roundoff( a_i_2d, v_i_2d, v_s_2d, sv_i_2d, oa_i_2d, a_ip_2d, v_ip_2d, lh_ip_2d, ze_s_2d, ze_i_2d ) 530 535 531 536 ! at_i must be <= rn_amax … … 555 560 CALL tab_2d_3d( npti, nptidx(1:npti), a_ip_2d(1:npti,1:jpl), a_ip ) 556 561 CALL tab_2d_3d( npti, nptidx(1:npti), v_ip_2d(1:npti,1:jpl), v_ip ) 562 CALL tab_2d_3d( npti, nptidx(1:npti), lh_ip_2d(1:npti,1:jpl), lh_ip ) 557 563 CALL tab_2d_3d( npti, nptidx(1:npti), t_su_2d(1:npti,1:jpl), t_su ) 558 564 DO jl = 1, jpl -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icerst.F90
r11715 r13080 132 132 CALL iom_rstput( iter, nitrst, numriw, 'a_ip' , a_ip ) 133 133 CALL iom_rstput( iter, nitrst, numriw, 'v_ip' , v_ip ) 134 CALL iom_rstput( iter, nitrst, numriw, 'lh_ip', lh_ip ) 134 135 ! Snow enthalpy 135 136 DO jk = 1, nlay_s … … 171 172 INTEGER :: jk 172 173 LOGICAL :: llok 173 INTEGER :: id0, id1, id2, id3, id4 ! local integer174 INTEGER :: id0, id1, id2, id3, id4, id5 ! local integer 174 175 CHARACTER(len=25) :: znam 175 176 CHARACTER(len=2) :: zchar, zchar1 … … 250 251 v_ip(:,:,:) = 0._wp 251 252 ENDIF 253 a_ip_eff(:,:,:) = a_ip(:,:,:) 254 ! melt pond lids 255 id5 = iom_varid( numrir, 'lh_ip' , ldstop = .FALSE. ) 256 IF( id5 > 0 ) THEN 257 CALL iom_get( numrir, jpdom_autoglo, 'lh_ip', lh_ip) 258 ELSE 259 lh_ip(:,:,:) = 0._wp 260 ENDIF 252 261 ! fields needed for Met Office (Jules) coupling 253 262 IF( ln_cpl ) THEN … … 274 283 CALL ice_istate( nit000 ) 275 284 ! 276 IF( .NOT.ln_iceini .OR. .NOT.ln_iceini_file) &277 & CALL ctl_stop('STOP', 'ice_rst_read: you need ln_ice_ini=T and ln_iceini_file=T')285 IF( .NOT.ln_iceini .OR. nn_iceini_file == 0 ) & 286 & CALL ctl_stop('STOP', 'ice_rst_read: you need ln_ice_ini=T and nn_iceini_file=0') 278 287 ! 279 288 ENDIF -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icesbc.F90
r11715 r13080 132 132 133 133 ! --- cloud-sky and overcast-sky ice albedos --- ! 134 CALL ice_alb( t_su, h_i, h_s, ln_pnd_alb, a_ip_ frac, h_ip, zalb_cs, zalb_os )134 CALL ice_alb( t_su, h_i, h_s, ln_pnd_alb, a_ip_eff, h_ip, zalb_cs, zalb_os ) 135 135 136 136 ! albedo depends on cloud fraction because of non-linear spectral effects -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icestp.F90
r11715 r13080 252 252 ! 253 253 ! ! Initial sea-ice state 254 IF( .NOT. ln_rstart ) THEN ! start from rest: sea-ice deduced from sst 255 CALL ice_istate_init 256 CALL ice_istate( nit000 ) 257 ELSE ! start from a restart file 258 CALL ice_rst_read 254 255 CALL ice_istate_init 256 IF ( ln_rstart .OR. nn_iceini_file == 2 ) THEN 257 CALL ice_rst_read ! start from a restart file 258 ELSE 259 CALL ice_istate( nit000 ) ! start from rest: sea-ice deduced from sst 259 260 ENDIF 260 261 CALL ice_var_glo2eqv -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icethd.F90
r11715 r13080 355 355 CALL tab_2d_1d( npti, nptidx(1:npti), h_ip_1d (1:npti), h_ip (:,:,kl) ) 356 356 CALL tab_2d_1d( npti, nptidx(1:npti), a_ip_frac_1d(1:npti), a_ip_frac(:,:,kl) ) 357 CALL tab_2d_1d( npti, nptidx(1:npti), a_ip_eff_1d (1:npti), a_ip_eff (:,:,kl) ) 358 CALL tab_2d_1d( npti, nptidx(1:npti), lh_ip_1d (1:npti), lh_ip (:,:,kl) ) 357 359 ! 358 360 CALL tab_2d_1d( npti, nptidx(1:npti), qprec_ice_1d (1:npti), qprec_ice ) … … 461 463 CALL tab_1d_2d( npti, nptidx(1:npti), h_ip_1d (1:npti), h_ip (:,:,kl) ) 462 464 CALL tab_1d_2d( npti, nptidx(1:npti), a_ip_frac_1d(1:npti), a_ip_frac(:,:,kl) ) 465 CALL tab_1d_2d( npti, nptidx(1:npti), a_ip_eff_1d (1:npti), a_ip_eff (:,:,kl) ) 466 CALL tab_1d_2d( npti, nptidx(1:npti), lh_ip_1d (1:npti), lh_ip (:,:,kl) ) 463 467 ! 464 468 CALL tab_1d_2d( npti, nptidx(1:npti), wfx_snw_sni_1d(1:npti), wfx_snw_sni ) -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icethd_pnd.F90
r11715 r13080 38 38 INTEGER, PARAMETER :: np_pndH12 = 2 ! Evolutive pond scheme (Holland et al. 2012) 39 39 40 REAL(wp), PARAMETER :: viscosity_dyn = 1.79e-3_wp 41 40 42 !! * Substitutions 41 43 # include "vectopt_loop_substitute.h90" … … 89 91 IF( a_i_1d(ji) > 0._wp .AND. t_su_1d(ji) >= rt0 ) THEN 90 92 a_ip_frac_1d(ji) = rn_apnd 93 a_ip_eff_1d(ji) = rn_apnd 91 94 h_ip_1d(ji) = rn_hpnd 92 95 a_ip_1d(ji) = a_ip_frac_1d(ji) * a_i_1d(ji) 93 96 ELSE 94 97 a_ip_frac_1d(ji) = 0._wp 98 a_ip_eff_1d(ji) = 0._wp 95 99 h_ip_1d(ji) = 0._wp 96 100 a_ip_1d(ji) = 0._wp … … 129 133 REAL(wp), PARAMETER :: zpnd_aspect = 0.8_wp ! pond aspect ratio 130 134 REAL(wp), PARAMETER :: zTp = -2._wp ! reference temperature 131 ! 135 REAL(wp), PARAMETER :: pnd_lid_max = 0.015_wp ! pond lid thickness above which the ponds disappear from the albedo calculation 136 REAL(wp), PARAMETER :: pnd_lid_min = 0.005_wp ! pond lid thickness below which the full pond area is used in the albedo calculation 137 ! 138 REAL(wp) :: tot_mlt ! Total ice and snow surface melt (some goes into ponds, some into the ocean) 132 139 REAL(wp) :: zfr_mlt ! fraction of available meltwater retained for melt ponding 133 REAL(wp) :: zdv_mlt ! available meltwater for melt ponding 140 REAL(wp) :: zdv_mlt ! available meltwater for melt ponding (equivalent volume change) 134 141 REAL(wp) :: z1_Tp ! inverse reference temperature 135 142 REAL(wp) :: z1_rhow ! inverse freshwater density 136 143 REAL(wp) :: z1_zpnd_aspect ! inverse pond aspect ratio 144 REAL(wp) :: z1_rhoi ! inverse ice density 137 145 REAL(wp) :: zfac, zdum 138 ! 139 INTEGER :: ji ! loop indices 146 REAL(wp) :: t_grad ! Temperature deficit for refreezing 147 REAL(wp) :: omega_dt ! Time independent accumulated variables used for freezing 148 REAL(wp) :: lh_ip_end ! Lid thickness at end of timestep (temporary variable) 149 REAL(wp) :: zdh_frz ! Amount of melt pond that freezes (m) 150 REAL(wp) :: dh_ip_over ! Pond thickness change due to leaking 151 REAL(wp) :: dh_i_pndleak ! Grid box mean change in water depth due to leaking back to the ocean 152 REAL(wp) :: h_gravity_head ! Height above sea level of the top of the melt pond 153 REAL(wp) :: h_percolation ! Distance between the base of the melt pond and the base of the sea ice 154 REAL(wp) :: Sbr ! Brine salinity 155 REAL(wp), DIMENSION(nlay_i) :: phi ! liquid fraction 156 REAL(wp) :: perm ! Permeability of the sea ice 157 REAL(wp) :: za_ip ! Temporary pond fraction 158 REAL(wp) :: max_h_diff_ts ! Maximum meltpond depth change due to leaking or overflow (m per ts) 159 REAL(wp) :: lfrac_pnd ! The fraction of the meltpond exposed (not inder a frozen lid) 160 161 ! 162 INTEGER :: ji, jk ! loop indices 140 163 !!------------------------------------------------------------------- 141 164 z1_rhow = 1._wp / rhow 142 165 z1_zpnd_aspect = 1._wp / zpnd_aspect 143 166 z1_Tp = 1._wp / zTp 167 z1_rhoi = 1._wp / rhoi 168 169 ! Define time-independent field for use in refreezing 170 omega_dt = 2.0_wp * rcnd_i * rdt_ice / (rLfus * rhow) 144 171 145 172 DO ji = 1, npti 146 ! !--------------------------------! 147 IF( h_i_1d(ji) < rn_himin) THEN ! Case ice thickness < rn_himin ! 148 ! !--------------------------------! 149 !--- Remove ponds on thin ice 173 174 ! !----------------------------------------------------! 175 IF( h_i_1d(ji) < rn_himin .OR. a_i_1d(ji) < epsi10 ) THEN ! Case ice thickness < rn_himin or tiny ice fraction ! 176 ! !----------------------------------------------------! 177 !--- Remove ponds on thin ice or tiny ice fractions 150 178 a_ip_1d(ji) = 0._wp 151 179 a_ip_frac_1d(ji) = 0._wp 152 180 h_ip_1d(ji) = 0._wp 181 lh_ip_1d(ji) = 0._wp 153 182 ! !--------------------------------! 154 183 ELSE ! Case ice thickness >= rn_himin ! 155 184 ! !--------------------------------! 156 185 v_ip_1d(ji) = h_ip_1d(ji) * a_ip_1d(ji) ! record pond volume at previous time step 157 ! 158 ! available meltwater for melt ponding [m, >0] and fraction 159 zdv_mlt = -( dh_i_sum(ji)*rhoi + dh_s_mlt(ji)*rhos ) * z1_rhow * a_i_1d(ji) 160 zfr_mlt = zrmin + ( zrmax - zrmin ) * a_i_1d(ji) ! from CICE doc 161 !zfr_mlt = zrmin + zrmax * a_i_1d(ji) ! from Holland paper 186 187 ! To avoid divide by zero errors in some calculations we will use a temporary pond fraction variable that has a minimum value of epsi06 188 za_ip = a_ip_1d(ji) 189 IF ( za_ip < epsi06 ) za_ip = epsi06 190 ! 191 ! available meltwater for melt ponding [m, >0] and fraction 192 tot_mlt = -( dh_i_sum(ji)*rhoi + dh_s_mlt(ji)*rhos ) * z1_rhow 193 IF ( ln_pnd_totfrac ) THEN 194 zfr_mlt = rn_pnd_min + ( rn_pnd_max - rn_pnd_min ) * at_i_1d(ji) ! Use total ice fraction 195 ELSE 196 zfr_mlt = rn_pnd_min + ( rn_pnd_max - rn_pnd_min ) * a_i_1d(ji) ! Use category ice fraction 197 ENDIF 198 zdv_mlt = zfr_mlt * tot_mlt 162 199 ! 163 200 !--- Pond gowth ---! 164 ! v_ip should never be negative, otherwise code crashes 165 v_ip_1d(ji) = MAX( 0._wp, v_ip_1d(ji) + zfr_mlt * zdv_mlt ) 201 v_ip_1d(ji) = v_ip_1d(ji) + zdv_mlt 202 ! 203 !--- Lid shrinking. ---! 204 IF ( ln_pnd_lids ) lh_ip_1d(ji) = lh_ip_1d(ji) - zdv_mlt / za_ip 166 205 ! 167 206 ! melt pond mass flux (<0) 168 IF( zdv_mlt > 0._wp ) THEN169 zfac = z fr_mlt * zdv_mlt * rhow * r1_rdtice207 IF( ln_use_pndmass .AND. zdv_mlt > 0._wp ) THEN 208 zfac = zdv_mlt * rhow * r1_rdtice 170 209 wfx_pnd_1d(ji) = wfx_pnd_1d(ji) - zfac 171 210 ! … … 177 216 ! 178 217 !--- Pond contraction (due to refreezing) ---! 179 v_ip_1d(ji) = v_ip_1d(ji) * EXP( 0.01_wp * MAX( zTp+rt0 - t_su_1d(ji), 0._wp ) * z1_Tp ) 218 IF ( ln_pnd_lids ) THEN 219 220 ! Code to use if using melt pond lids 221 IF ( t_su_1d(ji) < (zTp+rt0) .AND. v_ip_1d(ji) > 0._wp ) THEN 222 t_grad = (zTp+rt0) - t_su_1d(ji) 223 224 ! The following equation is a rearranged form of: 225 ! lid_thickness_end - lid_thickness_start = rcnd_i * t_grad * rdt_ice / (0.5*(lid_thickness_end + lid_thickness_start) * rLfus * rhow) 226 ! where: lid_thickness_start = lh_ip_1d(ji) 227 ! lid_thickness_end = lh_ip_end 228 ! omega_dt is a bunch of terms in the equation that do not change 229 ! Note the use of rhow instead of rhoi as we are working with volumes and it is mathematically easier 230 ! if the water and ice specific volumes (for the lid and the pond) are the same (have the same density). 231 232 lh_ip_end = SQRT(omega_dt * t_grad + lh_ip_1d(ji)**2) 233 zdh_frz = lh_ip_end - lh_ip_1d(ji) 234 235 ! Pond shrinking 236 v_ip_1d(ji) = v_ip_1d(ji) - zdh_frz * a_ip_1d(ji) 237 238 ! Lid growing 239 IF ( ln_pnd_lids ) lh_ip_1d(ji) = lh_ip_1d(ji) + zdh_frz 240 ELSE 241 zdh_frz = 0._wp 242 END IF 243 244 ELSE 245 ! Code to use if not using melt pond lids 246 v_ip_1d(ji) = v_ip_1d(ji) * EXP( 0.01_wp * MAX( zTp+rt0 - t_su_1d(ji), 0._wp ) * z1_Tp ) 247 ENDIF 248 249 ! 250 ! Make sure pond volume or lid thickness has not gone negative 251 IF ( v_ip_1d(ji) < 0._wp ) v_ip_1d(ji) = 0._wp 252 IF ( lh_ip_1d(ji) < 0._wp ) lh_ip_1d(ji) = 0._wp 180 253 ! 181 254 ! Set new pond area and depth assuming linear relation between h_ip and a_ip_frac … … 184 257 a_ip_frac_1d(ji) = a_ip_1d(ji) / a_i_1d(ji) 185 258 h_ip_1d(ji) = zpnd_aspect * a_ip_frac_1d(ji) 259 260 !--- Pond overflow and vertical flushing ---! 261 IF ( ln_pnd_overflow ) THEN 262 263 ! If pond area exceeds zfr_mlt * a_i_1d(ji) then reduce the pond volume 264 IF ( a_ip_1d(ji) > zfr_mlt * a_i_1d(ji) ) THEN 265 h_ip_1d(ji) = zpnd_aspect * zfr_mlt 266 a_ip_frac_1d(ji) = h_ip_1d(ji) / zpnd_aspect 267 a_ip_1d(ji) = a_ip_frac_1d(ji) * a_i_1d(ji) 268 v_ip_1d(ji) = h_ip_1d(ji) * a_ip_1d(ji) 269 ENDIF 270 271 ! If pond depth exceeds half the ice thickness then reduce the pond volume 272 IF ( h_ip_1d(ji) > 0.5_wp * h_i_1d(ji) ) THEN 273 h_ip_1d(ji) = 0.5_wp * h_i_1d(ji) 274 a_ip_frac_1d(ji) = h_ip_1d(ji) / zpnd_aspect 275 a_ip_1d(ji) = a_ip_frac_1d(ji) * a_i_1d(ji) 276 v_ip_1d(ji) = h_ip_1d(ji) * a_ip_1d(ji) 277 ENDIF 278 279 !-- Vertical flushing of pond water --! 280 ! The height above sea level of the top of the melt pond is the ratios of density of ice and water times the ice depth. 281 ! This assumes the pond is sitting on top of the ice. 282 h_gravity_head = h_i_1d(ji) * (rhow - rhoi) * z1_rhow 283 284 ! The depth through which water percolates is the distance under the melt pond 285 h_percolation = h_i_1d(ji) - h_ip_1d(ji) 286 287 ! Calculate the permeability of the ice (Assur 1958) 288 DO jk = 1, nlay_i 289 Sbr = - 1.2_wp & 290 - 21.8_wp * (t_i_1d(ji,jk)-rt0) & 291 - 0.919_wp * (t_i_1d(ji,jk)-rt0)**2 & 292 - 0.01878_wp * (t_i_1d(ji,jk)-rt0)**3 293 phi(jk) = sz_i_1d(ji,jk)/Sbr 294 END DO 295 perm = 3.0e-08_wp * (minval(phi))**3 296 297 ! Do the drainage using Darcy's law 298 IF ( perm > 0._wp ) THEN 299 dh_ip_over = -1.0_wp*perm*rhow*grav*h_gravity_head*rdt_ice / (viscosity_dyn*h_percolation) ! This should be a negative number 300 dh_ip_over = MIN(dh_ip_over, 0._wp) ! Make sure it is negative 301 302 h_ip_1d(ji) = MAX(0._wp, h_ip_1d(ji) + dh_ip_over) 303 a_ip_frac_1d(ji) = h_ip_1d(ji) / zpnd_aspect 304 a_ip_1d(ji) = a_ip_frac_1d(ji) * a_i_1d(ji) 305 v_ip_1d(ji) = h_ip_1d(ji) * a_ip_1d(ji) 306 ENDIF 307 ENDIF 308 309 ! If lid thickness is ten times greater than pond thickness then remove pond 310 IF ( ln_pnd_lids ) THEN 311 IF ( lh_ip_1d(ji) > h_ip_1d(ji) * 10._wp ) THEN 312 a_ip_1d(ji) = 0._wp 313 a_ip_frac_1d(ji) = 0._wp 314 h_ip_1d(ji) = 0._wp 315 lh_ip_1d(ji) = 0._wp 316 v_ip_1d(ji) = 0._wp 317 ENDIF 318 ENDIF 319 320 ! If any of the previous changes has removed all the ice thickness then remove ice area. 321 IF ( h_i_1d(ji) == 0._wp ) THEN 322 a_i_1d(ji) = 0._wp 323 h_s_1d(ji) = 0._wp 324 ENDIF 325 186 326 ! 187 327 ENDIF 328 329 ! Calculate how much meltpond is exposed to the atmosphere (not under a frozen lid) 330 IF ( lh_ip_1d(ji) < pnd_lid_min ) THEN ! Pond lid is very thin. Expose all the pond. 331 lfrac_pnd = 1.0 332 ELSE 333 IF ( lh_ip_1d(ji) > pnd_lid_max ) THEN ! Pond lid is very thick. Cover all the pond up with ice and nosw. 334 lfrac_pnd = 0.0 335 ELSE ! Pond lid is of intermediate thickness. Expose part of the melt pond. 336 lfrac_pnd = ( lh_ip_1d(ji) - pnd_lid_min ) / (pnd_lid_max - pnd_lid_min) 337 ENDIF 338 ENDIF 339 340 ! Calculate the melt pond effective area 341 a_ip_eff_1d(ji) = a_ip_frac_1d(ji) * lfrac_pnd 342 188 343 END DO 189 344 ! … … 205 360 INTEGER :: ios, ioptio ! Local integer 206 361 !! 207 NAMELIST/namthd_pnd/ ln_pnd, ln_pnd_H12, ln_pnd_CST, rn_apnd, rn_hpnd, ln_pnd_alb 362 NAMELIST/namthd_pnd/ ln_pnd, ln_pnd_H12, ln_pnd_CST, rn_apnd, rn_hpnd, ln_pnd_alb, & 363 rn_pnd_min, rn_pnd_max, ln_pnd_overflow, ln_pnd_lids, ln_pnd_totfrac, & 364 ln_use_pndmass 208 365 !!------------------------------------------------------------------- 209 366 ! … … 223 380 WRITE(numout,*) ' Melt ponds activated or not ln_pnd = ', ln_pnd 224 381 WRITE(numout,*) ' Evolutive melt pond fraction and depth (Holland et al 2012) ln_pnd_H12 = ', ln_pnd_H12 382 WRITE(numout,*) ' Minimum ice fraction that contributes to melt ponds rn_pnd_min = ', rn_pnd_min 383 WRITE(numout,*) ' Maximum ice fraction that contributes to melt ponds rn_pnd_max = ', rn_pnd_max 384 WRITE(numout,*) ' Use total ice fraction instead of category ice fraction ln_pnd_totfrac = ',ln_pnd_totfrac 385 WRITE(numout,*) ' Allow ponds to overflow and have vertical flushing ln_pnd_overflow = ',ln_pnd_overflow 386 WRITE(numout,*) ' Melt ponds can have frozen lids ln_pnd_lids = ',ln_pnd_lids 387 WRITE(numout,*) ' Use melt pond mass flux diagnostic, passing to ocean ln_use_pndmass = ',ln_use_pndmass 225 388 WRITE(numout,*) ' Prescribed melt pond fraction and depth ln_pnd_CST = ', ln_pnd_CST 226 389 WRITE(numout,*) ' Prescribed pond fraction rn_apnd = ', rn_apnd -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icethd_zdf_bl99.F90
r11715 r13080 31 31 !!---------------------------------------------------------------------- 32 32 !! NEMO/ICE 4.0 , NEMO Consortium (2018) 33 !! $Id $33 !! $Id: icethd_zdf_bl99.F90 10926 2019-05-03 12:32:10Z clem $ 34 34 !! Software governed by the CeCILL license (see ./LICENSE) 35 35 !!---------------------------------------------------------------------- … … 769 769 ! 770 770 ! --- calculate conduction fluxes (positive downward) 771 771 ! bottom ice conduction flux 772 772 DO ji = 1, npti 773 ! ! surface ice conduction flux 774 qcn_ice_top_1d(ji) = - isnow(ji) * zkappa_s(ji,0) * zg1s * ( t_s_1d(ji,1) - t_su_1d(ji) ) & 775 & - ( 1._wp - isnow(ji) ) * zkappa_i(ji,0) * zg1 * ( t_i_1d(ji,1) - t_su_1d(ji) ) 776 ! ! bottom ice conduction flux 777 qcn_ice_bot_1d(ji) = - zkappa_i(ji,nlay_i) * zg1 * ( t_bo_1d(ji ) - t_i_1d (ji,nlay_i) ) 773 qcn_ice_bot_1d(ji) = - zkappa_i(ji,nlay_i) * zg1 * ( t_bo_1d(ji ) - t_i_1d (ji,nlay_i) ) 778 774 END DO 779 775 ! surface ice conduction flux 776 IF( k_cnd == np_cnd_OFF .OR. k_cnd == np_cnd_EMU ) THEN 777 ! 778 DO ji = 1, npti 779 qcn_ice_top_1d(ji) = - isnow(ji) * zkappa_s(ji,0) * zg1s * ( t_s_1d(ji,1) - t_su_1d(ji) ) & 780 & - ( 1._wp - isnow(ji) ) * zkappa_i(ji,0) * zg1 * ( t_i_1d(ji,1) - t_su_1d(ji) ) 781 END DO 782 ! 783 ELSEIF( k_cnd == np_cnd_ON ) THEN 784 ! 785 DO ji = 1, npti 786 qcn_ice_top_1d(ji) = qcn_ice_1d(ji) 787 END DO 788 ! 789 ENDIF 790 ! surface ice temperature 791 IF( k_cnd == np_cnd_ON .AND. ln_cndemulate ) THEN 792 ! 793 DO ji = 1, npti 794 t_su_1d(ji) = ( qcn_ice_top_1d(ji) & ! calculate surface temperature 795 & + isnow(ji) * zkappa_s(ji,0) * zg1s * t_s_1d(ji,1) & 796 & + ( 1._wp - isnow(ji) ) * zkappa_i(ji,0) * zg1 * t_i_1d(ji,1) & 797 & ) / MAX( epsi10, isnow(ji) * zkappa_s(ji,0) * zg1s + ( 1._wp - isnow(ji) ) * zkappa_i(ji,0) * zg1 ) 798 t_su_1d(ji) = MAX( MIN( t_su_1d(ji), rt0 ), rt0 - 100._wp ) ! cap t_su 799 END DO 800 ! 801 ENDIF 780 802 ! 781 803 ! --- Diagnose the heat loss due to changing non-solar / conduction flux --- ! … … 785 807 DO ji = 1, npti 786 808 hfx_err_dif_1d(ji) = hfx_err_dif_1d(ji) - ( qns_ice_1d(ji) - zqns_ice_b(ji) ) * a_i_1d(ji) 787 END DO788 !789 ELSEIF( k_cnd == np_cnd_ON ) THEN790 !791 DO ji = 1, npti792 hfx_err_dif_1d(ji) = hfx_err_dif_1d(ji) - ( qcn_ice_top_1d(ji) - qcn_ice_1d(ji) ) * a_i_1d(ji)793 809 END DO 794 810 ! … … 856 872 t_i_1d (1:npti,:) = ztiold (1:npti,:) 857 873 qcn_ice_1d(1:npti) = qcn_ice_top_1d(1:npti) 858 859 !!clem860 ! remettre t_su_1d, qns_ice_1d et dqns_ice_1d comme avant puisqu'on devrait faire comme si on avant conduction = input861 !clem862 874 ENDIF 863 875 ! -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/iceupdate.F90
r11715 r13080 185 185 ! Snow/ice albedo (only if sent to coupler, useless in forced mode) 186 186 !------------------------------------------------------------------ 187 CALL ice_alb( t_su, h_i, h_s, ln_pnd_alb, a_ip_ frac, h_ip, zalb_cs, zalb_os ) ! cloud-sky and overcast-sky ice albedos187 CALL ice_alb( t_su, h_i, h_s, ln_pnd_alb, a_ip_eff, h_ip, zalb_cs, zalb_os ) ! cloud-sky and overcast-sky ice albedos 188 188 ! 189 189 alb_ice(:,:,:) = ( 1._wp - cldf_ice ) * zalb_cs(:,:,:) + cldf_ice * zalb_os(:,:,:) … … 279 279 IF( iom_use('hfxcndbot' ) ) CALL iom_put( 'hfxcndbot' , SUM( qcn_ice_bot * a_i_b, dim=3 ) ) ! Bottom conduction flux 280 280 IF( iom_use('hfxcndtop' ) ) CALL iom_put( 'hfxcndtop' , SUM( qcn_ice_top * a_i_b, dim=3 ) ) ! Surface conduction flux 281 IF( iom_use('hfxcndcpl' ) ) CALL iom_put( "hfxcndcpl" , SUM( qcn_ice * a_i_b, dim=3 ) ) ! Conduction flux we are giving it 281 282 282 283 ! controls -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icevar.F90
r11715 r13080 533 533 a_ip (ji,jj,jl) = a_ip (ji,jj,jl) * zswitch(ji,jj) 534 534 v_ip (ji,jj,jl) = v_ip (ji,jj,jl) * zswitch(ji,jj) 535 lh_ip (ji,jj,jl) = lh_ip (ji,jj,jl) * zswitch(ji,jj) 535 536 ! 536 537 END DO … … 555 556 556 557 557 SUBROUTINE ice_var_zapneg( pdt, pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, p e_s, pe_i )558 SUBROUTINE ice_var_zapneg( pdt, pato_i, pv_i, pv_s, psv_i, poa_i, pa_i, pa_ip, pv_ip, plh_ip, pe_s, pe_i ) 558 559 !!------------------------------------------------------------------- 559 560 !! *** ROUTINE ice_var_zapneg *** … … 570 571 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: pa_ip ! melt pond fraction 571 572 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: pv_ip ! melt pond volume 573 REAL(wp), DIMENSION(:,:,:) , INTENT(inout) :: plh_ip ! melt pond lid thickness 572 574 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pe_s ! snw heat content 573 575 REAL(wp), DIMENSION(:,:,:,:), INTENT(inout) :: pe_i ! ice heat content … … 637 639 WHERE( pv_ip (:,:,:) < 0._wp ) pv_ip (:,:,:) = 0._wp ! in theory one should change wfx_pnd(-) and wfx_sum(+) 638 640 ! but it does not change conservation, so keep it this way is ok 641 WHERE( plh_ip (:,:,:) < 0._wp ) plh_ip (:,:,:) = 0._wp 639 642 ! 640 643 END SUBROUTINE ice_var_zapneg 641 644 642 645 643 SUBROUTINE ice_var_roundoff( pa_i, pv_i, pv_s, psv_i, poa_i, pa_ip, pv_ip, p e_s, pe_i )646 SUBROUTINE ice_var_roundoff( pa_i, pv_i, pv_s, psv_i, poa_i, pa_ip, pv_ip, plh_ip, pe_s, pe_i ) 644 647 !!------------------------------------------------------------------- 645 648 !! *** ROUTINE ice_var_roundoff *** … … 654 657 REAL(wp), DIMENSION(:,:) , INTENT(inout) :: pa_ip ! melt pond fraction 655 658 REAL(wp), DIMENSION(:,:) , INTENT(inout) :: pv_ip ! melt pond volume 659 REAL(wp), DIMENSION(:,:) , INTENT(inout) :: plh_ip ! melt pond lid thickness 656 660 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: pe_s ! snw heat content 657 661 REAL(wp), DIMENSION(:,:,:), INTENT(inout) :: pe_i ! ice heat content … … 668 672 WHERE( pa_ip(1:npti,:) < 0._wp .AND. pa_ip(1:npti,:) > -epsi10 ) pa_ip(1:npti,:) = 0._wp ! a_ip must be >= 0 669 673 WHERE( pv_ip(1:npti,:) < 0._wp .AND. pv_ip(1:npti,:) > -epsi10 ) pv_ip(1:npti,:) = 0._wp ! v_ip must be >= 0 674 WHERE( plh_ip(1:npti,:) < 0._wp .AND. plh_ip(1:npti,:) > -epsi10 ) plh_ip(1:npti,:) = 0._wp ! lh_ip must be >= 0 670 675 ENDIF 671 676 ! -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/ICE/icewri.F90
r11715 r13080 149 149 150 150 ! --- category-dependent fields --- ! 151 IF( iom_use('icelhpnd_cat') ) CALL iom_put( 'icelhpnd_cat', lh_ip * zmsk00l ) ! melt pond lid thickness for categories 151 152 IF( iom_use('icemask_cat' ) ) CALL iom_put( 'icemask_cat' , zmsk00l ) ! ice mask 0% 152 153 IF( iom_use('iceconc_cat' ) ) CALL iom_put( 'iceconc_cat' , a_i * zmsk00l ) ! area for categories … … 165 166 IF( iom_use('iceafpnd_cat') ) CALL iom_put( 'iceafpnd_cat', a_ip_frac * zmsk00l ) ! melt pond frac for categories 166 167 IF( iom_use('icealb_cat' ) ) CALL iom_put( 'icealb_cat' , alb_ice * zmsk00l + zmiss_val * ( 1._wp - zmsk00l ) ) ! ice albedo for categories 168 169 IF( iom_use('icevol_cat' ) ) CALL iom_put( "icevol_cat" , v_i * zmsk00l ) ! volume for categories 170 IF( iom_use('iceaepnd_cat') ) CALL iom_put( 'iceaepnd_cat', a_ip_eff * zmsk00l ) ! melt pond effective frac for categories 167 171 168 172 !------------------ -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/OCE/BDY/bdyice.F90
r11715 r13080 98 98 & , v_i , 'T', 1., v_s , 'T', 1., sv_i, 'T', 1. & 99 99 & , kfillmode=jpfillnothing ,lsend=llsend1, lrecv=llrecv1 ) 100 CALL lbc_lnk_multi( 'bdyice', lh_ip,'T', 1., kfillmode=jpfillnothing ,lsend=llsend1, lrecv=llrecv1 ) 100 101 ! exchange 4d arrays : third dimension = 1 and then third dimension = jpk 101 102 CALL lbc_lnk_multi( 'bdyice', t_s , 'T', 1., e_s , 'T', 1., kfillmode=jpfillnothing ,lsend=llsend1, lrecv=llrecv1 ) … … 163 164 a_ip(ji,jj, jl) = ( a_ip(ji,jj, jl) * zwgt1 + dta%aip(i_bdy,jl) * zwgt ) * tmask(ji,jj,1) ! Ice pond concentration 164 165 h_ip(ji,jj, jl) = ( h_ip(ji,jj, jl) * zwgt1 + dta%hip(i_bdy,jl) * zwgt ) * tmask(ji,jj,1) ! Ice pond depth 166 lh_ip(ji,jj, jl) = ( lh_ip(ji,jj, jl) * zwgt1 + dta%hip(i_bdy,jl) * zwgt ) * tmask(ji,jj,1) ! Ice pond lid thickness 165 167 ! 166 168 sz_i(ji,jj,:,jl) = s_i(ji,jj,jl) … … 170 172 a_ip(ji,jj,jl) = 0._wp 171 173 h_ip(ji,jj,jl) = 0._wp 174 lh_ip(ji,jj,jl) = 0._wp 172 175 ENDIF 173 176 ! … … 231 234 a_ip(ji,jj, jl) = a_ip(ib,jb, jl) 232 235 h_ip(ji,jj, jl) = h_ip(ib,jb, jl) 236 lh_ip(ji,jj, jl) = lh_ip(ib,jb, jl) 233 237 ! 234 238 sz_i(ji,jj,:,jl) = sz_i(ib,jb,:,jl) … … 280 284 a_ip(ji,jj, jl) = 0._wp 281 285 v_ip(ji,jj, jl) = 0._wp 286 lh_ip(ji,jj, jl) = 0._wp 282 287 t_su(ji,jj, jl) = rt0 283 288 t_s (ji,jj,:,jl) = rt0 -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/OCE/SBC/sbc_ice.F90
r11715 r13080 70 70 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: wndm_ice !: wind speed module at T-point [m/s] 71 71 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:) :: sstfrz !: wind speed module at T-point [m/s] 72 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: tsfc_ice !: sea ice surface skin temperature (on categories)73 72 #endif 74 73 … … 93 92 94 93 ! already defined in ice.F90 for SI3 95 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: a_i 94 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: a_i ! Sea ice fraction on categories 95 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: a_i_last_couple ! Sea ice fraction on categories at the last coupling point 96 96 REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: h_i, h_s 97 97 … … 132 132 & qemp_ice(jpi,jpj) , qevap_ice(jpi,jpj,jpl) , qemp_oce (jpi,jpj) , & 133 133 & qns_oce (jpi,jpj) , qsr_oce (jpi,jpj) , emp_oce (jpi,jpj) , & 134 & emp_ice (jpi,jpj) , tsfc_ice (jpi,jpj,jpl) , sstfrz(jpi,jpj) , STAT= ierr(2) )134 & emp_ice (jpi,jpj) , sstfrz (jpi,jpj) , STAT= ierr(2) ) 135 135 #endif 136 136 -
NEMO/branches/UKMO/NEMO_4.0.1_fix_cpl_retain_melt/src/OCE/SBC/sbccpl.F90
r11715 r13080 48 48 USE lib_mpp ! distribued memory computing library 49 49 USE lbclnk ! ocean lateral boundary conditions (or mpp link) 50 51 #if defined key_oasis3 52 USE mod_oasis, ONLY : OASIS_Sent, OASIS_ToRest, OASIS_SentOut, OASIS_ToRestOut 53 #endif 50 54 51 55 IMPLICIT NONE … … 152 156 INTEGER, PARAMETER :: jps_wlev = 32 ! water level 153 157 INTEGER, PARAMETER :: jps_fice1 = 33 ! first-order ice concentration (for semi-implicit coupling of atmos-ice fluxes) 154 INTEGER, PARAMETER :: jps_a_p = 34 ! meltpond area 158 INTEGER, PARAMETER :: jps_a_p = 34 ! meltpond area fraction 155 159 INTEGER, PARAMETER :: jps_ht_p = 35 ! meltpond thickness 156 160 INTEGER, PARAMETER :: jps_kice = 36 ! sea ice effective conductivity … … 159 163 160 164 INTEGER, PARAMETER :: jpsnd = 38 ! total number of fields sent 165 166 #if ! defined key_oasis3 167 ! Dummy variables to enable compilation when oasis3 is not being used 168 INTEGER :: OASIS_Sent = -1 169 INTEGER :: OASIS_SentOut = -1 170 INTEGER :: OASIS_ToRest = -1 171 INTEGER :: OASIS_ToRestOut = -1 172 #endif 161 173 162 174 ! !!** namelist namsbc_cpl ** … … 184 196 LOGICAL :: ln_usecplmask ! use a coupling mask file to merge data received from several models 185 197 ! -> file cplmask.nc with the float variable called cplmask (jpi,jpj,nn_cplmodel) 198 LOGICAL :: ln_scale_ice_fluxes ! Scale sea ice fluxes by the sea ice fractions at the previous coupling point 186 199 TYPE :: DYNARR 187 200 REAL(wp), POINTER, DIMENSION(:,:,:) :: z3 … … 250 263 NAMELIST/namsbc_cpl/ sn_snd_temp , sn_snd_alb , sn_snd_thick, sn_snd_crt , sn_snd_co2 , & 251 264 & sn_snd_ttilyr, sn_snd_cond , sn_snd_mpnd , sn_snd_sstfrz, sn_snd_thick1, & 265 & ln_scale_ice_fluxes, & 252 266 & sn_snd_ifrac , sn_snd_crtw , sn_snd_wlev , sn_rcv_hsig , sn_rcv_phioc, & 253 267 & sn_rcv_w10m , sn_rcv_taumod, sn_rcv_tau , sn_rcv_dqnsdt, sn_rcv_qsr , & … … 279 293 ENDIF 280 294 IF( lwp .AND. ln_cpl ) THEN ! control print 295 WRITE(numout,*)' ln_scale_ice_fluxes = ', ln_scale_ice_fluxes 281 296 WRITE(numout,*)' received fields (mutiple ice categogies)' 282 297 WRITE(numout,*)' 10m wind module = ', TRIM(sn_rcv_w10m%cldes ), ' (', TRIM(sn_rcv_w10m%clcat ), ')' … … 815 830 END SELECT 816 831 832 ! Initialise ice fractions from last coupling time to zero 833 a_i_last_couple(:,:,:) = 0._wp 834 835 817 836 ! ! ------------------------- ! 818 837 ! ! Ice Meltponds ! … … 1243 1262 IF( srcv(jpr_co2)%laction ) atm_co2(:,:) = frcv(jpr_co2)%z3(:,:,1) 1244 1263 ! 1245 ! ! ================== !1246 ! ! ice skin temp. !1247 ! ! ================== !1248 #if defined key_si31249 ! needed by Met Office1250 IF( srcv(jpr_ts_ice)%laction ) THEN1251 WHERE ( frcv(jpr_ts_ice)%z3(:,:,:) > 0.0 ) ; tsfc_ice(:,:,:) = 0.01252 ELSEWHERE( frcv(jpr_ts_ice)%z3(:,:,:) < -60. ) ; tsfc_ice(:,:,:) = -60.1253 ELSEWHERE ; tsfc_ice(:,:,:) = frcv(jpr_ts_ice)%z3(:,:,:)1254 END WHERE1255 ENDIF1256 #endif1257 1264 ! ! ========================= ! 1258 1265 ! ! Mean Sea Level Pressure ! (taum) … … 1630 1637 !! sprecip solid precipitation over the ocean 1631 1638 !!---------------------------------------------------------------------- 1632 REAL(wp), INTENT(in) , DIMENSION(:,:) :: picefr ! ice fraction [0 to 1]1633 ! !! ! optional arguments, used only in 'mixed oce-ice' case1634 REAL(wp), INTENT(in) , DIMENSION(:,:,:), OPTIONAL :: palbi ! all skies ice albedo1635 REAL(wp), INTENT(in) , DIMENSION(:,: ), OPTIONAL :: psst ! sea surface temperature [Celsius]1636 REAL(wp), INTENT(in ), DIMENSION(:,:,:), OPTIONAL :: pist ! ice surface temperature [Kelvin]1637 REAL(wp), INTENT(in) , DIMENSION(:,:,:), OPTIONAL :: phs ! snow depth [m]1638 REAL(wp), INTENT(in) , DIMENSION(:,:,:), OPTIONAL :: phi ! ice thickness [m]1639 REAL(wp), INTENT(in) , DIMENSION(:,:) :: picefr ! ice fraction [0 to 1] 1640 ! !! ! optional arguments, used only in 'mixed oce-ice' case or for Met-Office coupling 1641 REAL(wp), INTENT(in) , DIMENSION(:,:,:), OPTIONAL :: palbi ! all skies ice albedo 1642 REAL(wp), INTENT(in) , DIMENSION(:,: ), OPTIONAL :: psst ! sea surface temperature [Celsius] 1643 REAL(wp), INTENT(inout), DIMENSION(:,:,:), OPTIONAL :: pist ! ice surface temperature [Kelvin] => inout for Met-Office 1644 REAL(wp), INTENT(in) , DIMENSION(:,:,:), OPTIONAL :: phs ! snow depth [m] 1645 REAL(wp), INTENT(in) , DIMENSION(:,:,:), OPTIONAL :: phi ! ice thickness [m] 1639 1646 ! 1640 1647 INTEGER :: ji, jj, jl ! dummy loop index … … 1643 1650 REAL(wp), DIMENSION(jpi,jpj) :: zemp_tot, zemp_ice, zemp_oce, ztprecip, zsprecip , zevap_oce, zdevap_ice 1644 1651 REAL(wp), DIMENSION(jpi,jpj) :: zqns_tot, zqns_oce, zqsr_tot, zqsr_oce, zqprec_ice, zqemp_oce, zqemp_ice 1645 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zqns_ice, zqsr_ice, zdqns_ice, zqevap_ice, zevap_ice !!gm , zfrqsr_tr_i 1652 REAL(wp), DIMENSION(jpi,jpj) :: zevap_ice_total 1653 REAL(wp), DIMENSION(jpi,jpj,jpl) :: zqns_ice, zqsr_ice, zdqns_ice, zqevap_ice, zevap_ice, zqtr_ice_top, ztsu 1646 1654 !!---------------------------------------------------------------------- 1647 1655 ! … … 1663 1671 ztprecip(:,:) = frcv(jpr_rain)%z3(:,:,1) + zsprecip(:,:) ! May need to ensure positive here 1664 1672 zemp_tot(:,:) = frcv(jpr_tevp)%z3(:,:,1) - ztprecip(:,:) 1665 zemp_ice(:,:) = ( frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_snow)%z3(:,:,1) ) * picefr(:,:)1666 1673 CASE( 'oce and ice' ) ! received fields: jpr_sbpr, jpr_semp, jpr_oemp, jpr_ievp 1667 1674 zemp_tot(:,:) = ziceld(:,:) * frcv(jpr_oemp)%z3(:,:,1) + picefr(:,:) * frcv(jpr_sbpr)%z3(:,:,1) 1668 1675 zemp_ice(:,:) = frcv(jpr_semp)%z3(:,:,1) * picefr(:,:) 1669 1676 zsprecip(:,:) = frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_semp)%z3(:,:,1) 1670 ztprecip(:,:) = frcv(jpr_semp)%z3(:,:,1) - frcv(jpr_sbpr)%z3(:,:,1) + zsprecip(:,:) 1677 ztprecip(:,:) = frcv(jpr_semp)%z3(:,:,1) - frcv(jpr_sbpr)%z3(:,:,1) + zsprecip(:,:) 1671 1678 CASE( 'none' ) ! Not available as for now: needs additional coding below when computing zevap_oce 1672 1679 ! ! since fields received are not defined with none option … … 1675 1682 1676 1683 #if defined key_si3 1684 1685 ! --- evaporation over ice (kg/m2/s) --- ! 1686 zevap_ice_total(:,:) = 0._wp 1687 IF (sn_rcv_emp%clcat == 'yes') THEN 1688 DO jl=1,jpl 1689 IF (ln_scale_ice_fluxes) THEN 1690 zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,jl) * a_i_last_couple(:,:,jl) 1691 ELSE 1692 zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,jl) 1693 ENDIF 1694 zevap_ice_total(:,:) = zevap_ice_total(:,:) + zevap_ice(:,:,jl) 1695 ENDDO 1696 ELSE 1697 zevap_ice(:,:,1) = frcv(jpr_ievp)%z3(:,:,1 ) 1698 zevap_ice_total(:,:) = zevap_ice(:,:,1) 1699 ENDIF 1700 1701 IF ( TRIM( sn_rcv_emp%cldes ) == 'conservative' ) THEN 1702 ! For conservative case zemp_ice has not been defined yet. Do it now. 1703 zemp_ice(:,:) = zevap_ice_total(:,:) - frcv(jpr_snow)%z3(:,:,1) * picefr(:,:) 1704 END IF 1705 1677 1706 ! zsnw = snow fraction over ice after wind blowing (=picefr if no blowing) 1678 1707 zsnw(:,:) = 0._wp ; CALL ice_thd_snwblow( ziceld, zsnw ) … … 1683 1712 1684 1713 ! --- evaporation over ocean (used later for qemp) --- ! 1685 zevap_oce(:,:) = frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) 1686 1687 ! --- evaporation over ice (kg/m2/s) --- ! 1688 DO jl=1,jpl 1689 IF (sn_rcv_emp%clcat == 'yes') THEN ; zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,jl) 1690 ELSE ; zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,1 ) ; ENDIF 1691 ENDDO 1714 zevap_oce(:,:) = frcv(jpr_tevp)%z3(:,:,1) - zevap_ice_total(:,:) 1715 1716 1717 1692 1718 1693 1719 ! since the sensitivity of evap to temperature (devap/dT) is not prescribed by the atmosphere, we set it to 0 … … 1727 1753 sprecip (:,:) = zsprecip (:,:) 1728 1754 tprecip (:,:) = ztprecip (:,:) 1729 evap_ice(:,:,:) = zevap_ice(:,:,:) 1755 IF (ln_scale_ice_fluxes) THEN 1756 ! Convert from grid box means to sea ice means 1757 WHERE( a_i(:,:,:) > 0.0_wp ) evap_ice(:,:,:) = zevap_ice(:,:,:) / a_i(:,:,:) 1758 WHERE( a_i(:,:,:) <= 0.0_wp ) evap_ice(:,:,:) = 0.0 1759 ELSE 1760 evap_ice(:,:,:) = zevap_ice(:,:,:) 1761 ENDIF 1730 1762 DO jl = 1, jpl 1731 1763 devap_ice(:,:,jl) = zdevap_ice(:,:) … … 1774 1806 IF( iom_use('snow_ao_cea') ) CALL iom_put( 'snow_ao_cea' , sprecip(:,:) * ( 1._wp - zsnw(:,:) ) ) ! Snow over ice-free ocean (cell average) 1775 1807 IF( iom_use('snow_ai_cea') ) CALL iom_put( 'snow_ai_cea' , sprecip(:,:) * zsnw(:,:) ) ! Snow over sea-ice (cell average) 1776 IF( iom_use('subl_ai_cea') ) CALL iom_put( 'subl_ai_cea' , frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) * tmask(:,:,1) ) ! Sublimation over sea-ice (cell average)1808 IF( iom_use('subl_ai_cea') ) CALL iom_put( 'subl_ai_cea' , zevap_ice_total(:,:) * tmask(:,:,1) ) ! Sublimation over sea-ice (cell average) 1777 1809 IF( iom_use('evap_ao_cea') ) CALL iom_put( 'evap_ao_cea' , ( frcv(jpr_tevp)%z3(:,:,1) & 1778 & - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) ) * tmask(:,:,1) ) ! ice-free oce evap (cell average)1810 & - zevap_ice_total(:,:) ) * tmask(:,:,1) ) ! ice-free oce evap (cell average) 1779 1811 ! note: runoff output is done in sbcrnf (which includes icebergs too) and iceshelf output is done in sbcisf 1780 1812 ! … … 1784 1816 CASE( 'oce only' ) ! the required field is directly provided 1785 1817 zqns_tot(:,:) = frcv(jpr_qnsoce)%z3(:,:,1) 1818 ! For Met Office sea ice non-solar fluxes are already delt with by JULES so setting to zero 1819 ! here so the only flux is the ocean only one. 1820 zqns_ice(:,:,:) = 0._wp 1786 1821 CASE( 'conservative' ) ! the required fields are directly provided 1787 1822 zqns_tot(:,:) = frcv(jpr_qnsmix)%z3(:,:,1) … … 1847 1882 & + zsprecip(:,:) * ( 1._wp - zsnw ) * ( zcptsnw (:,:) - rLfus ) ! solid precip over ocean + snow melting 1848 1883 zqemp_ice(:,:) = zsprecip(:,:) * zsnw * ( zcptsnw (:,:) - rLfus ) ! solid precip over ice (qevap_ice=0 since atm. does not take it into account) 1849 !! zqemp_ice(:,:) = - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) * zcptsnw (:,:) & ! ice evap1850 !! & + zsprecip(:,:) * zsnw * zqprec_ice(:,:) * r1_rhos ! solid precip over ice1851 1884 1852 1885 ! --- total non solar flux (including evap/precip) --- ! … … 1900 1933 IF ( srcv(jpr_icb)%laction ) CALL iom_put('hflx_icb_cea' , - frcv(jpr_icb)%z3(:,:,1) * rLfus ) ! latent heat from icebergs melting 1901 1934 IF ( iom_use('hflx_rain_cea') ) CALL iom_put('hflx_rain_cea' , ( tprecip(:,:) - sprecip(:,:) ) * zcptrain(:,:) ) ! heat flux from rain (cell average) 1902 IF ( iom_use('hflx_evap_cea') ) CALL iom_put('hflx_evap_cea' , ( frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1)&1903 & * picefr(:,:) )* zcptn(:,:) * tmask(:,:,1) ) ! heat flux from evap (cell average)1935 IF ( iom_use('hflx_evap_cea') ) CALL iom_put('hflx_evap_cea' , ( frcv(jpr_tevp)%z3(:,:,1) - zevap_ice_total(:,:) ) & 1936 * zcptn(:,:) * tmask(:,:,1) ) ! heat flux from evap (cell average) 1904 1937 IF ( iom_use('hflx_snow_cea') ) CALL iom_put('hflx_snow_cea' , sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) ) ! heat flux from snow (cell average) 1905 1938 IF ( iom_use('hflx_snow_ao_cea') ) CALL iom_put('hflx_snow_ao_cea', sprecip(:,:) * ( zcptsnw(:,:) - rLfus ) & … … 1914 1947 CASE( 'oce only' ) 1915 1948 zqsr_tot(:,: ) = MAX( 0._wp , frcv(jpr_qsroce)%z3(:,:,1) ) 1949 ! For Met Office sea ice solar fluxes are already delt with by JULES so setting to zero 1950 ! here so the only flux is the ocean only one. 1951 zqsr_ice(:,:,:) = 0._wp 1916 1952 CASE( 'conservative' ) 1917 1953 zqsr_tot(:,: ) = frcv(jpr_qsrmix)%z3(:,:,1) … … 1992 2028 ENDDO 1993 2029 ENDIF 2030 CASE( 'none' ) 2031 zdqns_ice(:,:,:) = 0._wp 1994 2032 END SELECT 1995 2033 … … 2007 2045 ! ! ========================= ! 2008 2046 CASE ('coupled') 2009 qml_ice(:,:,:) = frcv(jpr_topm)%z3(:,:,:) 2010 qcn_ice(:,:,:) = frcv(jpr_botm)%z3(:,:,:) 2047 IF (ln_scale_ice_fluxes) THEN 2048 WHERE( a_i(:,:,:) > 0.0_wp ) qml_ice(:,:,:) = frcv(jpr_topm)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:) 2049 WHERE( a_i(:,:,:) <= 0.0_wp ) qml_ice(:,:,:) = 0.0_wp 2050 WHERE( a_i(:,:,:) > 0.0_wp ) qcn_ice(:,:,:) = frcv(jpr_botm)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:) 2051 WHERE( a_i(:,:,:) <= 0.0_wp ) qcn_ice(:,:,:) = 0.0_wp 2052 ELSE 2053 qml_ice(:,:,:) = frcv(jpr_topm)%z3(:,:,:) 2054 qcn_ice(:,:,:) = frcv(jpr_botm)%z3(:,:,:) 2055 ENDIF 2011 2056 END SELECT 2012 2057 ! … … 2017 2062 ! 2018 2063 ! ! ===> used prescribed cloud fraction representative for polar oceans in summer (0.81) 2019 ztri = 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ! surface transmission parameter(Grenfell Maykut 77)2064 ztri = 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ! surface transmission when hi>10cm (Grenfell Maykut 77) 2020 2065 ! 2021 qtr_ice_top(:,:,:) = ztri * qsr_ice(:,:,:) 2022 WHERE( phs(:,:,:) >= 0.0_wp ) qtr_ice_top(:,:,:) = 0._wp ! snow fully opaque 2023 WHERE( phi(:,:,:) <= 0.1_wp ) qtr_ice_top(:,:,:) = qsr_ice(:,:,:) ! thin ice transmits all solar radiation 2066 WHERE ( phs(:,:,:) <= 0._wp .AND. phi(:,:,:) < 0.1_wp ) ! linear decrease from hi=0 to 10cm 2067 zqtr_ice_top(:,:,:) = qsr_ice(:,:,:) * ( ztri + ( 1._wp - ztri ) * ( 1._wp - phi(:,:,:) * 10._wp ) ) 2068 ELSEWHERE( phs(:,:,:) <= 0._wp .AND. phi(:,:,:) >= 0.1_wp ) ! constant (ztri) when hi>10cm 2069 zqtr_ice_top(:,:,:) = qsr_ice(:,:,:) * ztri 2070 ELSEWHERE ! zero when hs>0 2071 zqtr_ice_top(:,:,:) = 0._wp 2072 END WHERE 2024 2073 ! 2025 2074 ELSEIF( ln_cndflx .AND. .NOT.ln_cndemulate ) THEN !== conduction flux as surface forcing ==! … … 2027 2076 ! ! ===> here we must receive the qtr_ice_top array from the coupler 2028 2077 ! for now just assume zero (fully opaque ice) 2029 qtr_ice_top(:,:,:) = 0._wp 2078 zqtr_ice_top(:,:,:) = 0._wp 2079 ! 2080 ENDIF 2081 ! 2082 IF( ln_mixcpl ) THEN 2083 DO jl=1,jpl 2084 qtr_ice_top(:,:,jl) = qtr_ice_top(:,:,jl) * xcplmask(:,:,0) + zqtr_ice_top(:,:,jl) * zmsk(:,:) 2085 ENDDO 2086 ELSE 2087 qtr_ice_top(:,:,:) = zqtr_ice_top(:,:,:) 2088 ENDIF 2089 ! ! ================== ! 2090 ! ! ice skin temp. ! 2091 ! ! ================== ! 2092 ! needed by Met Office 2093 IF( srcv(jpr_ts_ice)%laction ) THEN 2094 WHERE ( frcv(jpr_ts_ice)%z3(:,:,:) > 0.0 ) ; ztsu(:,:,:) = 0.0 + rt0 2095 ELSEWHERE( frcv(jpr_ts_ice)%z3(:,:,:) < -60. ) ; ztsu(:,:,:) = -60. + rt0 2096 ELSEWHERE ; ztsu(:,:,:) = frcv(jpr_ts_ice)%z3(:,:,:) + rt0 2097 END WHERE 2098 ! 2099 IF( ln_mixcpl ) THEN 2100 DO jl=1,jpl 2101 pist(:,:,jl) = pist(:,:,jl) * xcplmask(:,:,0) + ztsu(:,:,jl) * zmsk(:,:) 2102 ENDDO 2103 ELSE 2104 pist(:,:,:) = ztsu(:,:,:) 2105 ENDIF 2030 2106 ! 2031 2107 ENDIF … … 2047 2123 INTEGER, INTENT(in) :: kt 2048 2124 ! 2125 ! 2049 2126 INTEGER :: ji, jj, jl ! dummy loop indices 2050 2127 INTEGER :: isec, info ! local integer 2051 2128 REAL(wp) :: zumax, zvmax 2052 2129 REAL(wp), DIMENSION(jpi,jpj) :: zfr_l, ztmp1, ztmp2, zotx1, zoty1, zotz1, zitx1, zity1, zitz1 2053 REAL(wp), DIMENSION(jpi,jpj,jpl) :: ztmp3, ztmp4 2130 REAL(wp), DIMENSION(jpi,jpj,jpl) :: ztmp3, ztmp4 2131 2054 2132 !!---------------------------------------------------------------------- 2055 2133 ! … … 2193 2271 ENDIF 2194 2272 2273 ! If this coupling was successful then save ice fraction for use between coupling points. 2274 ! This is needed for some calculations where the ice fraction at the last coupling point 2275 ! is needed. 2276 IF( info == OASIS_Sent .OR. info == OASIS_ToRest .OR. & 2277 & info == OASIS_SentOut .OR. info == OASIS_ToRestOut ) THEN 2278 IF ( sn_snd_thick%clcat == 'yes' ) THEN 2279 a_i_last_couple(:,:,1:jpl) = a_i(:,:,1:jpl) 2280 ENDIF 2281 ENDIF 2282 2195 2283 IF( ssnd(jps_fice1)%laction ) THEN 2196 2284 SELECT CASE( sn_snd_thick1%clcat ) … … 2250 2338 ! ! Ice melt ponds ! 2251 2339 ! ! ------------------------- ! 2252 ! needed by Met Office 2340 ! needed by Met Office - 1) fraction of ponded ice; 2) local/actual pond depth 2253 2341 IF( ssnd(jps_a_p)%laction .OR. ssnd(jps_ht_p)%laction ) THEN 2254 2342 SELECT CASE( sn_snd_mpnd%cldes) … … 2256 2344 SELECT CASE( sn_snd_mpnd%clcat ) 2257 2345 CASE( 'yes' ) 2258 ztmp3(:,:,1:jpl) = a_ip(:,:,1:jpl) 2259 ztmp4(:,:,1:jpl) = v_ip(:,:,1:jpl) 2346 2347 ztmp3(:,:,1:jpl) = a_ip_eff(:,:,1:jpl) 2348 ztmp4(:,:,1:jpl) = h_ip(:,:,1:jpl) 2349 2260 2350 CASE( 'no' ) 2261 2351 ztmp3(:,:,:) = 0.0 2262 2352 ztmp4(:,:,:) = 0.0 2263 2353 DO jl=1,jpl 2264 ztmp3(:,:,1) = ztmp3(:,:,1) + a_ip (:,:,jpl)2265 ztmp4(:,:,1) = ztmp4(:,:,1) + v_ip(:,:,jpl)2354 ztmp3(:,:,1) = ztmp3(:,:,1) + a_ip_frac(:,:,jpl) 2355 ztmp4(:,:,1) = ztmp4(:,:,1) + h_ip(:,:,jpl) 2266 2356 ENDDO 2267 2357 CASE default ; CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_mpnd%clcat' )
Note: See TracChangeset
for help on using the changeset viewer.