Changeset 12371
- Timestamp:
- 2020-02-12T13:13:26+01:00 (5 years ago)
- Location:
- NEMO/branches/UKMO/NEMO_4.0_add_pond_lids_prints/src/ICE
- Files:
-
- 5 edited
-
icethd_pnd.F90 (modified) (1 diff)
-
icethd_zdf_bl99.F90 (modified) (12 diffs)
-
iceupdate.F90 (modified) (1 diff)
-
icevar.F90 (modified) (3 diffs)
-
icewri.F90 (modified) (1 diff)
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/UKMO/NEMO_4.0_add_pond_lids_prints/src/ICE/icethd_pnd.F90
r12369 r12371 127 127 REAL(wp), PARAMETER :: zrmin = 0.15_wp ! minimum fraction of available meltwater retained for melt ponding 128 128 REAL(wp), PARAMETER :: zrmax = 0.70_wp ! maximum - - - - - 129 REAL(wp), PARAMETER :: zpnd_aspect = 0. 8_wp ! pond aspect ratio129 REAL(wp), PARAMETER :: zpnd_aspect = 0.174_wp ! pond aspect ratio 130 130 REAL(wp), PARAMETER :: zTp = -2._wp ! reference temperature 131 131 ! -
NEMO/branches/UKMO/NEMO_4.0_add_pond_lids_prints/src/ICE/icethd_zdf_bl99.F90
r12369 r12371 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 !!---------------------------------------------------------------------- … … 89 89 REAL(wp) :: zg1 = 2._wp ! 90 90 REAL(wp) :: zgamma = 18009._wp ! for specific heat 91 REAL(wp) :: zbeta = 0.1 17_wp ! for thermal conductivity (could be 0.13)91 REAL(wp) :: zbeta = 0.13_wp ! for thermal conductivity (could be 0.13) 92 92 REAL(wp) :: zraext_s = 10._wp ! extinction coefficient of radiation in the snow 93 93 REAL(wp) :: zkimin = 0.10_wp ! minimum ice thermal conductivity 94 94 REAL(wp) :: ztsu_err = 1.e-5_wp ! range around which t_su is considered at 0C 95 95 REAL(wp) :: zdti_bnd = 1.e-4_wp ! maximal authorized error on temperature 96 REAL(wp) :: zhs_min = 0. 01_wp! minimum snow thickness for conductivity calculation96 REAL(wp) :: zhs_min = 0.1_wp ! minimum snow thickness for conductivity calculation 97 97 REAL(wp) :: ztmelts ! ice melting temperature 98 98 REAL(wp) :: zdti_max ! current maximal error on temperature … … 101 101 REAL(wp) :: zfac ! dummy factor 102 102 ! 103 REAL(wp), DIMENSION(jpij) :: isnow ! switch for presence (1) or absence (0) of snow 103 REAL(wp), DIMENSION(jpij) :: isnow ! fraction of sea ice that is snow covered (for thermodynamic use only) 104 LOGICAL, DIMENSION(jpij) :: l_snow ! is snow above the threshold 104 105 REAL(wp), DIMENSION(jpij) :: ztsub ! surface temperature at previous iteration 105 106 REAL(wp), DIMENSION(jpij) :: zh_i, z1_h_i ! ice layer thickness … … 153 154 !------------------ 154 155 DO ji = 1, npti 155 isnow(ji) = 1._wp - MAX( 0._wp , SIGN(1._wp, - h_s_1d(ji) ) ) ! is there snow or not 156 157 ! If the snow thickness drops below zhs_min then reduce the snow fraction instead 158 IF( h_s_1d(ji) < zhs_min ) THEN 159 ! isnow(ji) = h_s_1d(ji) / zhs_min 160 isnow(ji) = 0.0_wp 161 zh_s(ji) = zhs_min * r1_nlay_s 162 l_snow(ji) = .false. 163 ELSE 164 isnow(ji) = 1.0_wp 165 zh_s(ji) = h_s_1d(ji) * r1_nlay_s 166 l_snow(ji) = .true. 167 END IF 168 156 169 ! layer thickness 157 170 zh_i(ji) = h_i_1d(ji) * r1_nlay_i 158 zh_s(ji) = h_s_1d(ji) * r1_nlay_s 171 159 172 IF( to_print(ji) == 10 ) THEN 160 173 write(numout,*)'icethd_zdf_bl99: v_i_1d(ji), a_i_1d(ji), h_i_1d(ji) = ',v_i_1d(ji), ' ', a_i_1d(ji), ' ', h_i_1d(ji) … … 166 179 END WHERE 167 180 ! 168 WHERE( zh_s(1:npti) > 0._wp ) zh_s(1:npti) = MAX( zhs_min * r1_nlay_s, zh_s(1:npti) )169 181 ! 170 182 WHERE( zh_s(1:npti) > 0._wp ) ; z1_h_s(1:npti) = 1._wp / zh_s(1:npti) … … 310 322 IF ( .NOT. l_T_converged(ji) ) & 311 323 zkappa_s(ji,jk) = zghe(ji) * rn_cnd_s * z1_h_s(ji) 324 IF ( .NOT. l_snow(ji) ) zkappa_s(ji,jk) = 0.0_wp 312 325 END DO 313 326 END DO … … 333 346 END DO 334 347 DO ji = 1, npti ! Snow-ice interface 335 IF ( .NOT. l_T_converged(ji) ) &336 zkappa_i(ji,0) = zkappa_s(ji,nlay_s) * isnow(ji) + zkappa_i(ji,0) * ( 1._wp - isnow(ji) )337 348 338 349 IF (to_print(ji) == 10) THEN … … 359 370 DO jk = 1, nlay_s 360 371 DO ji = 1, npti 361 zeta_s(ji,jk) = rdt_ice * r1_rhos * r1_rcpi * z1_h_s(ji) 372 IF (l_snow(ji)) then 373 zeta_s(ji,jk) = rdt_ice * r1_rhos * r1_rcpi * z1_h_s(ji) 374 ELSE 375 zeta_s(ji,jk) = 0.0_wp 376 END IF 362 377 END DO 363 378 END DO … … 707 722 DO ji = 1, npti 708 723 ! !---------------------! 709 IF( h_s_1d(ji) > 0._wp) THEN ! snow-covered cells !724 IF( l_snow(ji) ) THEN ! snow-covered cells ! 710 725 ! !---------------------! 711 726 ! snow interior terms (bottom equation has the same form as the others) … … 895 910 ! 896 911 ! --- calculate conduction fluxes (positive downward) 897 912 ! bottom ice conduction flux 898 913 DO ji = 1, npti 899 ! ! surface ice conduction flux900 qcn_ice_top_1d(ji) = - isnow(ji) * zkappa_s(ji,0) * zg1s * ( t_s_1d(ji,1) - t_su_1d(ji) ) & 901 & - ( 1._wp - isnow(ji) ) * zkappa_i(ji,0) * zg1 * ( t_i_1d(ji,1) - t_su_1d(ji) )902 ! ! bottom ice conduction flux903 qcn_ice_bot_1d(ji) = - zkappa_i(ji,nlay_i) * zg1 * ( t_bo_1d(ji ) - t_i_1d (ji,nlay_i) )914 qcn_ice_bot_1d(ji) = - zkappa_i(ji,nlay_i) * zg1 * ( t_bo_1d(ji ) - t_i_1d (ji,nlay_i) ) 915 916 IF (to_print(ji) == 1) THEN 917 write(numout,*) 'icethd_zdf_bl99: qcn_ice_bot_1d(ji), zkappa_i(ji,nlay_i), zg1, t_bo_1d(ji ), t_i_1d (ji,nlay_i) = ',qcn_ice_bot_1d(ji), zkappa_i(ji,nlay_i), zg1, t_bo_1d(ji ), t_i_1d (ji,nlay_i) 918 ENDIF 904 919 END DO 905 920 ! surface ice conduction flux 921 IF( k_cnd == np_cnd_OFF .OR. k_cnd == np_cnd_EMU ) THEN 922 ! 923 DO ji = 1, npti 924 qcn_ice_top_1d(ji) = - isnow(ji) * zkappa_s(ji,0) * zg1s * ( t_s_1d(ji,1) - t_su_1d(ji) ) & 925 & - ( 1._wp - isnow(ji) ) * zkappa_i(ji,0) * zg1 * ( t_i_1d(ji,1) - t_su_1d(ji) ) 926 END DO 927 ! 928 ELSEIF( k_cnd == np_cnd_ON ) THEN 929 ! 930 DO ji = 1, npti 931 qcn_ice_top_1d(ji) = qcn_ice_1d(ji) 932 ! 933 t_su_1d(ji) = ( qcn_ice_top_1d(ji) & ! calculate surface temperature 934 & + isnow(ji) * zkappa_s(ji,0) * zg1s * t_s_1d(ji,1) & 935 & + ( 1._wp - isnow(ji) ) * zkappa_i(ji,0) * zg1 * t_i_1d(ji,1) & 936 & ) / MAX( epsi10, isnow(ji) * zkappa_s(ji,0) * zg1s + ( 1._wp - isnow(ji) ) * zkappa_i(ji,0) * zg1 ) 937 t_su_1d(ji) = MAX( MIN( t_su_1d(ji), rt0 ), rt0 - 100._wp ) ! cap t_su 938 END DO 939 ! 940 ENDIF 906 941 ! 907 942 ! --- Diagnose the heat loss due to changing non-solar / conduction flux --- ! … … 911 946 DO ji = 1, npti 912 947 hfx_err_dif_1d(ji) = hfx_err_dif_1d(ji) - ( qns_ice_1d(ji) - zqns_ice_b(ji) ) * a_i_1d(ji) 913 END DO914 !915 ELSEIF( k_cnd == np_cnd_ON ) THEN916 !917 DO ji = 1, npti918 hfx_err_dif_1d(ji) = hfx_err_dif_1d(ji) - ( qcn_ice_top_1d(ji) - qcn_ice_1d(ji) ) * a_i_1d(ji)919 948 END DO 920 949 ! … … 985 1014 ! effective conductivity and 1st layer temperature (needed by Met Office) 986 1015 DO ji = 1, npti 987 IF( h_ s_1d(ji) > 0.1_wp ) THEN988 cnd_ice_1d(ji) = 2._wp * zkappa_ s(ji,0)1016 IF( h_i_1d(ji) > 0.1_wp ) THEN 1017 cnd_ice_1d(ji) = 2._wp * zkappa_i(ji,0) 989 1018 ELSE 990 IF( h_i_1d(ji) > 0.1_wp ) THEN 991 cnd_ice_1d(ji) = 2._wp * zkappa_i(ji,0) 992 ELSE 993 cnd_ice_1d(ji) = 2._wp * ztcond_i(ji,0) * 10._wp 994 ENDIF 1019 cnd_ice_1d(ji) = 2._wp * ztcond_i(ji,0) * 10._wp 995 1020 ENDIF 996 1021 t1_ice_1d(ji) = isnow(ji) * t_s_1d(ji,1) + ( 1._wp - isnow(ji) ) * t_i_1d(ji,1) -
NEMO/branches/UKMO/NEMO_4.0_add_pond_lids_prints/src/ICE/iceupdate.F90
r10888 r12371 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 ! diags -
NEMO/branches/UKMO/NEMO_4.0_add_pond_lids_prints/src/ICE/icevar.F90
r12369 r12371 317 317 REAL(wp), PARAMETER :: zsi0 = 3.5_wp 318 318 REAL(wp), PARAMETER :: zsi1 = 4.5_wp 319 320 REAL(wp), PARAMETER :: min_salin = 0.1_wp ! Minimum salinity in sea ice 321 REAL(wp) :: zn ! Fraction of ice depth 322 REAL(wp), PARAMETER :: msal = 0.573_wp 323 REAL(wp), PARAMETER :: nsal = 0.407_wp 324 REAL(wp), PARAMETER :: saltmax = 9.6_wp 325 ! REAL(wp), PARAMETER, DIMENSION(4) :: weights = (/0.06363419, 0.25545967, 0.33155539, 0.34935076/) 326 REAL(wp), PARAMETER, DIMENSION(4) :: weights = (/0.0703805, 0.25526203, 0.32860314, 0.34575433/) 327 319 328 !!------------------------------------------------------------------- 320 329 … … 350 359 DO jj = 1, jpj 351 360 DO ji = 1, jpi 352 zalpha(ji,jj,jl) = MAX( 0._wp , MIN( ( zsi1 - s_i(ji,jj,jl) ) * z1_dS , 1._wp ) )361 zalpha(ji,jj,jl) = 1._wp 353 362 ! ! force a constant profile when SSS too low (Baltic Sea) 354 363 IF( 2._wp * s_i(ji,jj,jl) >= sss_m(ji,jj) ) zalpha(ji,jj,jl) = 0._wp … … 363 372 DO jj = 1, jpj 364 373 DO ji = 1, jpi 365 ! ! linear profile with 0 surface value 366 zs0 = z_slope_s(ji,jj,jl) * ( REAL(jk,wp) - 0.5_wp ) * h_i(ji,jj,jl) * r1_nlay_i 374 375 ! Copied over from CICE 376 ! zn = (real(jk,kind=dbl_kind)-0.5_wp) / real(nlay_i,kind=dbl_kind) 377 ! sz_i(ji,jj,jk,jl)=(saltmax/2.0_wp)*(1.0_wp-cos(pi*zn**(nsal/(msal+zn)))) 378 ! sz_i(ji,jj,jk,jl) = max(sz_i(ji,jj,jk,jl), min_salin) 379 380 ! Weighting method to match CICE 381 zs0 = (s_i(ji,jj,jl)*4.0_wp) * weights(jk) 367 382 zs = zalpha(ji,jj,jl) * zs0 + ( 1._wp - zalpha(ji,jj,jl) ) * s_i(ji,jj,jl) ! weighting the profile 368 383 sz_i(ji,jj,jk,jl) = MIN( rn_simax, MAX( zs, rn_simin ) ) 384 385 386 ! zs0 = z_slope_s(ji,jj,jl) * ( REAL(jk,wp) - 0.5_wp ) * h_i(ji,jj,jl) * r1_nlay_i 387 ! zs = zalpha(ji,jj,jl) * zs0 + ( 1._wp - zalpha(ji,jj,jl) ) * s_i(ji,jj,jl) ! weighting the profile 388 ! sz_i(ji,jj,jk,jl) = MIN( rn_simax, MAX( zs, rn_simin ) ) 369 389 370 390 IF ( jl == 1 .AND. jj == 26 .AND. ji == 42 ) THEN 371 391 write(numout,*) 'icevar: jk, sz_i(ji,jj,jk,jl), s_i(ji,jj,jl), weights(jk) = ',jk, ' ',sz_i(ji,jj,jk,jl), ' ',s_i(ji,jj,jl), ' ',weights(jk) 372 392 ENDIF 393 373 394 END DO 374 395 END DO -
NEMO/branches/UKMO/NEMO_4.0_add_pond_lids_prints/src/ICE/icewri.F90
r11081 r12371 155 155 IF( iom_use('iceconc_cat' ) ) CALL iom_put( "iceconc_cat" , a_i * zmsk00l ) ! area for categories 156 156 IF( iom_use('icethic_cat' ) ) CALL iom_put( "icethic_cat" , h_i * zmsk00l ) ! thickness for categories 157 IF( iom_use('icevol_cat' ) ) CALL iom_put( "icevol_cat" , v_i * zmsk00l ) ! volume for categories 157 158 IF( iom_use('snwthic_cat' ) ) CALL iom_put( "snwthic_cat" , h_s * zmsksnl ) ! snow depth for categories 158 159 IF( iom_use('icesalt_cat' ) ) CALL iom_put( "icesalt_cat" , s_i * zmsk00l ) ! salinity for categories
Note: See TracChangeset
for help on using the changeset viewer.
