- Timestamp:
- 2013-12-11T15:38:42+01:00 (10 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limthd_dif.F90
r4220 r4332 31 31 PUBLIC lim_thd_dif ! called by lim_thd 32 32 33 REAL(wp) :: epsi20 = 1e-20 ! constant values 34 REAL(wp) :: epsi13 = 1e-13 ! constant values 35 33 REAL(wp) :: epsi10 = 1.e-10_wp ! 36 34 !!---------------------------------------------------------------------- 37 35 !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011) … … 107 105 INTEGER, DIMENSION(kiut) :: numeqmax ! reference number of bottom equation 108 106 INTEGER, DIMENSION(kiut) :: isnow ! switch for presence (1) or absence (0) of snow 109 REAL(wp) :: zeps = 1.e-10_wp !110 107 REAL(wp) :: zg1s = 2._wp ! for the tridiagonal system 111 108 REAL(wp) :: zg1 = 2._wp ! … … 114 111 REAL(wp) :: zraext_s = 1.e+8_wp ! extinction coefficient of radiation in the snow 115 112 REAL(wp) :: zkimin = 0.10_wp ! minimum ice thermal conductivity 116 REAL(wp) :: zht_smin = 1.e-4_wp ! minimum snow depth117 113 REAL(wp) :: ztmelt_i ! ice melting temperature 118 114 REAL(wp) :: zerritmax ! current maximal error on temperature … … 312 308 IF( thcon_i_swi == 0 ) THEN ! Untersteiner (1964) formula 313 309 DO ji = kideb , kiut 314 ztcond_i(ji,0) = rcdic + zbeta*s_i_b(ji,1) / MIN(- zeps,t_i_b(ji,1)-rtt)310 ztcond_i(ji,0) = rcdic + zbeta*s_i_b(ji,1) / MIN(-epsi10,t_i_b(ji,1)-rtt) 315 311 ztcond_i(ji,0) = MAX(ztcond_i(ji,0),zkimin) 316 312 END DO … … 318 314 DO ji = kideb , kiut 319 315 ztcond_i(ji,layer) = rcdic + zbeta*( s_i_b(ji,layer) + s_i_b(ji,layer+1) ) / & 320 MIN(-2.0_wp * zeps, t_i_b(ji,layer)+t_i_b(ji,layer+1) - 2.0_wp * rtt)316 MIN(-2.0_wp * epsi10, t_i_b(ji,layer)+t_i_b(ji,layer+1) - 2.0_wp * rtt) 321 317 ztcond_i(ji,layer) = MAX(ztcond_i(ji,layer),zkimin) 322 318 END DO … … 326 322 IF( thcon_i_swi == 1 ) THEN ! Pringle et al formula included: 2.11 + 0.09 S/T - 0.011.T 327 323 DO ji = kideb , kiut 328 ztcond_i(ji,0) = rcdic + 0.090_wp * s_i_b(ji,1) / MIN( - zeps, t_i_b(ji,1)-rtt ) &324 ztcond_i(ji,0) = rcdic + 0.090_wp * s_i_b(ji,1) / MIN( -epsi10, t_i_b(ji,1)-rtt ) & 329 325 & - 0.011_wp * ( t_i_b(ji,1) - rtt ) 330 326 ztcond_i(ji,0) = MAX( ztcond_i(ji,0), zkimin ) … … 333 329 DO ji = kideb , kiut 334 330 ztcond_i(ji,layer) = rcdic + 0.090_wp * ( s_i_b(ji,layer) + s_i_b(ji,layer+1) ) & 335 & / MIN(-2.0_wp * zeps, t_i_b(ji,layer)+t_i_b(ji,layer+1) - 2.0_wp * rtt) &331 & / MIN(-2.0_wp * epsi10, t_i_b(ji,layer)+t_i_b(ji,layer+1) - 2.0_wp * rtt) & 336 332 & - 0.0055_wp* ( t_i_b(ji,layer) + t_i_b(ji,layer+1) - 2.0*rtt ) 337 333 ztcond_i(ji,layer) = MAX( ztcond_i(ji,layer), zkimin ) … … 339 335 END DO 340 336 DO ji = kideb , kiut 341 ztcond_i(ji,nlay_i) = rcdic + 0.090_wp * s_i_b(ji,nlay_i) / MIN(- zeps,t_bo_b(ji)-rtt) &337 ztcond_i(ji,nlay_i) = rcdic + 0.090_wp * s_i_b(ji,nlay_i) / MIN(-epsi10,t_bo_b(ji)-rtt) & 342 338 & - 0.011_wp * ( t_bo_b(ji) - rtt ) 343 339 ztcond_i(ji,nlay_i) = MAX( ztcond_i(ji,nlay_i), zkimin ) … … 352 348 353 349 !-- Snow kappa factors 354 zkappa_s(ji,0) = rcdsn / MAX( zeps,zh_s(ji))355 zkappa_s(ji,nlay_s) = rcdsn / MAX( zeps,zh_s(ji))350 zkappa_s(ji,0) = rcdsn / MAX(epsi10,zh_s(ji)) 351 zkappa_s(ji,nlay_s) = rcdsn / MAX(epsi10,zh_s(ji)) 356 352 END DO 357 353 … … 359 355 DO ji = kideb , kiut 360 356 zkappa_s(ji,layer) = 2.0 * rcdsn / & 361 MAX( zeps,2.0*zh_s(ji))357 MAX(epsi10,2.0*zh_s(ji)) 362 358 END DO 363 359 END DO … … 367 363 !-- Ice kappa factors 368 364 zkappa_i(ji,layer) = 2.0*ztcond_i(ji,layer)/ & 369 MAX( zeps,2.0*zh_i(ji))370 END DO 371 END DO 372 373 DO ji = kideb , kiut 374 zkappa_i(ji,0) = ztcond_i(ji,0)/MAX( zeps,zh_i(ji))375 zkappa_i(ji,nlay_i) = ztcond_i(ji,nlay_i) / MAX( zeps,zh_i(ji))365 MAX(epsi10,2.0*zh_i(ji)) 366 END DO 367 END DO 368 369 DO ji = kideb , kiut 370 zkappa_i(ji,0) = ztcond_i(ji,0)/MAX(epsi10,zh_i(ji)) 371 zkappa_i(ji,nlay_i) = ztcond_i(ji,nlay_i) / MAX(epsi10,zh_i(ji)) 376 372 !-- Interface 377 zkappa_s(ji,nlay_s) = 2.0*rcdsn*ztcond_i(ji,0)/MAX( zeps, &373 zkappa_s(ji,nlay_s) = 2.0*rcdsn*ztcond_i(ji,0)/MAX(epsi10, & 378 374 (ztcond_i(ji,0)*zh_s(ji) + rcdsn*zh_i(ji))) 379 375 zkappa_i(ji,0) = zkappa_s(ji,nlay_s)*REAL( isnow(ji) ) & … … 389 385 ztitemp(ji,layer) = t_i_b(ji,layer) 390 386 zspeche_i(ji,layer) = cpic + zgamma*s_i_b(ji,layer)/ & 391 MAX((t_i_b(ji,layer)-rtt)*(ztiold(ji,layer)-rtt), zeps)387 MAX((t_i_b(ji,layer)-rtt)*(ztiold(ji,layer)-rtt),epsi10) 392 388 zeta_i(ji,layer) = rdt_ice / MAX(rhoic*zspeche_i(ji,layer)*zh_i(ji), & 393 zeps)389 epsi10) 394 390 END DO 395 391 END DO … … 398 394 DO ji = kideb , kiut 399 395 ztstemp(ji,layer) = t_s_b(ji,layer) 400 zeta_s(ji,layer) = rdt_ice / MAX(rhosn*cpic*zh_s(ji), zeps)396 zeta_s(ji,layer) = rdt_ice / MAX(rhosn*cpic*zh_s(ji),epsi10) 401 397 END DO 402 398 END DO … … 707 703 END DO ! End of the do while iterative procedure 708 704 709 IF( ln_nicep ) THEN705 IF( ln_nicep .AND. lwp ) THEN 710 706 WRITE(numout,*) ' zerritmax : ', zerritmax 711 707 WRITE(numout,*) ' nconv : ', nconv … … 732 728 ! Heat conservation ! 733 729 !-------------------------! 734 IF( con_i ) THEN730 IF( con_i .AND. jiindex_1d > 0 ) THEN 735 731 DO ji = kideb, kiut 736 732 ! Upper snow value
Note: See TracChangeset
for help on using the changeset viewer.