- Timestamp:
- 2015-06-04T16:12:19+02:00 (9 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2014/dev_r5134_UKMO4_CF_compliance/NEMOGCM/NEMO/LIM_SRC_3/limthd.F90
r5134 r5350 89 89 REAL(wp) :: zfric_u, zqld, zqfr 90 90 REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 91 REAL(wp), PARAMETER :: zfric_umin = 0._wp ! lower bound for the friction velocity (cice value=5.e-04)92 REAL(wp), PARAMETER :: zch = 0.0057_wp ! heat transfer coefficient91 REAL(wp), PARAMETER :: zfric_umin = 0._wp ! lower bound for the friction velocity (cice value=5.e-04) 92 REAL(wp), PARAMETER :: zch = 0.0057_wp ! heat transfer coefficient 93 93 ! 94 94 REAL(wp), POINTER, DIMENSION(:,:) :: zqsr, zqns 95 95 !!------------------------------------------------------------------- 96 CALL wrk_alloc( jpi, 96 CALL wrk_alloc( jpi,jpj, zqsr, zqns ) 97 97 98 98 IF( nn_timing == 1 ) CALL timing_start('limthd') … … 101 101 IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limthd', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 102 102 103 CALL lim_var_glo2eqv 103 104 !------------------------------------------------------------------------! 104 105 ! 1) Initialization of some variables ! … … 209 210 ! Net heat flux on top of ice-ocean [W.m-2] 210 211 ! ----------------------------------------- 211 ! First step here :heat flux at the ocean surface + precip212 ! Second step below : heat flux at the ice surface (after limthd_dif)212 ! heat flux at the ocean surface + precip 213 ! + heat flux at the ice surface 213 214 hfx_in(ji,jj) = hfx_in(ji,jj) & 214 215 ! heat flux above the ocean … … 216 217 ! latent heat of precip (note that precip is included in qns but not in qns_ice) 217 218 & + ( 1._wp - pfrld(ji,jj) ) * sprecip(ji,jj) * ( cpic * ( MIN( tatm_ice(ji,jj), rt0_snow ) - rt0 ) - lfus ) & 218 & + ( 1._wp - pfrld(ji,jj) ) * ( tprecip(ji,jj) - sprecip(ji,jj) ) * rcp * ( tatm_ice(ji,jj) - rt0 ) 219 & + ( 1._wp - pfrld(ji,jj) ) * ( tprecip(ji,jj) - sprecip(ji,jj) ) * rcp * ( tatm_ice(ji,jj) - rt0 ) & 220 ! heat flux above the ice 221 & + SUM( a_i_b(ji,jj,:) * ( qns_ice(ji,jj,:) + qsr_ice(ji,jj,:) ) ) 219 222 220 223 ! ----------------------------------------------------------------------------- … … 226 229 hfx_out(ji,jj) = hfx_out(ji,jj) & 227 230 ! Non solar heat flux received by the ocean 228 & + pfrld(ji,jj) * qns(ji,jj) &231 & + pfrld(ji,jj) * zqns(ji,jj) & 229 232 ! latent heat of precip (note that precip is included in qns but not in qns_ice) 230 233 & + ( pfrld(ji,jj)**rn_betas - pfrld(ji,jj) ) * sprecip(ji,jj) & … … 311 314 ! --- lateral melting if monocat --- ! 312 315 !------------------------------------! 313 IF ( ( ( nn_monocat == 1 ) .OR. ( nn_monocat == 4 ) ) .AND. ( jpl == 1 )) THEN316 IF ( ( nn_monocat == 1 .OR. nn_monocat == 4 ) .AND. jpl == 1 ) THEN 314 317 CALL lim_thd_lam( 1, nbpb ) 315 318 END IF … … 324 327 ENDIF 325 328 ! 326 END DO 329 END DO !jl 327 330 328 331 !------------------------------------------------------------------------------! … … 350 353 END DO 351 354 352 !------------------------353 ! Ice natural aging354 !------------------------355 oa_i(:,:,:) = oa_i(:,:,:) + a_i(:,:,:) * rdt_ice /rday356 357 355 !---------------------------------- 358 356 ! Change thickness to volume 359 357 !---------------------------------- 360 CALL lim_var_eqv2glo 358 v_i(:,:,:) = ht_i(:,:,:) * a_i(:,:,:) 359 v_s(:,:,:) = ht_s(:,:,:) * a_i(:,:,:) 360 smv_i(:,:,:) = sm_i(:,:,:) * v_i(:,:,:) 361 362 ! update ice age (in case a_i changed, i.e. becomes 0 or lateral melting in monocat) 363 DO jl = 1, jpl 364 DO jj = 1, jpj 365 DO ji = 1, jpi 366 rswitch = MAX( 0._wp , SIGN( 1._wp, a_i_b(ji,jj,jl) - epsi10 ) ) 367 oa_i(ji,jj,jl) = rswitch * oa_i(ji,jj,jl) * a_i(ji,jj,jl) / MAX( a_i_b(ji,jj,jl), epsi10 ) 368 END DO 369 END DO 370 END DO 361 371 362 372 CALL lim_var_zapsmall 373 363 374 !-------------------------------------------- 364 375 ! Diagnostic thermodynamic growth rates … … 399 410 ! 400 411 ! 401 CALL wrk_dealloc( jpi, jpj, zqsr, zqns )402 403 412 IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limthd', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 413 414 CALL wrk_dealloc( jpi,jpj, zqsr, zqns ) 415 404 416 !------------------------------------------------------------------------------| 405 417 ! 6) Transport of ice between thickness categories. | 406 418 !------------------------------------------------------------------------------| 419 ! Given thermodynamic growth rates, transport ice between thickness categories. 407 420 IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limitd_th_rem', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 408 421 409 ! Given thermodynamic growth rates, transport ice between thickness categories. 410 IF( jpl > 1 ) CALL lim_itd_th_rem( 1, jpl, kt ) 411 ! 412 CALL lim_var_glo2eqv ! only for info 413 CALL lim_var_agg(1) 422 IF( jpl > 1 ) CALL lim_itd_th_rem( 1, jpl, kt ) 414 423 415 424 IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limitd_th_rem', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 425 416 426 !------------------------------------------------------------------------------| 417 427 ! 7) Add frazil ice growing in leads. 418 428 !------------------------------------------------------------------------------| 419 429 IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limthd_lac', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 430 420 431 CALL lim_thd_lac 421 CALL lim_var_glo2eqv ! only for info422 432 423 ! conservation test424 433 IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limthd_lac', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 425 434 426 IF(ln_ctl) THEN ! Control print 435 ! Control print 436 IF(ln_ctl) THEN 437 CALL lim_var_glo2eqv 438 427 439 CALL prt_ctl_info(' ') 428 440 CALL prt_ctl_info(' - Cell values : ') … … 503 515 REAL(wp) :: zhi_bef ! ice thickness before thermo 504 516 REAL(wp) :: zdh_mel, zda_mel ! net melting 505 REAL(wp) :: zv ! ice volume517 REAL(wp) :: zvi, zvs ! ice/snow volumes 506 518 507 519 DO ji = kideb, kiut 508 520 zdh_mel = MIN( 0._wp, dh_i_surf(ji) + dh_i_bott(ji) + dh_snowice(ji) ) 509 IF( zdh_mel < 0._wp ) THEN 510 zv = a_i_1d(ji) * ht_i_1d(ji) 521 IF( zdh_mel < 0._wp .AND. a_i_1d(ji) > 0._wp ) THEN 522 zvi = a_i_1d(ji) * ht_i_1d(ji) 523 zvs = a_i_1d(ji) * ht_s_1d(ji) 511 524 ! lateral melting = concentration change 512 525 zhi_bef = ht_i_1d(ji) - zdh_mel 513 zda_mel = a_i_1d(ji) * zdh_mel / ( 2._wp * MAX( zhi_bef, epsi10 ) ) 514 a_i_1d(ji) = MAX( 0._wp, a_i_1d(ji) + zda_mel ) 515 ! adjust thickness 516 rswitch = MAX( 0._wp , SIGN( 1._wp , a_i_1d(ji) - epsi20 ) ) 517 ht_i_1d(ji) = rswitch * zv / MAX( a_i_1d(ji), epsi20 ) 526 rswitch = MAX( 0._wp , SIGN( 1._wp , zhi_bef - epsi20 ) ) 527 zda_mel = rswitch * a_i_1d(ji) * zdh_mel / ( 2._wp * MAX( zhi_bef, epsi20 ) ) 528 a_i_1d(ji) = MAX( epsi20, a_i_1d(ji) + zda_mel ) 529 ! adjust thickness 530 ht_i_1d(ji) = zvi / a_i_1d(ji) 531 ht_s_1d(ji) = zvs / a_i_1d(ji) 518 532 ! retrieve total concentration 519 533 at_i_1d(ji) = a_i_1d(ji) … … 601 615 CALL tab_2d_1d( nbpb, hfx_err_1d (1:nbpb), hfx_err , jpi, jpj, npb(1:nbpb) ) 602 616 CALL tab_2d_1d( nbpb, hfx_res_1d (1:nbpb), hfx_res , jpi, jpj, npb(1:nbpb) ) 617 CALL tab_2d_1d( nbpb, hfx_err_dif_1d (1:nbpb), hfx_err_dif , jpi, jpj, npb(1:nbpb) ) 603 618 CALL tab_2d_1d( nbpb, hfx_err_rem_1d (1:nbpb), hfx_err_rem , jpi, jpj, npb(1:nbpb) ) 604 619 … … 651 666 CALL tab_1d_2d( nbpb, hfx_res , npb, hfx_res_1d(1:nbpb) , jpi, jpj ) 652 667 CALL tab_1d_2d( nbpb, hfx_err_rem , npb, hfx_err_rem_1d(1:nbpb), jpi, jpj ) 668 CALL tab_1d_2d( nbpb, hfx_err_dif , npb, hfx_err_dif_1d(1:nbpb), jpi, jpj ) 653 669 ! 654 670 CALL tab_1d_2d( nbpb, qns_ice(:,:,jl), npb, qns_ice_1d(1:nbpb) , jpi, jpj) … … 674 690 INTEGER :: ios ! Local integer output status for namelist read 675 691 NAMELIST/namicethd/ rn_hnewice, ln_frazil, rn_maxfrazb, rn_vfrazb, rn_Cfrazb, & 676 & rn_himin, parsub,rn_betas, rn_kappa_i, nn_conv_dif, rn_terr_dif, nn_ice_thcon, &677 & nn_monocat 692 & rn_himin, rn_betas, rn_kappa_i, nn_conv_dif, rn_terr_dif, nn_ice_thcon, & 693 & nn_monocat, ln_it_qnsice 678 694 !!------------------------------------------------------------------- 679 695 ! … … 698 714 ENDIF 699 715 700 IF( lk_cpl .AND. parsub /= 0.0 ) CALL ctl_stop( 'In coupled mode, use parsub = 0. or send dqla' )701 716 ! 702 717 IF(lwp) THEN ! control print … … 710 725 WRITE(numout,*)' minimum ice thickness rn_himin = ', rn_himin 711 726 WRITE(numout,*)' numerical carac. of the scheme for diffusion in ice ' 712 WRITE(numout,*)' switch for snow sublimation (=1) or not (=0) parsub = ', parsub713 727 WRITE(numout,*)' coefficient for ice-lead partition of snowfall rn_betas = ', rn_betas 714 728 WRITE(numout,*)' extinction radiation parameter in sea ice rn_kappa_i = ', rn_kappa_i … … 718 732 WRITE(numout,*)' check heat conservation in the ice/snow con_i = ', con_i 719 733 WRITE(numout,*)' virtual ITD mono-category parameterizations (1) or not nn_monocat = ', nn_monocat 734 WRITE(numout,*)' iterate the surface non-solar flux (T) or not (F) ln_it_qnsice = ', ln_it_qnsice 720 735 ENDIF 721 736 !
Note: See TracChangeset
for help on using the changeset viewer.