- Timestamp:
- 2015-06-04T20:39:20+02:00 (9 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
branches/2015/dev_r5218_CNRS17_coupling/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90
r5065 r5357 46 46 USE sbc_ice ! Surface boundary condition: ice fields 47 47 USE lib_fortran ! to use key_nosignedzero 48 #if defined key_lim3 49 USE ice, ONLY : u_ice, v_ice, jpl, pfrld, a_i_b 50 USE limthd_dh ! for CALL lim_thd_snwblow 51 #elif defined key_lim2 52 USE ice_2, ONLY : u_ice, v_ice 53 USE par_ice_2 54 #endif 48 55 49 56 IMPLICIT NONE … … 51 58 52 59 PUBLIC sbc_blk_core ! routine called in sbcmod module 53 PUBLIC blk_ice_core ! routine called in sbc_ice_lim module 60 #if defined key_lim2 || defined key_lim3 61 PUBLIC blk_ice_core_tau ! routine called in sbc_ice_lim module 62 PUBLIC blk_ice_core_flx ! routine called in sbc_ice_lim module 63 #endif 54 64 PUBLIC blk_ice_meanqsr ! routine called in sbc_ice_lim module 55 65 PUBLIC turb_core_2z ! routine calles in sbcblk_mfs module … … 376 386 emp (:,:) = ( zevap(:,:) & ! mass flux (evap. - precip.) 377 387 & - sf(jp_prec)%fnow(:,:,1) * rn_pfac ) * tmask(:,:,1) 378 qns(:,:) = zqlw(:,:) - zqsb(:,:) - zqla(:,:) & ! Downward Non Solar flux 388 ! 389 qns(:,:) = zqlw(:,:) - zqsb(:,:) - zqla(:,:) & ! Downward Non Solar 379 390 & - sf(jp_snow)%fnow(:,:,1) * rn_pfac * lfus & ! remove latent melting heat for solid precip 380 391 & - zevap(:,:) * pst(:,:) * rcp & ! remove evap heat content at SST … … 383 394 & + sf(jp_snow)%fnow(:,:,1) * rn_pfac & ! add solid precip heat content at min(Tair,Tsnow) 384 395 & * ( MIN( sf(jp_tair)%fnow(:,:,1), rt0_snow ) - rt0 ) * cpic * tmask(:,:,1) 396 ! 397 #if defined key_lim3 398 qns_oce(:,:) = zqlw(:,:) - zqsb(:,:) - zqla(:,:) ! non solar without emp (only needed by LIM3) 399 qsr_oce(:,:) = qsr(:,:) 400 #endif 385 401 ! 386 402 CALL iom_put( "qlw_oce", zqlw ) ! output downward longwave heat over the ocean … … 406 422 407 423 408 SUBROUTINE blk_ice_core( pst , pui , pvi , palb , & 409 & p_taui, p_tauj, p_qns , p_qsr, & 410 & p_qla , p_dqns, p_dqla, & 411 & p_tpr , p_spr , & 412 & p_fr1 , p_fr2 , cd_grid, pdim ) 413 !!--------------------------------------------------------------------- 414 !! *** ROUTINE blk_ice_core *** 424 #if defined key_lim2 || defined key_lim3 425 SUBROUTINE blk_ice_core_tau 426 !!--------------------------------------------------------------------- 427 !! *** ROUTINE blk_ice_core_tau *** 415 428 !! 416 429 !! ** Purpose : provide the surface boundary condition over sea-ice 417 430 !! 418 !! ** Method : compute momentum, heat and freshwater exchanged 419 !! between atmosphere and sea-ice using CORE bulk 420 !! formulea, ice variables and read atmmospheric fields. 431 !! ** Method : compute momentum using CORE bulk 432 !! formulea, ice variables and read atmospheric fields. 421 433 !! NB: ice drag coefficient is assumed to be a constant 422 !! 423 !! caution : the net upward water flux has with mm/day unit 424 !!--------------------------------------------------------------------- 425 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: pst ! ice surface temperature (>0, =rt0 over land) [Kelvin] 426 REAL(wp), DIMENSION(:,:) , INTENT(in ) :: pui ! ice surface velocity (i- and i- components [m/s] 427 REAL(wp), DIMENSION(:,:) , INTENT(in ) :: pvi ! at I-point (B-grid) or U & V-point (C-grid) 428 REAL(wp), DIMENSION(:,:,:), INTENT(in ) :: palb ! ice albedo (all skies) [%] 429 REAL(wp), DIMENSION(:,:) , INTENT( out) :: p_taui ! i- & j-components of surface ice stress [N/m2] 430 REAL(wp), DIMENSION(:,:) , INTENT( out) :: p_tauj ! at I-point (B-grid) or U & V-point (C-grid) 431 REAL(wp), DIMENSION(:,:,:), INTENT( out) :: p_qns ! non solar heat flux over ice (T-point) [W/m2] 432 REAL(wp), DIMENSION(:,:,:), INTENT( out) :: p_qsr ! solar heat flux over ice (T-point) [W/m2] 433 REAL(wp), DIMENSION(:,:,:), INTENT( out) :: p_qla ! latent heat flux over ice (T-point) [W/m2] 434 REAL(wp), DIMENSION(:,:,:), INTENT( out) :: p_dqns ! non solar heat sensistivity (T-point) [W/m2] 435 REAL(wp), DIMENSION(:,:,:), INTENT( out) :: p_dqla ! latent heat sensistivity (T-point) [W/m2] 436 REAL(wp), DIMENSION(:,:) , INTENT( out) :: p_tpr ! total precipitation (T-point) [Kg/m2/s] 437 REAL(wp), DIMENSION(:,:) , INTENT( out) :: p_spr ! solid precipitation (T-point) [Kg/m2/s] 438 REAL(wp), DIMENSION(:,:) , INTENT( out) :: p_fr1 ! 1sr fraction of qsr penetration in ice (T-point) [%] 439 REAL(wp), DIMENSION(:,:) , INTENT( out) :: p_fr2 ! 2nd fraction of qsr penetration in ice (T-point) [%] 440 CHARACTER(len=1) , INTENT(in ) :: cd_grid ! ice grid ( C or B-grid) 441 INTEGER , INTENT(in ) :: pdim ! number of ice categories 442 !! 443 INTEGER :: ji, jj, jl ! dummy loop indices 444 INTEGER :: ijpl ! number of ice categories (size of 3rd dim of input arrays) 445 REAL(wp) :: zst2, zst3 446 REAL(wp) :: zcoef_wnorm, zcoef_wnorm2, zcoef_dqlw, zcoef_dqla, zcoef_dqsb 447 REAL(wp) :: zztmp ! temporary variable 448 REAL(wp) :: zwnorm_f, zwndi_f , zwndj_f ! relative wind module and components at F-point 449 REAL(wp) :: zwndi_t , zwndj_t ! relative wind components at T-point 450 !! 451 REAL(wp), DIMENSION(:,:) , POINTER :: z_wnds_t ! wind speed ( = | U10m - U_ice | ) at T-point 452 REAL(wp), DIMENSION(:,:,:), POINTER :: z_qlw ! long wave heat flux over ice 453 REAL(wp), DIMENSION(:,:,:), POINTER :: z_qsb ! sensible heat flux over ice 454 REAL(wp), DIMENSION(:,:,:), POINTER :: z_dqlw ! long wave heat sensitivity over ice 455 REAL(wp), DIMENSION(:,:,:), POINTER :: z_dqsb ! sensible heat sensitivity over ice 456 !!--------------------------------------------------------------------- 457 ! 458 IF( nn_timing == 1 ) CALL timing_start('blk_ice_core') 459 ! 460 CALL wrk_alloc( jpi,jpj, z_wnds_t ) 461 CALL wrk_alloc( jpi,jpj,pdim, z_qlw, z_qsb, z_dqlw, z_dqsb ) 462 463 ijpl = pdim ! number of ice categories 464 434 !!--------------------------------------------------------------------- 435 INTEGER :: ji, jj ! dummy loop indices 436 REAL(wp) :: zcoef_wnorm, zcoef_wnorm2 437 REAL(wp) :: zwnorm_f, zwndi_f , zwndj_f ! relative wind module and components at F-point 438 REAL(wp) :: zwndi_t , zwndj_t ! relative wind components at T-point 439 !!--------------------------------------------------------------------- 440 ! 441 IF( nn_timing == 1 ) CALL timing_start('blk_ice_core_tau') 442 ! 465 443 ! local scalars ( place there for vector optimisation purposes) 466 444 zcoef_wnorm = rhoa * Cice 467 445 zcoef_wnorm2 = rhoa * Cice * 0.5 468 zcoef_dqlw = 4.0 * 0.95 * Stef469 zcoef_dqla = -Ls * Cice * 11637800. * (-5897.8)470 zcoef_dqsb = rhoa * cpa * Cice471 446 472 447 !!gm brutal.... 473 z_wnds_t(:,:) = 0.e0474 p_taui (:,:) = 0.e0475 p_tauj (:,:) = 0.e0448 utau_ice (:,:) = 0._wp 449 vtau_ice (:,:) = 0._wp 450 wndm_ice (:,:) = 0._wp 476 451 !!gm end 477 452 478 #if defined key_lim3479 tatm_ice(:,:) = sf(jp_tair)%fnow(:,:,1) ! LIM3: make Tair available in sea-ice. WARNING allocated after call to ice_init480 #endif481 453 ! ----------------------------------------------------------------------------- ! 482 454 ! Wind components and module relative to the moving ocean ( U10m - U_ice ) ! 483 455 ! ----------------------------------------------------------------------------- ! 484 SELECT CASE( c d_grid)456 SELECT CASE( cp_ice_msh ) 485 457 CASE( 'I' ) ! B-grid ice dynamics : I-point (i.e. F-point with sea-ice indexation) 486 458 ! and scalar wind at T-point ( = | U10m - U_ice | ) (masked) … … 489 461 ! ... scalar wind at I-point (fld being at T-point) 490 462 zwndi_f = 0.25 * ( sf(jp_wndi)%fnow(ji-1,jj ,1) + sf(jp_wndi)%fnow(ji ,jj ,1) & 491 & + sf(jp_wndi)%fnow(ji-1,jj-1,1) + sf(jp_wndi)%fnow(ji ,jj-1,1) ) - rn_vfac * pui(ji,jj)463 & + sf(jp_wndi)%fnow(ji-1,jj-1,1) + sf(jp_wndi)%fnow(ji ,jj-1,1) ) - rn_vfac * u_ice(ji,jj) 492 464 zwndj_f = 0.25 * ( sf(jp_wndj)%fnow(ji-1,jj ,1) + sf(jp_wndj)%fnow(ji ,jj ,1) & 493 & + sf(jp_wndj)%fnow(ji-1,jj-1,1) + sf(jp_wndj)%fnow(ji ,jj-1,1) ) - rn_vfac * pvi(ji,jj)465 & + sf(jp_wndj)%fnow(ji-1,jj-1,1) + sf(jp_wndj)%fnow(ji ,jj-1,1) ) - rn_vfac * v_ice(ji,jj) 494 466 zwnorm_f = zcoef_wnorm * SQRT( zwndi_f * zwndi_f + zwndj_f * zwndj_f ) 495 467 ! ... ice stress at I-point 496 p_taui(ji,jj) = zwnorm_f * zwndi_f497 p_tauj(ji,jj) = zwnorm_f * zwndj_f468 utau_ice(ji,jj) = zwnorm_f * zwndi_f 469 vtau_ice(ji,jj) = zwnorm_f * zwndj_f 498 470 ! ... scalar wind at T-point (fld being at T-point) 499 zwndi_t = sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac * 0.25 * ( pui(ji,jj+1) + pui(ji+1,jj+1) &500 & + pui(ji,jj ) + pui(ji+1,jj ) )501 zwndj_t = sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac * 0.25 * ( pvi(ji,jj+1) + pvi(ji+1,jj+1) &502 & + pvi(ji,jj ) + pvi(ji+1,jj ) )503 z_wnds_t(ji,jj) = SQRT( zwndi_t * zwndi_t + zwndj_t * zwndj_t ) * tmask(ji,jj,1)471 zwndi_t = sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac * 0.25 * ( u_ice(ji,jj+1) + u_ice(ji+1,jj+1) & 472 & + u_ice(ji,jj ) + u_ice(ji+1,jj ) ) 473 zwndj_t = sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac * 0.25 * ( v_ice(ji,jj+1) + v_ice(ji+1,jj+1) & 474 & + v_ice(ji,jj ) + v_ice(ji+1,jj ) ) 475 wndm_ice(ji,jj) = SQRT( zwndi_t * zwndi_t + zwndj_t * zwndj_t ) * tmask(ji,jj,1) 504 476 END DO 505 477 END DO 506 CALL lbc_lnk( p_taui, 'I', -1. )507 CALL lbc_lnk( p_tauj, 'I', -1. )508 CALL lbc_lnk( z_wnds_t, 'T', 1. )478 CALL lbc_lnk( utau_ice, 'I', -1. ) 479 CALL lbc_lnk( vtau_ice, 'I', -1. ) 480 CALL lbc_lnk( wndm_ice, 'T', 1. ) 509 481 ! 510 482 CASE( 'C' ) ! C-grid ice dynamics : U & V-points (same as ocean) 511 483 DO jj = 2, jpj 512 484 DO ji = fs_2, jpi ! vect. opt. 513 zwndi_t = ( sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac * 0.5 * ( pui(ji-1,jj ) + pui(ji,jj) ) )514 zwndj_t = ( sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac * 0.5 * ( pvi(ji ,jj-1) + pvi(ji,jj) ) )515 z_wnds_t(ji,jj)= SQRT( zwndi_t * zwndi_t + zwndj_t * zwndj_t ) * tmask(ji,jj,1)485 zwndi_t = ( sf(jp_wndi)%fnow(ji,jj,1) - rn_vfac * 0.5 * ( u_ice(ji-1,jj ) + u_ice(ji,jj) ) ) 486 zwndj_t = ( sf(jp_wndj)%fnow(ji,jj,1) - rn_vfac * 0.5 * ( v_ice(ji ,jj-1) + v_ice(ji,jj) ) ) 487 wndm_ice(ji,jj) = SQRT( zwndi_t * zwndi_t + zwndj_t * zwndj_t ) * tmask(ji,jj,1) 516 488 END DO 517 489 END DO 518 490 DO jj = 2, jpjm1 519 491 DO ji = fs_2, fs_jpim1 ! vect. opt. 520 p_taui(ji,jj) = zcoef_wnorm2 * ( z_wnds_t(ji+1,jj ) + z_wnds_t(ji,jj) ) &521 & * ( 0.5 * (sf(jp_wndi)%fnow(ji+1,jj,1) + sf(jp_wndi)%fnow(ji,jj,1) ) - rn_vfac * pui(ji,jj) )522 p_tauj(ji,jj) = zcoef_wnorm2 * ( z_wnds_t(ji,jj+1 ) + z_wnds_t(ji,jj) ) &523 & * ( 0.5 * (sf(jp_wndj)%fnow(ji,jj+1,1) + sf(jp_wndj)%fnow(ji,jj,1) ) - rn_vfac * pvi(ji,jj) )492 utau_ice(ji,jj) = zcoef_wnorm2 * ( wndm_ice(ji+1,jj ) + wndm_ice(ji,jj) ) & 493 & * ( 0.5 * (sf(jp_wndi)%fnow(ji+1,jj,1) + sf(jp_wndi)%fnow(ji,jj,1) ) - rn_vfac * u_ice(ji,jj) ) 494 vtau_ice(ji,jj) = zcoef_wnorm2 * ( wndm_ice(ji,jj+1 ) + wndm_ice(ji,jj) ) & 495 & * ( 0.5 * (sf(jp_wndj)%fnow(ji,jj+1,1) + sf(jp_wndj)%fnow(ji,jj,1) ) - rn_vfac * v_ice(ji,jj) ) 524 496 END DO 525 497 END DO 526 CALL lbc_lnk( p_taui, 'U', -1. )527 CALL lbc_lnk( p_tauj, 'V', -1. )528 CALL lbc_lnk( z_wnds_t, 'T', 1. )498 CALL lbc_lnk( utau_ice, 'U', -1. ) 499 CALL lbc_lnk( vtau_ice, 'V', -1. ) 500 CALL lbc_lnk( wndm_ice, 'T', 1. ) 529 501 ! 530 502 END SELECT 503 504 IF(ln_ctl) THEN 505 CALL prt_ctl(tab2d_1=utau_ice , clinfo1=' blk_ice_core: utau_ice : ', tab2d_2=vtau_ice , clinfo2=' vtau_ice : ') 506 CALL prt_ctl(tab2d_1=wndm_ice , clinfo1=' blk_ice_core: wndm_ice : ') 507 ENDIF 508 509 IF( nn_timing == 1 ) CALL timing_stop('blk_ice_core_tau') 510 511 END SUBROUTINE blk_ice_core_tau 512 513 514 SUBROUTINE blk_ice_core_flx( ptsu, palb ) 515 !!--------------------------------------------------------------------- 516 !! *** ROUTINE blk_ice_core_flx *** 517 !! 518 !! ** Purpose : provide the surface boundary condition over sea-ice 519 !! 520 !! ** Method : compute heat and freshwater exchanged 521 !! between atmosphere and sea-ice using CORE bulk 522 !! formulea, ice variables and read atmmospheric fields. 523 !! 524 !! caution : the net upward water flux has with mm/day unit 525 !!--------------------------------------------------------------------- 526 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: ptsu ! sea ice surface temperature 527 REAL(wp), DIMENSION(:,:,:), INTENT(in) :: palb ! ice albedo (all skies) 528 !! 529 INTEGER :: ji, jj, jl ! dummy loop indices 530 REAL(wp) :: zst2, zst3 531 REAL(wp) :: zcoef_dqlw, zcoef_dqla, zcoef_dqsb 532 REAL(wp) :: zztmp, z1_lsub ! temporary variable 533 !! 534 REAL(wp), DIMENSION(:,:,:), POINTER :: z_qlw ! long wave heat flux over ice 535 REAL(wp), DIMENSION(:,:,:), POINTER :: z_qsb ! sensible heat flux over ice 536 REAL(wp), DIMENSION(:,:,:), POINTER :: z_dqlw ! long wave heat sensitivity over ice 537 REAL(wp), DIMENSION(:,:,:), POINTER :: z_dqsb ! sensible heat sensitivity over ice 538 REAL(wp), DIMENSION(:,:) , POINTER :: zevap, zsnw ! evaporation and snw distribution after wind blowing (LIM3) 539 !!--------------------------------------------------------------------- 540 ! 541 IF( nn_timing == 1 ) CALL timing_start('blk_ice_core_flx') 542 ! 543 CALL wrk_alloc( jpi,jpj,jpl, z_qlw, z_qsb, z_dqlw, z_dqsb ) 544 545 ! local scalars ( place there for vector optimisation purposes) 546 zcoef_dqlw = 4.0 * 0.95 * Stef 547 zcoef_dqla = -Ls * Cice * 11637800. * (-5897.8) 548 zcoef_dqsb = rhoa * cpa * Cice 531 549 532 550 zztmp = 1. / ( 1. - albo ) 533 551 ! ! ========================== ! 534 DO jl = 1, ijpl! Loop over ice categories !552 DO jl = 1, jpl ! Loop over ice categories ! 535 553 ! ! ========================== ! 536 554 DO jj = 1 , jpj … … 539 557 ! I Radiative FLUXES ! 540 558 ! ----------------------------! 541 zst2 = p st(ji,jj,jl) * pst(ji,jj,jl)542 zst3 = p st(ji,jj,jl) * zst2559 zst2 = ptsu(ji,jj,jl) * ptsu(ji,jj,jl) 560 zst3 = ptsu(ji,jj,jl) * zst2 543 561 ! Short Wave (sw) 544 p_qsr(ji,jj,jl) = zztmp * ( 1. - palb(ji,jj,jl) ) * qsr(ji,jj)562 qsr_ice(ji,jj,jl) = zztmp * ( 1. - palb(ji,jj,jl) ) * qsr(ji,jj) 545 563 ! Long Wave (lw) 546 z_qlw(ji,jj,jl) = 0.95 * ( sf(jp_qlw)%fnow(ji,jj,1) - Stef * p st(ji,jj,jl) * zst3 ) * tmask(ji,jj,1)564 z_qlw(ji,jj,jl) = 0.95 * ( sf(jp_qlw)%fnow(ji,jj,1) - Stef * ptsu(ji,jj,jl) * zst3 ) * tmask(ji,jj,1) 547 565 ! lw sensitivity 548 566 z_dqlw(ji,jj,jl) = zcoef_dqlw * zst3 … … 554 572 ! ... turbulent heat fluxes 555 573 ! Sensible Heat 556 z_qsb(ji,jj,jl) = rhoa * cpa * Cice * z_wnds_t(ji,jj) * ( pst(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1) )574 z_qsb(ji,jj,jl) = rhoa * cpa * Cice * wndm_ice(ji,jj) * ( ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1) ) 557 575 ! Latent Heat 558 p_qla(ji,jj,jl) = rn_efac * MAX( 0.e0, rhoa * Ls * Cice * z_wnds_t(ji,jj) &559 & * ( 11637800. * EXP( -5897.8 / p st(ji,jj,jl) ) / rhoa - sf(jp_humi)%fnow(ji,jj,1) ) )560 561 IF( p_qla(ji,jj,jl) > 0._wp ) THEN562 p_dqla(ji,jj,jl) = rn_efac * zcoef_dqla * z_wnds_t(ji,jj) / ( zst2 ) * EXP( -5897.8 / pst(ji,jj,jl) )576 qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, rhoa * Ls * Cice * wndm_ice(ji,jj) & 577 & * ( 11637800. * EXP( -5897.8 / ptsu(ji,jj,jl) ) / rhoa - sf(jp_humi)%fnow(ji,jj,1) ) ) 578 ! Latent heat sensitivity for ice (Dqla/Dt) 579 IF( qla_ice(ji,jj,jl) > 0._wp ) THEN 580 dqla_ice(ji,jj,jl) = rn_efac * zcoef_dqla * wndm_ice(ji,jj) / ( zst2 ) * EXP( -5897.8 / ptsu(ji,jj,jl) ) 563 581 ELSE 564 p_dqla(ji,jj,jl) = 0._wp582 dqla_ice(ji,jj,jl) = 0._wp 565 583 ENDIF 566 584 567 585 ! Sensible heat sensitivity (Dqsb_ice/Dtn_ice) 568 z_dqsb(ji,jj,jl) = zcoef_dqsb * z_wnds_t(ji,jj)586 z_dqsb(ji,jj,jl) = zcoef_dqsb * wndm_ice(ji,jj) 569 587 570 588 ! ----------------------------! … … 572 590 ! ----------------------------! 573 591 ! Downward Non Solar flux 574 p_qns (ji,jj,jl) = z_qlw (ji,jj,jl) - z_qsb (ji,jj,jl) - p_qla(ji,jj,jl)592 qns_ice (ji,jj,jl) = z_qlw (ji,jj,jl) - z_qsb (ji,jj,jl) - qla_ice (ji,jj,jl) 575 593 ! Total non solar heat flux sensitivity for ice 576 p_dqns(ji,jj,jl) = - ( z_dqlw(ji,jj,jl) + z_dqsb(ji,jj,jl) + p_dqla(ji,jj,jl) )594 dqns_ice(ji,jj,jl) = - ( z_dqlw(ji,jj,jl) + z_dqsb(ji,jj,jl) + dqla_ice(ji,jj,jl) ) 577 595 END DO 578 596 ! … … 581 599 END DO 582 600 ! 601 tprecip(:,:) = sf(jp_prec)%fnow(:,:,1) * rn_pfac ! total precipitation [kg/m2/s] 602 sprecip(:,:) = sf(jp_snow)%fnow(:,:,1) * rn_pfac ! solid precipitation [kg/m2/s] 603 CALL iom_put( 'snowpre', sprecip * 86400. ) ! Snow precipitation 604 CALL iom_put( 'precip' , tprecip * 86400. ) ! Total precipitation 605 606 #if defined key_lim3 607 CALL wrk_alloc( jpi,jpj, zevap, zsnw ) 608 609 ! --- evaporation --- ! 610 z1_lsub = 1._wp / Lsub 611 evap_ice (:,:,:) = qla_ice (:,:,:) * z1_lsub ! sublimation 612 devap_ice(:,:,:) = dqla_ice(:,:,:) * z1_lsub 613 zevap (:,:) = emp(:,:) + tprecip(:,:) ! evaporation over ocean 614 615 ! --- evaporation minus precipitation --- ! 616 CALL lim_thd_snwblow( pfrld, zsnw ) ! snow distribution over ice after wind blowing 617 emp_oce(:,:) = pfrld(:,:) * zevap(:,:) - ( tprecip(:,:) - sprecip(:,:) ) - sprecip(:,:) * (1._wp - zsnw ) 618 emp_ice(:,:) = SUM( a_i_b(:,:,:) * evap_ice(:,:,:), dim=3 ) - sprecip(:,:) * zsnw 619 emp_tot(:,:) = emp_oce(:,:) + emp_ice(:,:) 620 621 ! --- heat flux associated with emp --- ! 622 qemp_oce(:,:) = - pfrld(:,:) * zevap(:,:) * sst_m(:,:) * rcp & ! evap 623 & + ( tprecip(:,:) - sprecip(:,:) ) * ( sf(jp_tair)%fnow(:,:,1) - rt0 ) * rcp & ! liquid precip 624 & + sprecip(:,:) * ( 1._wp - zsnw ) * & ! solid precip 625 & ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0_snow ) - rt0 ) * cpic * tmask(:,:,1) - lfus ) 626 qemp_ice(:,:) = sprecip(:,:) * zsnw * & ! solid precip (only) 627 & ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0_snow ) - rt0 ) * cpic * tmask(:,:,1) - lfus ) 628 629 ! --- total solar and non solar fluxes --- ! 630 qns_tot(:,:) = pfrld(:,:) * qns_oce(:,:) + SUM( a_i_b(:,:,:) * qns_ice(:,:,:), dim=3 ) + qemp_ice(:,:) + qemp_oce(:,:) 631 qsr_tot(:,:) = pfrld(:,:) * qsr_oce(:,:) + SUM( a_i_b(:,:,:) * qsr_ice(:,:,:), dim=3 ) 632 633 ! --- heat content of precip over ice in J/m3 (to be used in 1D-thermo) --- ! 634 qprec_ice(:,:) = rhosn * ( ( MIN( sf(jp_tair)%fnow(:,:,1), rt0_snow ) - rt0 ) * cpic * tmask(:,:,1) - lfus ) 635 636 CALL wrk_dealloc( jpi,jpj, zevap, zsnw ) 637 #endif 638 583 639 !-------------------------------------------------------------------- 584 640 ! FRACTIONs of net shortwave radiation which is not absorbed in the … … 586 642 ! ( Maykut and Untersteiner, 1971 ; Ebert and Curry, 1993 ) 587 643 ! 588 p_fr1(:,:) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ) 589 p_fr2(:,:) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice ) 590 ! 591 p_tpr(:,:) = sf(jp_prec)%fnow(:,:,1) * rn_pfac ! total precipitation [kg/m2/s] 592 p_spr(:,:) = sf(jp_snow)%fnow(:,:,1) * rn_pfac ! solid precipitation [kg/m2/s] 593 CALL iom_put( 'snowpre', p_spr * 86400. ) ! Snow precipitation 594 CALL iom_put( 'precip' , p_tpr * 86400. ) ! Total precipitation 644 fr1_i0(:,:) = ( 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice ) 645 fr2_i0(:,:) = ( 0.82 * ( 1.0 - cldf_ice ) + 0.65 * cldf_ice ) 646 ! 595 647 ! 596 648 IF(ln_ctl) THEN 597 CALL prt_ctl(tab3d_1=p_qla , clinfo1=' blk_ice_core: p_qla : ', tab3d_2=z_qsb , clinfo2=' z_qsb : ', kdim=ijpl) 598 CALL prt_ctl(tab3d_1=z_qlw , clinfo1=' blk_ice_core: z_qlw : ', tab3d_2=p_dqla , clinfo2=' p_dqla : ', kdim=ijpl) 599 CALL prt_ctl(tab3d_1=z_dqsb , clinfo1=' blk_ice_core: z_dqsb : ', tab3d_2=z_dqlw , clinfo2=' z_dqlw : ', kdim=ijpl) 600 CALL prt_ctl(tab3d_1=p_dqns , clinfo1=' blk_ice_core: p_dqns : ', tab3d_2=p_qsr , clinfo2=' p_qsr : ', kdim=ijpl) 601 CALL prt_ctl(tab3d_1=pst , clinfo1=' blk_ice_core: pst : ', tab3d_2=p_qns , clinfo2=' p_qns : ', kdim=ijpl) 602 CALL prt_ctl(tab2d_1=p_tpr , clinfo1=' blk_ice_core: p_tpr : ', tab2d_2=p_spr , clinfo2=' p_spr : ') 603 CALL prt_ctl(tab2d_1=p_taui , clinfo1=' blk_ice_core: p_taui : ', tab2d_2=p_tauj , clinfo2=' p_tauj : ') 604 CALL prt_ctl(tab2d_1=z_wnds_t, clinfo1=' blk_ice_core: z_wnds_t : ') 649 CALL prt_ctl(tab3d_1=qla_ice , clinfo1=' blk_ice_core: qla_ice : ', tab3d_2=z_qsb , clinfo2=' z_qsb : ', kdim=jpl) 650 CALL prt_ctl(tab3d_1=z_qlw , clinfo1=' blk_ice_core: z_qlw : ', tab3d_2=dqla_ice, clinfo2=' dqla_ice : ', kdim=jpl) 651 CALL prt_ctl(tab3d_1=z_dqsb , clinfo1=' blk_ice_core: z_dqsb : ', tab3d_2=z_dqlw , clinfo2=' z_dqlw : ', kdim=jpl) 652 CALL prt_ctl(tab3d_1=dqns_ice, clinfo1=' blk_ice_core: dqns_ice : ', tab3d_2=qsr_ice , clinfo2=' qsr_ice : ', kdim=jpl) 653 CALL prt_ctl(tab3d_1=ptsu , clinfo1=' blk_ice_core: ptsu : ', tab3d_2=qns_ice , clinfo2=' qns_ice : ', kdim=jpl) 654 CALL prt_ctl(tab2d_1=tprecip , clinfo1=' blk_ice_core: tprecip : ', tab2d_2=sprecip , clinfo2=' sprecip : ') 605 655 ENDIF 606 656 607 CALL wrk_dealloc( jpi,jpj, z_wnds_t)608 CALL wrk_dealloc( jpi,jpj, pdim, z_qlw, z_qsb, z_dqlw, z_dqsb )609 !610 IF( nn_timing == 1 ) CALL timing_stop('blk_ice_core')611 !612 END SUBROUTINE blk_ice_core 657 CALL wrk_dealloc( jpi,jpj,jpl, z_qlw, z_qsb, z_dqlw, z_dqsb ) 658 ! 659 IF( nn_timing == 1 ) CALL timing_stop('blk_ice_core_flx') 660 661 END SUBROUTINE blk_ice_core_flx 662 #endif 613 663 614 664
Note: See TracChangeset
for help on using the changeset viewer.