Changeset 15549
- Timestamp:
- 2021-11-28T21:00:36+01:00 (3 years ago)
- Location:
- NEMO/trunk
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/trunk/cfgs/SHARED/namelist_ice_ref
r14247 r15549 73 73 !------------------------------------------------------------------------------ 74 74 ! -- ice_rdgrft_strength -- ! 75 ln_str_H79 = .true. ! ice strength param.: Hibler_79 => P = pstar*<h>*exp(-c_rhg*A) 75 ln_str_H79 = .true. ! ice strength param.: Hibler_79 => P = pstar*<h>*exp(-c_rhg*A) 76 76 rn_pstar = 2.0e+04 ! ice strength thickness parameter [N/m2] 77 77 rn_crhg = 20.0 ! ice strength conc. parameter (-) 78 ln_str_R75 = .false. ! ice strength param.: Rothrock_75 => P = fn of potential energy 79 rn_pe_rdg = 17.0 ! coef accouting for frictional dissipation 80 ln_str_CST = .false. ! ice strength param.: Constant 81 rn_str = 0.0 ! ice strength value 82 ln_str_smooth = .true. ! spatial smoothing of the ice strength 78 83 ! -- ice_rdgrft -- ! 84 ln_distf_lin = .true. ! redistribution function of ridged ice: linear (Hibler 1980) 85 ln_distf_exp = .false. ! redistribution function of ridged ice: exponential => not coded yet 86 rn_murdg = 3.0 ! e-folding scale of ridged ice (m**.5) 79 87 rn_csrdg = 0.5 ! fraction of shearing energy contributing to ridging 80 88 ! -- ice_rdgrft_prep -- ! … … 139 147 ln_cndflx = .false. ! Use conduction flux as surface boundary conditions (i.e. for Jules coupling) 140 148 ln_cndemulate = .false. ! emulate conduction flux (if not provided in the inputs) 141 nn_qtrice = 1! Solar flux transmitted thru the surface scattering layer:149 nn_qtrice = 0 ! Solar flux transmitted thru the surface scattering layer: 142 150 ! = 0 Grenfell and Maykut 1977 (depends on cloudiness and is 0 when there is snow) 143 151 ! = 1 Lebrun 2019 (equals 0.3 anytime with different melting/dry snw conductivities) … … 268 276 rn_alb_imlt = 0.50 ! bare puddled ice albedo : 0.49 -- 0.58 269 277 rn_alb_dpnd = 0.27 ! ponded ice albedo : 0.10 -- 0.30 278 rn_alb_hpiv = 1.00 ! pivotal ice thickness in m (above which albedo is constant) 270 279 / 271 280 !------------------------------------------------------------------------------ -
NEMO/trunk/src/ICE/icealb.F90
r15334 r15549 38 38 REAL(wp) :: rn_alb_imlt ! bare puddled ice albedo 39 39 REAL(wp) :: rn_alb_dpnd ! ponded ice albedo 40 REAL(wp) :: rn_alb_hpiv ! pivotal ice thickness in meters (above which albedo is constant) 40 41 41 42 !! * Substitutions … … 59 60 !! which are an update of Allison et al. (JGR 1993) ; Brandt et al. 1999 60 61 !! 0-5cm : linear function of ice thickness 61 !! 5-1 50cm: log function of ice thickness62 !! > 1 50cm: constant62 !! 5-100cm: log function of ice thickness 63 !! > 100cm: constant 63 64 !! 2) Albedo dependency on snow thickness follows the findings from Grenfell & Perovich (2004) 64 65 !! i.e. it increases as -EXP(-snw_thick/0.02) during freezing and -EXP(-snw_thick/0.03) during melting … … 111 112 REAL(wp) :: zalb_snw, zafrac_snw ! snow-covered sea ice albedo & relative snow fraction 112 113 REAL(wp) :: zalb_cs, zalb_os ! albedo of ice under clear/overcast sky 113 !! clem114 REAL(wp), PARAMETER :: zhi_albcst = 1.5_wp ! pivotal thickness (should be in the namelist)115 114 !!--------------------------------------------------------------------- 116 115 ! … … 118 117 ! 119 118 z1_href_pnd = 1._wp / 0.05_wp 120 z1_c1 = 1._wp / ( LOG( zhi_albcst) - LOG(0.05_wp) )119 z1_c1 = 1._wp / ( LOG(rn_alb_hpiv) - LOG(0.05_wp) ) 121 120 z1_c2 = 1._wp / 0.05_wp 122 121 z1_c3 = 1._wp / 0.02_wp … … 142 141 !--- Albedos ---! 143 142 !---------------! 144 ! !--- Bare ice albedo (for hi > 1 50cm)143 ! !--- Bare ice albedo (for hi > 100cm) 145 144 IF( ld_pnd_alb ) THEN 146 145 zalb_ice = rn_alb_idry … … 149 148 ELSE ; zalb_ice = rn_alb_idry ; ENDIF 150 149 ENDIF 151 ! !--- Bare ice albedo (for hi < 1 50cm)152 IF( 0.05 < ph_ice(ji,jj,jl) .AND. ph_ice(ji,jj,jl) <= zhi_albcst ) THEN ! 5cm < hi < 150cm153 zalb_ice = zalb_ice + ( 0.18_wp - zalb_ice ) * z1_c1 * ( LOG( zhi_albcst) - LOG(ph_ice(ji,jj,jl)) )154 ELSEIF( ph_ice(ji,jj,jl) <= 0.05_wp ) THEN ! 0cm < hi < 5cm150 ! !--- Bare ice albedo (for hi < 100cm) 151 IF( 0.05 < ph_ice(ji,jj,jl) .AND. ph_ice(ji,jj,jl) <= rn_alb_hpiv ) THEN ! 5cm < hi < 100cm 152 zalb_ice = zalb_ice + ( 0.18_wp - zalb_ice ) * z1_c1 * ( LOG(rn_alb_hpiv) - LOG(ph_ice(ji,jj,jl)) ) 153 ELSEIF( ph_ice(ji,jj,jl) <= 0.05_wp ) THEN ! 0cm < hi < 5cm 155 154 zalb_ice = rn_alb_oce + ( 0.18_wp - rn_alb_oce ) * z1_c2 * ph_ice(ji,jj,jl) 156 155 ENDIF … … 193 192 INTEGER :: ios ! Local integer output status for namelist read 194 193 !! 195 NAMELIST/namalb/ rn_alb_sdry, rn_alb_smlt, rn_alb_idry, rn_alb_imlt, rn_alb_dpnd 194 NAMELIST/namalb/ rn_alb_sdry, rn_alb_smlt, rn_alb_idry, rn_alb_imlt, rn_alb_dpnd, rn_alb_hpiv 196 195 !!---------------------------------------------------------------------- 197 196 ! … … 212 211 WRITE(numout,*) ' albedo of bare puddled ice rn_alb_imlt = ', rn_alb_imlt 213 212 WRITE(numout,*) ' albedo of ponded ice rn_alb_dpnd = ', rn_alb_dpnd 213 WRITE(numout,*) ' pivotal ice thickness (m) rn_alb_hpiv = ', rn_alb_hpiv 214 214 ENDIF 215 215 ! -
NEMO/trunk/src/ICE/icedyn_rdgrft.F90
r15340 r15549 38 38 PUBLIC ice_strength ! called by icedyn_rhg_evp 39 39 40 INTEGER :: nice_str ! choice of the type of strength 41 ! ! associated indices: 42 INTEGER, PARAMETER :: np_strh79 = 1 ! Hibler 79 43 INTEGER, PARAMETER :: np_strr75 = 2 ! Rothrock 75 44 INTEGER, PARAMETER :: np_strcst = 3 ! Constant value 45 40 46 ! Variables shared among ridging subroutines 41 47 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: closing_net ! net rate at which area is removed (1/s) … … 57 63 ! 58 64 ! ** namelist (namdyn_rdgrft) ** 65 LOGICAL :: ln_str_R75 ! ice strength parameterization: Rothrock 75 66 REAL(wp) :: rn_pe_rdg ! coef accounting for frictional dissipation 67 LOGICAL :: ln_str_CST ! ice strength parameterization: Constant 68 REAL(wp) :: rn_str ! constant value of ice strength 69 LOGICAL :: ln_str_smooth ! ice strength spatial smoothing 70 LOGICAL :: ln_distf_lin ! redistribution of ridged ice: linear (Hibler 1980) 71 LOGICAL :: ln_distf_exp ! redistribution of ridged ice: exponential 72 REAL(wp) :: rn_murdg ! gives e-folding scale of ridged ice (m^.5) 59 73 REAL(wp) :: rn_csrdg ! fraction of shearing energy contributing to ridging 60 74 LOGICAL :: ln_partf_lin ! participation function linear (Thorndike et al. (1975)) … … 774 788 !!---------------------------------------------------------------------- 775 789 INTEGER :: ji, jj, jl ! dummy loop indices 776 INTEGER :: ismooth ! smoothing the resistance to deformation 777 INTEGER :: itframe ! number of time steps for the P smoothing 778 REAL(wp) :: zp, z1_3 ! local scalars 790 REAL(wp) :: z1_3 ! local scalars 779 791 REAL(wp), DIMENSION(jpi,jpj) :: zworka ! temporary array used here 780 REAL(wp), DIMENSION(jpi,jpj) :: zstrp1, zstrp2 ! strength at previous time steps 792 !!clem 793 LOGICAL :: ln_str_R75 794 REAL(wp) :: zhi, zcp, rn_pe_rdg 795 REAL(wp), DIMENSION(jpij) :: zstrength, zaksum ! strength in 1D 781 796 !!---------------------------------------------------------------------- 782 ! !--------------------------------------------------! 783 IF( ln_str_H79 ) THEN ! Ice strength => Hibler (1979) method ! 784 ! !--------------------------------------------------! 797 ! 798 SELECT CASE( nice_str ) !--- Set which ice strength is chosen 799 800 CASE ( np_strr75 ) !== Rothrock(1975)'s method ==! 801 802 ! these 2 param should be defined once for all at the 1st time step 803 zcp = 0.5_wp * grav * (rho0-rhoi) * rhoi * r1_rho0 ! proport const for PE 804 ! 805 strength(:,:) = 0._wp 806 ! 807 ! Identify grid cells with ice 808 at_i(:,:) = SUM( a_i, dim=3 ) 809 npti = 0 ; nptidx(:) = 0 810 DO_2D( nn_hls, nn_hls, nn_hls, nn_hls ) 811 IF ( at_i(ji,jj) > epsi10 ) THEN 812 npti = npti + 1 813 nptidx( npti ) = (jj - 1) * jpi + ji 814 ENDIF 815 END_2D 816 817 IF( npti > 0 ) THEN 818 CALL tab_3d_2d( npti, nptidx(1:npti), a_i_2d (1:npti,1:jpl), a_i ) 819 CALL tab_3d_2d( npti, nptidx(1:npti), v_i_2d (1:npti,1:jpl), v_i ) 820 CALL tab_2d_1d( npti, nptidx(1:npti), ato_i_1d(1:npti) , ato_i ) 821 CALL tab_2d_1d( npti, nptidx(1:npti), zstrength(1:npti) , strength ) 822 823 CALL rdgrft_prep( a_i_2d, v_i_2d, ato_i_1d, closing_net ) 824 ! 825 zaksum(1:npti) = apartf(1:npti,0) !clem: aksum should be defined in the header => local to module 826 DO jl = 1, jpl 827 DO ji = 1, npti 828 IF ( apartf(ji,jl) > 0._wp ) THEN 829 zaksum(ji) = zaksum(ji) + aridge(ji,jl) * ( 1._wp - hi_hrdg(ji,jl) ) & 830 & + araft (ji,jl) * ( 1._wp - hi_hrft ) 831 ENDIF 832 END DO 833 END DO 834 ! 835 z1_3 = 1._wp / 3._wp 836 DO jl = 1, jpl 837 DO ji = 1, npti 838 ! 839 IF( apartf(ji,jl) > 0._wp ) THEN 840 ! 841 IF( a_i_2d(ji,jl) > epsi10 ) THEN ; zhi = v_i_2d(ji,jl) / a_i_2d(ji,jl) 842 ELSE ; zhi = 0._wp 843 ENDIF 844 zstrength(ji) = zstrength(ji) - apartf(ji,jl) * zhi * zhi ! PE loss from deforming ice 845 zstrength(ji) = zstrength(ji) + 2._wp * araft (ji,jl) * zhi * zhi ! PE gain from rafting ice 846 zstrength(ji) = zstrength(ji) + aridge(ji,jl) * hi_hrdg(ji,jl) * z1_3 * & ! PE gain from ridging ice 847 & ( hrmax(ji,jl) * hrmax(ji,jl) + & ! (a**3-b**3)/(a-b) = a*a+ab+b*b 848 & hrmin(ji,jl) * hrmin(ji,jl) + & 849 & hrmax(ji,jl) * hrmin(ji,jl) ) 850 ENDIF 851 ! 852 END DO 853 END DO 854 ! 855 zstrength(1:npti) = rn_pe_rdg * zcp * zstrength(1:npti) / zaksum(1:npti) 856 ! 857 CALL tab_1d_2d( npti, nptidx(1:npti), zstrength(1:npti), strength ) 858 ! 859 ENDIF 860 ! 861 CASE ( np_strh79 ) !== Hibler(1979)'s method ==! 785 862 strength(:,:) = rn_pstar * SUM( v_i(:,:,:), dim=3 ) * EXP( -rn_crhg * ( 1._wp - SUM( a_i(:,:,:), dim=3 ) ) ) 786 ismooth = 1 ! original code 787 ! ismooth = 0 ! try for EAP stability 788 ! !--------------------------------------------------! 789 ELSE ! Zero strength ! 790 ! !--------------------------------------------------! 791 strength(:,:) = 0._wp 792 ismooth = 0 793 ENDIF 794 ! !--------------------------------------------------! 795 SELECT CASE( ismooth ) ! Smoothing ice strength ! 796 ! !--------------------------------------------------! 797 CASE( 1 ) !--- Spatial smoothing 863 ! 864 CASE ( np_strcst ) !== Constant strength ==! 865 strength(:,:) = rn_str 866 ! 867 END SELECT 868 ! 869 IF( ln_str_smooth ) THEN !--- Spatial smoothing 798 870 DO_2D( 0, 0, 0, 0 ) 799 871 IF ( SUM( a_i(ji,jj,:) ) > 0._wp ) THEN 800 zworka(ji,jj) = ( 4. 0* strength(ji,jj) &801 & + strength(ji-1,jj) * tmask(ji-1,jj,1) + strength(ji+1,jj) * tmask(ji+1,jj,1) &802 & + strength(ji,jj-1) * tmask(ji,jj-1,1) + strength(ji,jj+1) * tmask(ji,jj+1,1) &803 & ) / ( 4. 0+ tmask(ji-1,jj,1) + tmask(ji+1,jj,1) + tmask(ji,jj-1,1) + tmask(ji,jj+1,1) )872 zworka(ji,jj) = ( 4._wp * strength(ji,jj) & 873 & + strength(ji-1,jj) * tmask(ji-1,jj,1) + strength(ji+1,jj) * tmask(ji+1,jj,1) & 874 & + strength(ji,jj-1) * tmask(ji,jj-1,1) + strength(ji,jj+1) * tmask(ji,jj+1,1) & 875 & ) / ( 4._wp + tmask(ji-1,jj,1) + tmask(ji+1,jj,1) + tmask(ji,jj-1,1) + tmask(ji,jj+1,1) ) 804 876 ELSE 805 877 zworka(ji,jj) = 0._wp … … 812 884 CALL lbc_lnk( 'icedyn_rdgrft', strength, 'T', 1.0_wp ) 813 885 ! 814 CASE( 2 ) !--- Temporal smoothing 815 IF ( kt_ice == nit000 ) THEN 816 zstrp1(:,:) = 0._wp 817 zstrp2(:,:) = 0._wp 818 ENDIF 819 ! 820 DO_2D( 0, 0, 0, 0 ) 821 IF ( SUM( a_i(ji,jj,:) ) > 0._wp ) THEN 822 itframe = 1 ! number of time steps for the running mean 823 IF ( zstrp1(ji,jj) > 0._wp ) itframe = itframe + 1 824 IF ( zstrp2(ji,jj) > 0._wp ) itframe = itframe + 1 825 zp = ( strength(ji,jj) + zstrp1(ji,jj) + zstrp2(ji,jj) ) / itframe 826 zstrp2 (ji,jj) = zstrp1 (ji,jj) 827 zstrp1 (ji,jj) = strength(ji,jj) 828 strength(ji,jj) = zp 829 ENDIF 830 END_2D 831 CALL lbc_lnk( 'icedyn_rdgrft', strength, 'T', 1.0_wp ) 832 ! 833 END SELECT 886 ENDIF 834 887 ! 835 888 END SUBROUTINE ice_strength … … 924 977 !! ** input : Namelist namdyn_rdgrft 925 978 !!------------------------------------------------------------------- 926 INTEGER :: ios ! Local integer output status for namelist read 927 !! 928 NAMELIST/namdyn_rdgrft/ ln_str_H79, rn_pstar, rn_crhg, & 929 & rn_csrdg , & 930 & ln_partf_lin, rn_gstar, & 931 & ln_partf_exp, rn_astar, & 979 INTEGER :: ios, ioptio ! Local integer output status for namelist read 980 !! 981 NAMELIST/namdyn_rdgrft/ ln_str_H79, rn_pstar, rn_crhg, ln_str_R75, rn_pe_rdg, ln_str_CST, rn_str, ln_str_smooth, & 982 & ln_distf_lin, ln_distf_exp, rn_murdg, rn_csrdg, & 983 & ln_partf_lin, rn_gstar, ln_partf_exp, rn_astar, & 932 984 & ln_ridging, rn_hstar, rn_porordg, rn_fsnwrdg, rn_fpndrdg, & 933 985 & ln_rafting, rn_hraft, rn_craft , rn_fsnwrft, rn_fpndrft … … 948 1000 WRITE(numout,*) ' 1st bulk-rheology parameter rn_pstar = ', rn_pstar 949 1001 WRITE(numout,*) ' 2nd bulk-rhelogy parameter rn_crhg = ', rn_crhg 1002 WRITE(numout,*) ' ice strength parameterization Rothrock (1975) ln_str_R75 = ', ln_str_R75 1003 WRITE(numout,*) ' coef accounting for frictional dissipation rn_pe_rdg = ', rn_pe_rdg 1004 WRITE(numout,*) ' ice strength parameterization Constant ln_str_CST = ', ln_str_CST 1005 WRITE(numout,*) ' ice strength value rn_str = ', rn_str 1006 WRITE(numout,*) ' spatial smoothing of the strength ln_str_smooth= ', ln_str_smooth 1007 WRITE(numout,*) ' redistribution of ridged ice: linear (Hibler 1980) ln_distf_lin = ', ln_distf_lin 1008 WRITE(numout,*) ' redistribution of ridged ice: exponential ln_distf_exp = ', ln_distf_exp 1009 WRITE(numout,*) ' e-folding scale of ridged ice rn_murdg = ', rn_murdg 950 1010 WRITE(numout,*) ' Fraction of shear energy contributing to ridging rn_csrdg = ', rn_csrdg 951 1011 WRITE(numout,*) ' linear ridging participation function ln_partf_lin = ', ln_partf_lin … … 965 1025 ENDIF 966 1026 ! 1027 ioptio = 0 1028 IF( ln_str_H79 ) THEN ; ioptio = ioptio + 1 ; nice_str = np_strh79 ; ENDIF 1029 IF( ln_str_R75 ) THEN ; ioptio = ioptio + 1 ; nice_str = np_strr75 ; ENDIF 1030 IF( ln_str_CST ) THEN ; ioptio = ioptio + 1 ; nice_str = np_strcst ; ENDIF 1031 IF( ioptio /= 1 ) CALL ctl_stop( 'ice_dyn_rdgrft_init: one and only one ice strength option has to be defined ' ) 1032 ! 1033 IF ( ( ln_distf_lin .AND. ln_distf_exp ) .OR. ( .NOT.ln_distf_lin .AND. .NOT.ln_distf_exp ) ) THEN 1034 CALL ctl_stop( 'ice_dyn_rdgrft_init: choose one and only one redistribution function (ln_distf_lin or ln_distf_exp)' ) 1035 ENDIF 1036 !!clem 1037 IF( ln_distf_exp ) CALL ctl_stop( 'ice_dyn_rdgrft_init: exponential redistribution function not coded yet (ln_distf_exp)' ) 1038 ! 967 1039 IF ( ( ln_partf_lin .AND. ln_partf_exp ) .OR. ( .NOT.ln_partf_lin .AND. .NOT.ln_partf_exp ) ) THEN 968 1040 CALL ctl_stop( 'ice_dyn_rdgrft_init: choose one and only one participation function (ln_partf_lin or ln_partf_exp)' )
Note: See TracChangeset
for help on using the changeset viewer.