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