Changeset 1572 for trunk/NEMO/LIM_SRC_3/limthd.F90
- Timestamp:
- 2009-08-03T16:53:15+02:00 (15 years ago)
- File:
-
- 1 edited
-
trunk/NEMO/LIM_SRC_3/limthd.F90 (modified) (32 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/NEMO/LIM_SRC_3/limthd.F90
r1571 r1572 2 2 !!====================================================================== 3 3 !! *** MODULE limthd *** 4 !! LIM thermo ice model :ice thermodynamic4 !! LIM-3 : ice thermodynamic 5 5 !!====================================================================== 6 !! History : LIM ! 2000-01 (M.A. Morales Maqueda, H. Goosse, T. Fichefet) LIM-1 7 !! 2.0 ! 2002-07 (C. Ethe, G. Madec) LIM-2 (F90 rewriting) 8 !! 3.0 ! 2005-11 (M. Vancoppenolle) LIM-3 : Multi-layer thermodynamics + salinity variations 9 !! - ! 2007-04 (M. Vancoppenolle) add lim_thd_glohec and lim_thd_con_dif 10 !! 3.2 ! 2009-07 (M. Vancoppenolle, Y. Aksenov, G. Madec) bug correction in rdmsnif 11 !!---------------------------------------------------------------------- 6 12 #if defined key_lim3 7 13 !!---------------------------------------------------------------------- … … 11 17 !! lim_thd_init : initialisation of sea-ice thermodynamic 12 18 !!---------------------------------------------------------------------- 13 !! * Modules used14 19 USE phycst ! physical constants 15 20 USE dom_oce ! ocean space and time domain variables 16 USE domvvl ! Variable volume17 USE lbclnk18 USE in_out_manager ! I/O manager19 21 USE ice ! LIM sea-ice variables 22 USE par_ice 20 23 USE sbc_oce ! Surface boundary condition: ocean fields 21 24 USE sbc_ice ! Surface boundary condition: ice fields 22 25 USE thd_ice ! LIM thermodynamic sea-ice variables 23 26 USE dom_ice ! LIM sea-ice domain 27 USE domvvl ! Variable volume 24 28 USE iceini 25 29 USE limthd_dif … … 28 32 USE limthd_ent 29 33 USE limtab 30 USE par_ice31 34 USE limvar 35 USE in_out_manager ! I/O manager 32 36 USE prtctl ! Print control 37 USE lbclnk 33 38 USE lib_mpp 34 39 … … 36 41 PRIVATE 37 42 38 !! * Routine accessibility 39 PUBLIC lim_thd ! called by lim_step 40 41 !! * Module variables 42 REAL(wp) :: & ! constant values 43 epsi20 = 1e-20 , & 44 epsi16 = 1e-16 , & 45 epsi06 = 1e-06 , & 46 epsi04 = 1e-04 , & 47 zzero = 0.e0 , & 48 zone = 1.e0 43 PUBLIC lim_thd ! called by lim_step 44 45 REAL(wp) :: epsi20 = 1e-20 ! constant values 46 REAL(wp) :: epsi16 = 1e-16 ! 47 REAL(wp) :: epsi06 = 1e-06 ! 48 REAL(wp) :: epsi04 = 1e-04 ! 49 REAL(wp) :: zzero = 0.e0 ! 50 REAL(wp) :: zone = 1.e0 ! 49 51 50 52 !! * Substitutions … … 52 54 # include "vectopt_loop_substitute.h90" 53 55 !!---------------------------------------------------------------------- 54 !! LIM 3.0, UCL-LOCEAN-IPSL (2005)56 !! NEMO/LIM 3.2 , UCL-LOCEAN-IPSL (2009) 55 57 !! $Id$ 56 58 !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) … … 75 77 !! - back to the geographic grid 76 78 !! 77 !! ** References : 78 !! H. Goosse et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90 79 !! ** References : H. Goosse et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90 80 !!--------------------------------------------------------------------- 81 INTEGER, INTENT(in) :: kt ! number of iteration 79 82 !! 80 !! History : 81 !! 1.0 ! 00-01 (M.A. Morales Maqueda, H. Goosse, T. Fichefet) 82 !! 2.0 ! 02-07 (C. Ethe, G. Madec) F90 83 !! 3.0 ! 05-11 (M. Vancoppenolle ) Multi-layer thermodynamics, 84 !! salinity variations 85 !!--------------------------------------------------------------------- 86 INTEGER, INTENT(in) :: kt ! number of iteration 87 !! * Local variables 88 INTEGER :: ji, jj, jk, jl, nbpb ! nb of icy pts for thermo. cal. 89 90 REAL(wp) :: & 91 zfric_umin = 5e-03 , & ! lower bound for the friction velocity 92 zfric_umax = 2e-02 ! upper bound for the friction velocity 93 94 REAL(wp) :: & 95 zinda , & ! switch for test. the val. of concen. 96 zindb, & ! switches for test. the val of arg 97 zthsnice , & 98 zfric_u , & ! friction velocity 99 zfnsol , & ! total non solar heat 100 zfontn , & ! heat flux from snow thickness 101 zfntlat, zpareff , & ! test. the val. of lead heat budget 102 zeps 103 104 REAL(wp) :: & 105 zareamin 106 83 INTEGER :: ji, jj, jk, jl ! dummy loop indices 84 INTEGER :: nbpb ! nb of icy pts for thermo. cal. 85 REAL(wp) :: zfric_umin = 5e-03 ! lower bound for the friction velocity 86 REAL(wp) :: zfric_umax = 2e-02 ! upper bound for the friction velocity 87 REAL(wp) :: zinda, zindb, zthsnice, zfric_u ! temporary scalar 88 REAL(wp) :: zfnsol, zfontn, zfntlat, zpareff ! - - 89 REAL(wp) :: zeps, zareamin, zcoef 107 90 REAL(wp), DIMENSION(jpi,jpj) :: zqlbsbq ! link with lead energy budget qldif 108 109 91 !!------------------------------------------------------------------- 110 92 111 IF( numit == nstart ) CALL lim_thd_init ! Initialization (first time-step only) 112 113 IF( kt == nit000 .AND. lwp ) THEN 114 WRITE(numout,*) 'limthd : Ice Thermodynamics' 115 WRITE(numout,*) '~~~~~~' 116 ENDIF 117 118 IF( numit == nstart ) CALL lim_thd_sal_init ! Initialization (first time-step only) 93 IF( numit == nstart ) CALL lim_thd_init ! Initialization (first time-step only) 94 95 IF( numit == nstart ) CALL lim_thd_sal_init ! Initialization (first time-step only) 96 119 97 !------------------------------------------------------------------------------! 120 98 ! 1) Initialization of diagnostic variables ! … … 125 103 ! 1.2) Heat content 126 104 !-------------------- 127 ! Change the units of heat content; from global units to 128 ! J.m3 105 ! Change the units of heat content; from global units to J.m3 129 106 130 107 DO jl = 1, jpl … … 133 110 DO ji = 1, jpi 134 111 !Energy of melting q(S,T) [J.m-3] 135 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / area(ji,jj) / & 136 MAX( v_i(ji,jj,jl) , epsi06 ) * nlay_i 112 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / ( area(ji,jj) * MAX( v_i(ji,jj,jl) , epsi06 ) ) * nlay_i 137 113 !0 if no ice and 1 if yes 138 zindb 114 zindb = 1.0 - MAX ( 0.0 , SIGN ( 1.0 , - ht_i(ji,jj,jl) ) ) 139 115 !convert units ! very important that this line is here 140 116 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * unit_fac * zindb … … 142 118 END DO 143 119 END DO 144 END DO145 146 DO jl = 1, jpl147 120 DO jk = 1, nlay_s 148 121 DO jj = 1, jpj 149 122 DO ji = 1, jpi 150 123 !Energy of melting q(S,T) [J.m-3] 151 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) / area(ji,jj) / & 152 MAX( v_s(ji,jj,jl) , epsi06 ) * nlay_s 124 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) / ( area(ji,jj) * MAX( v_s(ji,jj,jl) , epsi06 ) ) * nlay_s 153 125 !0 if no ice and 1 if yes 154 zindb 126 zindb = 1.0 - MAX ( 0.0 , SIGN ( 1.0 , - ht_s(ji,jj,jl) ) ) 155 127 !convert units ! very important that this line is here 156 128 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * unit_fac * zindb … … 180 152 ! 1.4) Compute global heat content 181 153 !----------------------------------- 182 qt_i_in (:,:)= 0.e0183 qt_s_in (:,:)= 0.e0184 qt_i_fin (:,:)= 0.e0185 qt_s_fin (:,:)= 0.e0154 qt_i_in (:,:) = 0.e0 155 qt_s_in (:,:) = 0.e0 156 qt_i_fin (:,:) = 0.e0 157 qt_s_fin (:,:) = 0.e0 186 158 sum_fluxq(:,:) = 0.e0 187 fatm (:,:) = 0.e0159 fatm (:,:) = 0.e0 188 160 189 161 ! 2) Partial computation of forcing for the thermodynamic sea ice model. ! … … 220 192 zfontn = sprecip(ji,jj) * lfus ! energy of melting 221 193 zfnsol = qns(ji,jj) ! total non solar flux 222 qldif(ji,jj) = tms(ji,jj) * ( qsr(ji,jj) &223 & + zfnsol + fdtcn(ji,jj) - zfontn &224 & + ( 1.0 - zindb ) * fsbbq(ji,jj) ) &194 qldif(ji,jj) = tms(ji,jj) * ( qsr(ji,jj) & 195 & + zfnsol + fdtcn(ji,jj) - zfontn & 196 & + ( 1.0 - zindb ) * fsbbq(ji,jj) ) & 225 197 & * ( 1.0 - at_i(ji,jj) ) * rdt_ice 226 198 … … 229 201 != 1 if positive heat budget 230 202 zpareff = 1.0 - zinda * zfntlat 231 != 0 if ice and positive heat budget and 1 if one of those two is 232 !false 233 zqlbsbq(ji,jj) = qldif(ji,jj) * ( 1.0 - zpareff ) / & 234 MAX( at_i(ji,jj) * rdt_ice , epsi16 ) 203 != 0 if ice and positive heat budget and 1 if one of those two is false 204 zqlbsbq(ji,jj) = qldif(ji,jj) * ( 1.0 - zpareff ) / MAX( at_i(ji,jj) * rdt_ice , epsi16 ) 235 205 236 206 ! Heat budget of the lead, energy transferred from ice to ocean … … 238 208 qdtcn (ji,jj) = zpareff * qdtcn(ji,jj) 239 209 240 ! Energy needed to bring ocean surface layer until its freezing 241 ! qcmif, limflx 242 qcmif (ji,jj) = rau0 * rcp * fse3t_m(ji,jj,1) & 243 & * ( t_bo(ji,jj) - (sst_m(ji,jj) + rt0) ) * ( 1. - zinda ) 244 245 ! calculate oceanic heat flux (limthd_dh) 210 ! Energy needed to bring ocean surface layer until its freezing (qcmif, limflx) 211 qcmif (ji,jj) = rau0 * rcp * fse3t(ji,jj,1) * ( t_bo(ji,jj) - (sst_m(ji,jj) + rt0) ) * ( 1. - zinda ) 212 213 ! oceanic heat flux (limthd_dh) 246 214 fbif (ji,jj) = zindb * ( fsbbq(ji,jj) / MAX( at_i(ji,jj) , epsi20 ) + fdtcn(ji,jj) ) 247 215 ! … … 279 247 !------------------------------------------------------------------------------! 280 248 281 IF( lk_mpp ) CALL mpp_ini_ice(nbpb)282 283 IF (nbpb > 0) THEN ! If there is no ice, do nothing.249 IF( lk_mpp ) CALL mpp_ini_ice( nbpb ) 250 251 IF( nbpb > 0 ) THEN ! If there is no ice, do nothing. 284 252 285 253 !------------------------- … … 287 255 !------------------------- 288 256 289 CALL tab_2d_1d( nbpb, at_i_b (1:nbpb) 290 CALL tab_2d_1d( nbpb, a_i_b (1:nbpb) 291 CALL tab_2d_1d( nbpb, ht_i_b (1:nbpb) 292 CALL tab_2d_1d( nbpb, ht_s_b (1:nbpb) 293 294 CALL tab_2d_1d( nbpb, t_su_b (1:nbpb) , t_su(:,:,jl), jpi, jpj, npb(1:nbpb) )295 CALL tab_2d_1d( nbpb, sm_i_b (1:nbpb) , sm_i(:,:,jl), jpi, jpj, npb(1:nbpb) )257 CALL tab_2d_1d( nbpb, at_i_b (1:nbpb), at_i , jpi, jpj, npb(1:nbpb) ) 258 CALL tab_2d_1d( nbpb, a_i_b (1:nbpb), a_i(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 259 CALL tab_2d_1d( nbpb, ht_i_b (1:nbpb), ht_i(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 260 CALL tab_2d_1d( nbpb, ht_s_b (1:nbpb), ht_s(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 261 262 CALL tab_2d_1d( nbpb, t_su_b (1:nbpb), t_su(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 263 CALL tab_2d_1d( nbpb, sm_i_b (1:nbpb), sm_i(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 296 264 DO jk = 1, nlay_s 297 CALL tab_2d_1d( nbpb, t_s_b(1:nbpb,jk) , t_s(:,:,jk,jl), jpi, jpj, npb(1:nbpb) )298 CALL tab_2d_1d( nbpb, q_s_b(1:nbpb,jk) , e_s(:,:,jk,jl), jpi, jpj, npb(1:nbpb) )265 CALL tab_2d_1d( nbpb, t_s_b(1:nbpb,jk), t_s(:,:,jk,jl) , jpi, jpj, npb(1:nbpb) ) 266 CALL tab_2d_1d( nbpb, q_s_b(1:nbpb,jk), e_s(:,:,jk,jl) , jpi, jpj, npb(1:nbpb) ) 299 267 END DO 300 268 DO jk = 1, nlay_i 301 CALL tab_2d_1d( nbpb, t_i_b(1:nbpb,jk) , t_i(:,:,jk,jl), jpi, jpj, npb(1:nbpb) )302 CALL tab_2d_1d( nbpb, q_i_b(1:nbpb,jk) , e_i(:,:,jk,jl), jpi, jpj, npb(1:nbpb) )303 CALL tab_2d_1d( nbpb, s_i_b(1:nbpb,jk) , s_i(:,:,jk,jl), jpi, jpj, npb(1:nbpb) )269 CALL tab_2d_1d( nbpb, t_i_b(1:nbpb,jk), t_i(:,:,jk,jl) , jpi, jpj, npb(1:nbpb) ) 270 CALL tab_2d_1d( nbpb, q_i_b(1:nbpb,jk), e_i(:,:,jk,jl) , jpi, jpj, npb(1:nbpb) ) 271 CALL tab_2d_1d( nbpb, s_i_b(1:nbpb,jk), s_i(:,:,jk,jl) , jpi, jpj, npb(1:nbpb) ) 304 272 END DO 305 273 306 CALL tab_2d_1d( nbpb, tatm_ice_1d(1:nbpb) 307 CALL tab_2d_1d( nbpb, qsr_ice_1d (1:nbpb) 308 CALL tab_2d_1d( nbpb, fr1_i0_1d (1:nbpb) , fr1_i0, jpi, jpj, npb(1:nbpb) )309 CALL tab_2d_1d( nbpb, fr2_i0_1d (1:nbpb) , fr2_i0, jpi, jpj, npb(1:nbpb) )310 CALL tab_2d_1d( nbpb, qnsr_ice_1d(1:nbpb) 274 CALL tab_2d_1d( nbpb, tatm_ice_1d(1:nbpb), tatm_ice(:,:) , jpi, jpj, npb(1:nbpb) ) 275 CALL tab_2d_1d( nbpb, qsr_ice_1d (1:nbpb), qsr_ice(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 276 CALL tab_2d_1d( nbpb, fr1_i0_1d (1:nbpb), fr1_i0 , jpi, jpj, npb(1:nbpb) ) 277 CALL tab_2d_1d( nbpb, fr2_i0_1d (1:nbpb), fr2_i0 , jpi, jpj, npb(1:nbpb) ) 278 CALL tab_2d_1d( nbpb, qnsr_ice_1d(1:nbpb), qns_ice(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 311 279 312 280 #if ! defined key_coupled 313 CALL tab_2d_1d( nbpb, qla_ice_1d (1:nbpb) 314 CALL tab_2d_1d( nbpb, dqla_ice_1d(1:nbpb) 281 CALL tab_2d_1d( nbpb, qla_ice_1d (1:nbpb), qla_ice(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 282 CALL tab_2d_1d( nbpb, dqla_ice_1d(1:nbpb), dqla_ice(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 315 283 #endif 316 284 317 CALL tab_2d_1d( nbpb, dqns_ice_1d(1:nbpb) 318 CALL tab_2d_1d( nbpb, t_bo_b (1:nbpb) 319 CALL tab_2d_1d( nbpb, sprecip_1d (1:nbpb) 320 CALL tab_2d_1d( nbpb, fbif_1d (1:nbpb) 321 CALL tab_2d_1d( nbpb, qldif_1d (1:nbpb) 322 CALL tab_2d_1d( nbpb, rdmicif_1d (1:nbpb) 323 CALL tab_2d_1d( nbpb, rdmsnif_1d (1:nbpb) 324 CALL tab_2d_1d( nbpb, dmgwi_1d (1:nbpb) 325 CALL tab_2d_1d( nbpb, qlbbq_1d (1:nbpb) 326 327 CALL tab_2d_1d( nbpb, fseqv_1d (1:nbpb) 328 CALL tab_2d_1d( nbpb, fsbri_1d (1:nbpb) 329 CALL tab_2d_1d( nbpb, fhbri_1d (1:nbpb) 330 CALL tab_2d_1d( nbpb, fstbif_1d (1:nbpb) 331 CALL tab_2d_1d( nbpb, qfvbq_1d (1:nbpb) 285 CALL tab_2d_1d( nbpb, dqns_ice_1d(1:nbpb), dqns_ice(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 286 CALL tab_2d_1d( nbpb, t_bo_b (1:nbpb), t_bo , jpi, jpj, npb(1:nbpb) ) 287 CALL tab_2d_1d( nbpb, sprecip_1d (1:nbpb), sprecip , jpi, jpj, npb(1:nbpb) ) 288 CALL tab_2d_1d( nbpb, fbif_1d (1:nbpb), fbif , jpi, jpj, npb(1:nbpb) ) 289 CALL tab_2d_1d( nbpb, qldif_1d (1:nbpb), qldif , jpi, jpj, npb(1:nbpb) ) 290 CALL tab_2d_1d( nbpb, rdmicif_1d (1:nbpb), rdmicif , jpi, jpj, npb(1:nbpb) ) 291 CALL tab_2d_1d( nbpb, rdmsnif_1d (1:nbpb), rdmsnif , jpi, jpj, npb(1:nbpb) ) 292 CALL tab_2d_1d( nbpb, dmgwi_1d (1:nbpb), dmgwi , jpi, jpj, npb(1:nbpb) ) 293 CALL tab_2d_1d( nbpb, qlbbq_1d (1:nbpb), zqlbsbq , jpi, jpj, npb(1:nbpb) ) 294 295 CALL tab_2d_1d( nbpb, fseqv_1d (1:nbpb), fseqv , jpi, jpj, npb(1:nbpb) ) 296 CALL tab_2d_1d( nbpb, fsbri_1d (1:nbpb), fsbri , jpi, jpj, npb(1:nbpb) ) 297 CALL tab_2d_1d( nbpb, fhbri_1d (1:nbpb), fhbri , jpi, jpj, npb(1:nbpb) ) 298 CALL tab_2d_1d( nbpb, fstbif_1d (1:nbpb), fstric , jpi, jpj, npb(1:nbpb) ) 299 CALL tab_2d_1d( nbpb, qfvbq_1d (1:nbpb), qfvbq , jpi, jpj, npb(1:nbpb) ) 332 300 333 301 !-------------------------------- … … 335 303 !-------------------------------- 336 304 337 IF ( con_i ) CALL lim_thd_enmelt(1,nbpb) ! computes sea ice energy of melting 338 IF ( con_i ) CALL lim_thd_glohec( qt_i_in , qt_s_in , & 339 q_i_layer_in , 1 , nbpb , jl ) 340 341 !---------------------------------! 342 CALL lim_thd_dif(1,nbpb,jl) ! Ice/Snow Temperature profile ! 343 !---------------------------------! 344 345 CALL lim_thd_enmelt(1,nbpb) ! computes sea ice energy of melting 346 ! compulsory for limthd_dh 347 348 IF ( con_i ) CALL lim_thd_glohec( qt_i_fin , qt_s_fin , & 349 q_i_layer_fin , 1 , nbpb , jl ) 350 IF ( con_i ) CALL lim_thd_con_dif( 1 , nbpb , jl ) 351 352 !---------------------------------! 353 CALL lim_thd_dh(1,nbpb,jl) ! Ice/Snow thickness ! 354 !---------------------------------! 355 356 !---------------------------------! 357 CALL lim_thd_ent(1,nbpb,jl) ! Ice/Snow enthalpy remapping ! 358 !---------------------------------! 359 360 !---------------------------------! 361 CALL lim_thd_sal(1,nbpb) ! Ice salinity computation ! 362 !---------------------------------! 305 IF( con_i ) CALL lim_thd_enmelt( 1, nbpb ) ! computes sea ice energy of melting 306 IF( con_i ) CALL lim_thd_glohec( qt_i_in, qt_s_in, q_i_layer_in, 1, nbpb, jl ) 307 308 ! !---------------------------------! 309 CALL lim_thd_dif( 1, nbpb, jl ) ! Ice/Snow Temperature profile ! 310 ! !---------------------------------! 311 312 CALL lim_thd_enmelt( 1, nbpb ) ! computes sea ice energy of melting compulsory for limthd_dh 313 314 IF( con_i ) CALL lim_thd_glohec ( qt_i_fin, qt_s_fin, q_i_layer_fin, 1, nbpb, jl ) 315 IF( con_i ) CALL lim_thd_con_dif( 1 , nbpb , jl ) 316 317 ! !---------------------------------! 318 CALL lim_thd_dh( 1, nbpb, jl ) ! Ice/Snow thickness ! 319 ! !---------------------------------! 320 321 ! !---------------------------------! 322 CALL lim_thd_ent( 1, nbpb, jl ) ! Ice/Snow enthalpy remapping ! 323 ! !---------------------------------! 324 325 ! !---------------------------------! 326 CALL lim_thd_sal( 1, nbpb ) ! Ice salinity computation ! 327 ! !---------------------------------! 363 328 364 329 ! CALL lim_thd_enmelt(1,nbpb) ! computes sea ice energy of melting 365 IF ( con_i ) CALL lim_thd_glohec( qt_i_fin, qt_s_fin, & 366 q_i_layer_fin , 1 , nbpb , jl ) 367 IF ( con_i ) CALL lim_thd_con_dh ( 1 , nbpb , jl ) 330 IF( con_i ) CALL lim_thd_glohec( qt_i_fin, qt_s_fin, q_i_layer_fin, 1, nbpb, jl ) 331 IF( con_i ) CALL lim_thd_con_dh ( 1 , nbpb , jl ) 368 332 369 333 !-------------------------------- … … 371 335 !-------------------------------- 372 336 373 CALL tab_1d_2d( nbpb, at_i , npb, at_i_b 337 CALL tab_1d_2d( nbpb, at_i , npb, at_i_b(1:nbpb), jpi, jpj ) 374 338 CALL tab_1d_2d( nbpb, ht_i(:,:,jl), npb, ht_i_b(1:nbpb), jpi, jpj ) 375 339 CALL tab_1d_2d( nbpb, ht_s(:,:,jl), npb, ht_s_b(1:nbpb), jpi, jpj ) … … 389 353 END DO 390 354 391 CALL tab_1d_2d( nbpb, fstric , npb, fstbif_1d (1:nbpb), jpi, jpj )392 CALL tab_1d_2d( nbpb, qldif , npb, qldif_1d (1:nbpb), jpi, jpj )393 CALL tab_1d_2d( nbpb, qfvbq , npb, qfvbq_1d (1:nbpb), jpi, jpj )394 CALL tab_1d_2d( nbpb, rdmicif , npb, rdmicif_1d(1:nbpb), jpi, jpj )395 CALL tab_1d_2d( nbpb, rdmsnif , npb, rdmsnif_1d(1:nbpb), jpi, jpj )396 CALL tab_1d_2d( nbpb, dmgwi , npb, dmgwi_1d (1:nbpb), jpi, jpj )397 CALL tab_1d_2d( nbpb, rdvosif , npb, dvsbq_1d (1:nbpb), jpi, jpj )398 CALL tab_1d_2d( nbpb, rdvobif , npb, dvbbq_1d (1:nbpb), jpi, jpj )399 CALL tab_1d_2d( nbpb, fdvolif , npb, dvlbq_1d (1:nbpb), jpi, jpj )400 CALL tab_1d_2d( nbpb, rdvonif , npb, dvnbq_1d (1:nbpb), jpi, jpj )401 CALL tab_1d_2d( nbpb, fseqv , npb, fseqv_1d (1:nbpb), jpi, jpj )402 403 IF ( num_sal .EQ.2 ) THEN404 CALL tab_1d_2d( nbpb, fsbri , npb, fsbri_1d (1:nbpb), jpi, jpj )405 CALL tab_1d_2d( nbpb, fhbri , npb, fhbri_1d (1:nbpb), jpi, jpj )355 CALL tab_1d_2d( nbpb, fstric , npb, fstbif_1d (1:nbpb), jpi, jpj ) 356 CALL tab_1d_2d( nbpb, qldif , npb, qldif_1d (1:nbpb), jpi, jpj ) 357 CALL tab_1d_2d( nbpb, qfvbq , npb, qfvbq_1d (1:nbpb), jpi, jpj ) 358 CALL tab_1d_2d( nbpb, rdmicif, npb, rdmicif_1d(1:nbpb), jpi, jpj ) 359 CALL tab_1d_2d( nbpb, rdmsnif, npb, rdmsnif_1d(1:nbpb), jpi, jpj ) 360 CALL tab_1d_2d( nbpb, dmgwi , npb, dmgwi_1d (1:nbpb), jpi, jpj ) 361 CALL tab_1d_2d( nbpb, rdvosif, npb, dvsbq_1d (1:nbpb), jpi, jpj ) 362 CALL tab_1d_2d( nbpb, rdvobif, npb, dvbbq_1d (1:nbpb), jpi, jpj ) 363 CALL tab_1d_2d( nbpb, fdvolif, npb, dvlbq_1d (1:nbpb), jpi, jpj ) 364 CALL tab_1d_2d( nbpb, rdvonif, npb, dvnbq_1d (1:nbpb), jpi, jpj ) 365 CALL tab_1d_2d( nbpb, fseqv , npb, fseqv_1d (1:nbpb), jpi, jpj ) 366 367 IF( num_sal == 2 ) THEN 368 CALL tab_1d_2d( nbpb, fsbri, npb, fsbri_1d(1:nbpb), jpi, jpj ) 369 CALL tab_1d_2d( nbpb, fhbri, npb, fhbri_1d(1:nbpb), jpi, jpj ) 406 370 ENDIF 407 371 408 372 !+++++ 409 !temporary stuff for a dummy version373 !temporary stuff for a dummy version 410 374 CALL tab_1d_2d( nbpb, dh_i_surf2D, npb, dh_i_surf(1:nbpb) , jpi, jpj ) 411 375 CALL tab_1d_2d( nbpb, dh_i_bott2D, npb, dh_i_bott(1:nbpb) , jpi, jpj ) … … 417 381 !+++++ 418 382 419 IF( lk_mpp ) CALL mpp_comm_free(ncomm_ice) !RB necessary ??420 ENDIF ! nbpb421 422 END DO ! jl383 IF( lk_mpp ) CALL mpp_comm_free( ncomm_ice ) !RB necessary ?? 384 ENDIF 385 ! 386 END DO 423 387 424 388 !------------------------------------------------------------------------------! … … 431 395 432 396 ! Enthalpies are global variables we have to readjust the units 397 zcoef = 1.e0 / ( unit_fac * REAL(nlay_i) ) 433 398 DO jl = 1, jpl 434 399 DO jk = 1, nlay_i 435 DO jj = 1, jpj 436 DO ji = 1, jpi 437 ! Change dimensions 438 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / unit_fac 439 440 ! Multiply by volume, divide by nlayers so that heat content in 10^9 Joules 441 e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * & 442 area(ji,jj) * a_i(ji,jj,jl) * & 443 ht_i(ji,jj,jl) / nlay_i 444 END DO !ji 445 END DO !jj 446 END DO !jk 447 END DO !jl 400 ! Multiply by volume, divide by nlayers so that heat content in 10^9 Joules 401 e_i(:,:,jk,jl) = e_i(:,:,jk,jl) * area(:,:) * a_i(:,:,jl) * ht_i(:,:,jl) * zcoef 402 END DO 403 END DO 448 404 449 405 !------------------------ … … 452 408 453 409 ! Enthalpies are global variables we have to readjust the units 410 zcoef = 1.e0 / ( unit_fac * REAL(nlay_s) ) 454 411 DO jl = 1, jpl 455 412 DO jk = 1, nlay_s 456 DO jj = 1, jpj 457 DO ji = 1, jpi 458 ! Change dimensions 459 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) / unit_fac 460 ! Multiply by volume, so that heat content in 10^9 Joules 461 e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * area(ji,jj) * & 462 a_i(ji,jj,jl) * ht_s(ji,jj,jl) / nlay_s 463 END DO !ji 464 END DO !jj 465 END DO !jk 466 END DO !jl 413 ! Multiply by volume, so that heat content in 10^9 Joules 414 e_s(:,:,jk,jl) = e_s(:,:,jk,jl) * area(:,:) * a_i(:,:,jl) * ht_s(:,:,jl) * zcoef 415 END DO 416 END DO 467 417 468 418 !---------------------------------- … … 474 424 ! 5.4) Diagnostic thermodynamic growth rates 475 425 !-------------------------------------------- 476 d_v_i_thd (:,:,:) = v_i(:,:,:)- old_v_i(:,:,:) ! ice volumes477 dv_dt_thd(:,:,:) 478 479 IF ( con_i )fbif(:,:) = fbif(:,:) + zqlbsbq(:,:)426 d_v_i_thd(:,:,:) = v_i (:,:,:) - old_v_i(:,:,:) ! ice volumes 427 dv_dt_thd(:,:,:) = d_v_i_thd(:,:,:) / rdt_ice * 86400.0 428 429 IF( con_i ) fbif(:,:) = fbif(:,:) + zqlbsbq(:,:) 480 430 481 431 IF(ln_ctl) THEN ! Control print … … 514 464 END SUBROUTINE lim_thd 515 465 516 !=============================================================================== 517 518 SUBROUTINE lim_thd_glohec(eti,ets,etilayer,kideb,kiut,jl) 466 467 SUBROUTINE lim_thd_glohec( eti, ets, etilayer, kideb, kiut, jl ) 519 468 !!----------------------------------------------------------------------- 520 469 !! *** ROUTINE lim_thd_glohec *** … … 522 471 !! ** Purpose : Compute total heat content for each category 523 472 !! Works with 1d vectors only 473 !!----------------------------------------------------------------------- 474 INTEGER , INTENT(in ) :: kideb, kiut ! bounds for the spatial loop 475 INTEGER , INTENT(in ) :: jl ! category number 476 REAL(wp), INTENT( out), DIMENSION (jpij,jpl ) :: eti, ets ! vertically-summed heat content for ice & snow 477 REAL(wp), INTENT( out), DIMENSION (jpij,jkmax) :: etilayer ! heat content for ice layers 524 478 !! 525 !! history : 526 !! 9.9 ! 07-04 (M.Vancoppenolle) original code 479 INTEGER :: ji,jk ! loop indices 480 REAL(wp) :: zdes ! snow heat content increment (dummy) 481 REAL(wp) :: zeps ! very small value (1.e-10) 527 482 !!----------------------------------------------------------------------- 528 !! * Local variables 529 INTEGER, INTENT(in) :: & 530 kideb, kiut, & ! bounds for the spatial loop 531 jl ! category number 532 533 REAL(wp), DIMENSION (jpij,jpl), INTENT(out) :: & 534 eti, ets ! vertically-summed heat content for ice /snow 535 536 REAL(wp), DIMENSION (jpij,jkmax), INTENT(out) :: & 537 etilayer ! heat content for ice layers 538 539 REAL(wp) :: & 540 zdes, & ! snow heat content increment (dummy) 541 zeps ! very small value (1.e-10) 542 543 INTEGER :: & 544 ji,jk ! loop indices 545 546 !!----------------------------------------------------------------------- 547 eti(:,:) = 0.0 548 ets(:,:) = 0.0 483 eti(:,:) = 0.e0 484 ets(:,:) = 0.e0 549 485 zeps = 1.0e-10 550 486 551 ! total q over all layers, ice [J.m-2] 552 DO jk = 1, nlay_i 487 DO jk = 1, nlay_i ! total q over all layers, ice [J.m-2] 553 488 DO ji = kideb, kiut 554 etilayer(ji,jk) = q_i_b(ji,jk) & 555 * ht_i_b(ji) / nlay_i 556 eti(ji,jl) = eti(ji,jl) + etilayer(ji,jk) 489 etilayer(ji,jk) = q_i_b(ji,jk) * ht_i_b(ji) / nlay_i 490 eti (ji,jl) = eti(ji,jl) + etilayer(ji,jk) 557 491 END DO 558 492 END DO 559 560 ! total q over all layers, snow [J.m-2] 561 DO ji = kideb, kiut 562 zdes = q_s_b(ji,1) * ht_s_b(ji) / nlay_s 563 ets(ji,jl) = ets(ji,jl) + zdes 564 END DO 565 566 WRITE(numout,*) ' lim_thd_glohec ' 567 WRITE(numout,*) ' qt_i_in : ', eti(jiindex_1d,jl) / rdt_ice 568 WRITE(numout,*) ' qt_s_in : ', ets(jiindex_1d,jl) / rdt_ice 569 WRITE(numout,*) ' qt_in : ', ( eti(jiindex_1d,jl) + & 570 ets(jiindex_1d,jl) ) / rdt_ice 571 493 DO ji = kideb, kiut ! total q over all layers, snow [J.m-2] 494 ets(ji,jl) = ets(ji,jl) + q_s_b(ji,1) * ht_s_b(ji) / nlay_s 495 END DO 496 497 IF(lwp) WRITE(numout,*) ' lim_thd_glohec ' 498 IF(lwp) WRITE(numout,*) ' qt_i_in : ', eti(jiindex_1d,jl) / rdt_ice 499 IF(lwp) WRITE(numout,*) ' qt_s_in : ', ets(jiindex_1d,jl) / rdt_ice 500 IF(lwp) WRITE(numout,*) ' qt_in : ', ( eti(jiindex_1d,jl) + ets(jiindex_1d,jl) ) / rdt_ice 501 ! 572 502 END SUBROUTINE lim_thd_glohec 573 503 574 !=============================================================================== 575 576 SUBROUTINE lim_thd_con_dif(kideb,kiut,jl) 504 505 SUBROUTINE lim_thd_con_dif( kideb, kiut, jl ) 577 506 !!----------------------------------------------------------------------- 578 507 !! *** ROUTINE lim_thd_con_dif *** 579 508 !! 580 509 !! ** Purpose : Test energy conservation after heat diffusion 581 !!582 !! history :583 !! 9.9 ! 07-04 (M.Vancoppenolle) original code584 510 !!------------------------------------------------------------------- 585 !! * Local variables 586 INTEGER, INTENT(in) :: & 587 kideb, kiut, & !: bounds for the spatial loop 588 jl !: category number 589 590 REAL(wp) :: & !: ! goes to trash 591 meance, & !: mean conservation error 592 max_cons_err, & !: maximum tolerated conservation error 593 max_surf_err !: maximum tolerated surface error 594 595 INTEGER :: & 596 numce !: number of points for which conservation 597 ! is violated 598 INTEGER :: & 599 ji,jk, & !: loop indices 600 zji, zjj 511 INTEGER , INTENT(in ) :: kideb, kiut ! bounds for the spatial loop 512 INTEGER , INTENT(in ) :: jl ! category number 513 514 INTEGER :: ji, jk ! loop indices 515 INTEGER :: zji, zjj 516 INTEGER :: numce ! number of points for which conservation is violated 517 REAL(wp) :: meance ! mean conservation error 518 REAL(wp) :: max_cons_err, max_surf_err 601 519 !!--------------------------------------------------------------------- 602 520 603 max_cons_err = 1.0 604 max_surf_err = 0.001 521 max_cons_err = 1.0 ! maximum tolerated conservation error 522 max_surf_err = 0.001 ! maximum tolerated surface error 605 523 606 524 !-------------------------- … … 609 527 ! global 610 528 DO ji = kideb, kiut 611 dq_i(ji,jl) = qt_i_fin(ji,jl) - qt_i_in(ji,jl) & 612 + qt_s_fin(ji,jl) - qt_s_in(ji,jl) 529 dq_i(ji,jl) = qt_i_fin(ji,jl) - qt_i_in(ji,jl) + qt_s_fin(ji,jl) - qt_s_in(ji,jl) 613 530 END DO 614 531 ! layer by layer … … 620 537 621 538 DO ji = kideb, kiut 622 zji = MOD( npb(ji) - 1, jpi ) + 1 623 zjj = ( npb(ji) - 1 ) / jpi + 1 624 625 fatm(ji,jl) = & 626 qnsr_ice_1d(ji) + & ! atm non solar 627 (1.0-i0(ji))*qsr_ice_1d(ji) ! atm solar 628 629 sum_fluxq(ji,jl) = fc_su(ji) - fc_bo_i(ji) + qsr_ice_1d(ji)*i0(ji) & 630 - fstroc(zji,zjj,jl) 539 zji = MOD( npb(ji) - 1, jpi ) + 1 540 zjj = ( npb(ji) - 1 ) / jpi + 1 541 542 fatm(ji,jl) = qnsr_ice_1d(ji) + (1.0-i0(ji))*qsr_ice_1d(ji) 543 544 sum_fluxq(ji,jl) = fc_su(ji) - fc_bo_i(ji) + qsr_ice_1d(ji)*i0(ji) - fstroc(zji,zjj,jl) 631 545 END DO 632 546 … … 651 565 WRITE(numout,*) ' Maximum tolerated conservation error : ', max_cons_err 652 566 WRITE(numout,*) ' After lim_thd_dif, category : ', jl 653 WRITE(numout,*) ' Mean conservation error on big error points ', meance, & 654 numit 567 WRITE(numout,*) ' Mean conservation error on big error points ', meance, numit 655 568 WRITE(numout,*) ' Number of points where there is a cons err gt than c.e. : ', numce, numit 656 569 … … 751 664 752 665 END DO 753 666 ! 754 667 END SUBROUTINE lim_thd_con_dif 755 668 756 !==============================================================================757 669 758 670 SUBROUTINE lim_thd_con_dh(kideb,kiut,jl) … … 765 677 !! 9.9 ! 07-04 (M.Vancoppenolle) original code 766 678 !!----------------------------------------------------------------------- 767 !! * Local variables768 679 INTEGER, INTENT(in) :: & 769 680 kideb, kiut, & !: bounds for the spatial loop … … 882 793 883 794 ENDIF 884 885 END DO 886 795 ! 796 END DO 797 ! 887 798 END SUBROUTINE lim_thd_con_dh 888 !============================================================================== 889 890 SUBROUTINE lim_thd_enmelt( kideb,kiut)799 800 801 SUBROUTINE lim_thd_enmelt( kideb, kiut ) 891 802 !!----------------------------------------------------------------------- 892 803 !! *** ROUTINE lim_thd_enmelt *** … … 896 807 !! ** Method : Formula (Bitz and Lipscomb, 1999) 897 808 !! 898 !! history : Martin Vancoppenolle, May 2007899 809 !!------------------------------------------------------------------- 900 INTEGER, INTENT(in) :: & 901 kideb, kiut !: bounds for the spatial loop 902 903 REAL(wp) :: & !: goes to trash 904 ztmelts , & !: sea ice freezing point in K 905 zeps 906 907 INTEGER :: & 908 ji, & !: spatial loop index 909 jk !: vertical index 910 810 INTEGER, INTENT(in) :: kideb, kiut ! bounds for the spatial loop 811 !! 812 INTEGER :: ji, jk !dummy loop indices 813 REAL(wp) :: ztmelts, zeps ! temporary scalar 911 814 !!------------------------------------------------------------------- 912 zeps = 1.0e-10 913 914 ! Sea ice energy of melting 915 DO jk = 1, nlay_i 815 zeps = 1.e-10 816 ! 817 DO jk = 1, nlay_i ! Sea ice energy of melting 916 818 DO ji = kideb, kiut 917 ztmelts = - tmut * s_i_b(ji,jk) + rtt 918 q_i_b(ji,jk) = rhoic*( cpic * ( ztmelts - t_i_b(ji,jk) ) & 919 + lfus * ( 1.0 - (ztmelts-rtt)/MIN((t_i_b(ji,jk)-rtt),-zeps) ) & 920 - rcp * ( ztmelts-rtt ) ) 921 END DO !ji 922 END DO !jk 923 924 ! Snow energy of melting 925 DO jk = 1, nlay_s 819 ztmelts = - tmut * s_i_b(ji,jk) + rtt 820 q_i_b(ji,jk) = rhoic * ( cpic * ( ztmelts - t_i_b(ji,jk) ) & 821 & + lfus * ( 1.0 - (ztmelts-rtt) / MIN( t_i_b(ji,jk)-rtt, -zeps ) ) & 822 & - rcp * ( ztmelts-rtt ) ) 823 END DO 824 END DO 825 DO jk = 1, nlay_s ! Snow energy of melting 926 826 DO ji = kideb,kiut 927 827 q_s_b(ji,jk) = rhosn * ( cpic * ( rtt - t_s_b(ji,jk) ) + lfus ) 928 END DO !ji929 END DO !jk930 828 END DO 829 END DO 830 ! 931 831 END SUBROUTINE lim_thd_enmelt 932 832 933 !==============================================================================934 833 935 834 SUBROUTINE lim_thd_init … … 939 838 !! 940 839 !! ** Purpose : Physical constants and parameters linked to the ice 941 !! thermodynamics840 !! thermodynamics 942 841 !! 943 842 !! ** Method : Read the namicethd namelist and check the ice-thermo 944 !! parameter values called at the first timestep (nit000)843 !! parameter values called at the first timestep (nit000) 945 844 !! 946 845 !! ** input : Namelist namicether 947 !! 948 !! history : 949 !! 8.5 ! 03-08 (C. Ethe) original code 950 !!------------------------------------------------------------------- 951 NAMELIST/namicethd/ hmelt , hiccrit, fraz_swi, maxfrazb, vfrazb, Cfrazb, & 952 & hicmin, hiclim, amax , & 953 & sbeta , parlat, hakspl, hibspl, exld, & 954 & hakdif, hnzst , thth , parsub, alphs, betas, & 846 !!------------------------------------------------------------------- 847 NAMELIST/namicethd/ hmelt , hiccrit, fraz_swi, maxfrazb, vfrazb, Cfrazb, & 848 & hicmin, hiclim, amax , & 849 & sbeta , parlat, hakspl, hibspl, exld, & 850 & hakdif, hnzst , thth , parsub, alphs, betas, & 955 851 & kappa_i, nconv_i_thd, maxer_i_thd, thcon_i_swi 956 852 !!------------------------------------------------------------------- 957 853 958 ! Define the initial parameters 959 ! ------------------------- 960 REWIND( numnam_ice ) 854 IF(lwp) THEN 855 WRITE(numout,*) 856 WRITE(numout,*) 'lim_thd : Ice Thermodynamics' 857 WRITE(numout,*) '~~~~~~~' 858 ENDIF 859 860 REWIND( numnam_ice ) ! read Namelist numnam_ice 961 861 READ ( numnam_ice , namicethd ) 962 IF (lwp) THEN 862 863 IF(lwp) THEN ! control print 963 864 WRITE(numout,*) 964 WRITE(numout,*)'lim_thd_init : ice parameters for ice thermodynamic computation ' 965 WRITE(numout,*)'~~~~~~~~~~~~' 966 WRITE(numout,*)' maximum melting at the bottom hmelt = ', hmelt 967 WRITE(numout,*)' ice thick. for lateral accretion in NH (SH) hiccrit(1/2) = ', hiccrit 968 WRITE(numout,*)' Frazil ice thickness as a function of wind or not fraz_swi = ', fraz_swi 969 WRITE(numout,*)' Maximum proportion of frazil ice collecting at bottom maxfrazb = ', maxfrazb 970 WRITE(numout,*)' Thresold relative drift speed for collection of frazil vfrazb = ', vfrazb 971 WRITE(numout,*)' Squeezing coefficient for collection of frazil Cfrazb = ', Cfrazb 972 WRITE(numout,*)' ice thick. corr. to max. energy stored in brine pocket hicmin = ', hicmin 973 WRITE(numout,*)' minimum ice thickness hiclim = ', hiclim 974 WRITE(numout,*)' maximum lead fraction amax = ', amax 975 WRITE(numout,*)' numerical carac. of the scheme for diffusion in ice ' 976 WRITE(numout,*)' Cranck-Nicholson (=0.5), implicit (=1), explicit (=0) sbeta = ', sbeta 977 WRITE(numout,*)' percentage of energy used for lateral ablation parlat = ', parlat 978 WRITE(numout,*)' slope of distr. for Hakkinen-Mellor lateral melting hakspl = ', hakspl 979 WRITE(numout,*)' slope of distribution for Hibler lateral melting hibspl = ', hibspl 980 WRITE(numout,*)' exponent for leads-closure rate exld = ', exld 981 WRITE(numout,*)' coefficient for diffusions of ice and snow hakdif = ', hakdif 982 WRITE(numout,*)' threshold thick. for comp. of eq. thermal conductivity zhth = ', thth 983 WRITE(numout,*)' thickness of the surf. layer in temp. computation hnzst = ', hnzst 984 WRITE(numout,*)' switch for snow sublimation (=1) or not (=0) parsub = ', parsub 985 WRITE(numout,*)' coefficient for snow density when snow ice formation alphs = ', alphs 986 WRITE(numout,*)' coefficient for ice-lead partition of snowfall betas = ', betas 987 WRITE(numout,*)' extinction radiation parameter in sea ice (1.0) kappa_i = ', kappa_i 988 WRITE(numout,*)' maximal n. of iter. for heat diffusion computation nconv_i_thd = ', nconv_i_thd 989 WRITE(numout,*)' maximal err. on T for heat diffusion computation maxer_i_thd = ', maxer_i_thd 990 WRITE(numout,*)' switch for comp. of thermal conductivity in the ice thcon_i_swi = ', thcon_i_swi 991 WRITE(numout,*) 865 WRITE(numout,*)' Namelist of ice parameters for ice thermodynamic computation ' 866 WRITE(numout,*)' maximum melting at the bottom hmelt = ', hmelt 867 WRITE(numout,*)' ice thick. for lateral accretion in NH (SH) hiccrit(1/2) = ', hiccrit 868 WRITE(numout,*)' Frazil ice thickness as a function of wind or not fraz_swi = ', fraz_swi 869 WRITE(numout,*)' Maximum proportion of frazil ice collecting at bottom maxfrazb = ', maxfrazb 870 WRITE(numout,*)' Thresold relative drift speed for collection of frazil vfrazb = ', vfrazb 871 WRITE(numout,*)' Squeezing coefficient for collection of frazil Cfrazb = ', Cfrazb 872 WRITE(numout,*)' ice thick. corr. to max. energy stored in brine pocket hicmin = ', hicmin 873 WRITE(numout,*)' minimum ice thickness hiclim = ', hiclim 874 WRITE(numout,*)' maximum lead fraction amax = ', amax 875 WRITE(numout,*)' numerical carac. of the scheme for diffusion in ice ' 876 WRITE(numout,*)' Cranck-Nicholson (=0.5), implicit (=1), explicit (=0) sbeta = ', sbeta 877 WRITE(numout,*)' percentage of energy used for lateral ablation parlat = ', parlat 878 WRITE(numout,*)' slope of distr. for Hakkinen-Mellor lateral melting hakspl = ', hakspl 879 WRITE(numout,*)' slope of distribution for Hibler lateral melting hibspl = ', hibspl 880 WRITE(numout,*)' exponent for leads-closure rate exld = ', exld 881 WRITE(numout,*)' coefficient for diffusions of ice and snow hakdif = ', hakdif 882 WRITE(numout,*)' threshold thick. for comp. of eq. thermal conductivity zhth = ', thth 883 WRITE(numout,*)' thickness of the surf. layer in temp. computation hnzst = ', hnzst 884 WRITE(numout,*)' switch for snow sublimation (=1) or not (=0) parsub = ', parsub 885 WRITE(numout,*)' coefficient for snow density when snow ice formation alphs = ', alphs 886 WRITE(numout,*)' coefficient for ice-lead partition of snowfall betas = ', betas 887 WRITE(numout,*)' extinction radiation parameter in sea ice (1.0) kappa_i = ', kappa_i 888 WRITE(numout,*)' maximal n. of iter. for heat diffusion computation nconv_i_thd = ', nconv_i_thd 889 WRITE(numout,*)' maximal err. on T for heat diffusion computation maxer_i_thd = ', maxer_i_thd 890 WRITE(numout,*)' switch for comp. of thermal conductivity in the ice thcon_i_swi = ', thcon_i_swi 992 891 ENDIF 993 892 ! 994 893 rcdsn = hakdif * rcdsn 995 894 rcdic = hakdif * rcdic 996 997 895 ! 998 896 END SUBROUTINE lim_thd_init 999 897 1000 898 #else 1001 !!====================================================================== 1002 !! *** MODULE limthd *** 1003 !! No sea ice model 1004 !!====================================================================== 899 !!---------------------------------------------------------------------- 900 !! Default option NO LIM3 sea-ice model 901 !!---------------------------------------------------------------------- 1005 902 CONTAINS 1006 903 SUBROUTINE lim_thd ! Empty routine
Note: See TracChangeset
for help on using the changeset viewer.
