Changeset 11111 for NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk/src/OCE/SBC/sbcblk.F90
- Timestamp:
- 2019-06-14T15:55:28+02:00 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk/src/OCE/SBC/sbcblk.F90
r10535 r11111 15 15 !! 3.7 ! 2014-06 (L. Brodeau) simplification and optimization of CORE bulk 16 16 !! 4.0 ! 2016-06 (L. Brodeau) sbcblk_core becomes sbcblk and is not restricted to the CORE algorithm anymore 17 !! ! ==> based on AeroBulk (http ://aerobulk.sourceforge.net/)17 !! ! ==> based on AeroBulk (https://github.com/brodeau/aerobulk/) 18 18 !! 4.0 ! 2016-10 (G. Madec) introduce a sbc_blk_init routine 19 19 !! 4.0 ! 2016-10 (M. Vancoppenolle) Introduce conduction flux emulator (M. Vancoppenolle) … … 24 24 !! sbc_blk : bulk formulation as ocean surface boundary condition 25 25 !! blk_oce : computes momentum, heat and freshwater fluxes over ocean 26 !! rho_air : density of (moist) air (depends on T_air, q_air and SLP27 !! cp_air : specific heat of (moist) air (depends spec. hum. q_air)28 !! q_sat : saturation humidity as a function of SLP and temperature29 !! L_vap : latent heat of vaporization of water as a function of temperature30 26 !! sea-ice case only : 31 27 !! blk_ice_tau : provide the air-ice stress … … 60 56 USE prtctl ! Print control 61 57 58 USE sbcblk_phymbl !LB: all thermodynamics functions, rho_air, q_sat, etc... 59 60 62 61 IMPLICIT NONE 63 62 PRIVATE … … 70 69 PUBLIC blk_ice_qcn ! routine called in icesbc 71 70 #endif 72 73 !!Lolo: should ultimately be moved in the module with all physical constants ?74 !!gm : In principle, yes.75 REAL(wp), PARAMETER :: Cp_dry = 1005.0 !: Specic heat of dry air, constant pressure [J/K/kg]76 REAL(wp), PARAMETER :: Cp_vap = 1860.0 !: Specic heat of water vapor, constant pressure [J/K/kg]77 REAL(wp), PARAMETER :: R_dry = 287.05_wp !: Specific gas constant for dry air [J/K/kg]78 REAL(wp), PARAMETER :: R_vap = 461.495_wp !: Specific gas constant for water vapor [J/K/kg]79 REAL(wp), PARAMETER :: reps0 = R_dry/R_vap !: ratio of gas constant for dry air and water vapor => ~ 0.62280 REAL(wp), PARAMETER :: rctv0 = R_vap/R_dry !: for virtual temperature (== (1-eps)/eps) => ~ 0.60881 71 82 72 INTEGER , PARAMETER :: jpfld =10 ! maximum number of files to read … … 94 84 TYPE(FLD), ALLOCATABLE, DIMENSION(:) :: sf ! structure of input fields (file informations, fields read) 95 85 96 ! !!! Bulk parameters97 REAL(wp), PARAMETER :: cpa = 1000.5 ! specific heat of air (only used for ice fluxes now...)98 REAL(wp), PARAMETER :: Ls = 2.839e6 ! latent heat of sublimation99 REAL(wp), PARAMETER :: Stef = 5.67e-8 ! Stefan Boltzmann constant100 REAL(wp), PARAMETER :: Cd_ice = 1.4e-3 ! transfer coefficient over ice101 REAL(wp), PARAMETER :: albo = 0.066 ! ocean albedo assumed to be constant102 !103 86 ! !!* Namelist namsbc_blk : bulk parameters 104 87 LOGICAL :: ln_NCAR ! "NCAR" algorithm (Large and Yeager 2008) … … 124 107 REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: cdn_oce, chn_oce, cen_oce ! needed by Lupkes 2015 bulk scheme 125 108 109 !LB: 110 LOGICAL :: ln_skin ! use the cool-skin / warm-layer parameterization (only available in ECMWF and COARE algorithms) !LB 111 !LB: => sbc_oce.F90: REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: tsk ! sea-surface skin temperature out of the cool-skin/warm-layer parameterization [Celsius] 112 !LB. 113 126 114 INTEGER :: nblk ! choice of the bulk algorithm 127 115 ! ! associated indices: … … 150 138 IF( sbc_blk_alloc /= 0 ) CALL ctl_stop( 'STOP', 'sbc_blk_alloc: failed to allocate arrays' ) 151 139 END FUNCTION sbc_blk_alloc 140 141 !LB: 142 INTEGER FUNCTION sbc_blk_cswl_alloc() 143 !!------------------------------------------------------------------- 144 !! *** ROUTINE sbc_blk_cswl_alloc *** 145 !!------------------------------------------------------------------- 146 !PRINT *, '*** LB: allocating tsk!' 147 ALLOCATE( tsk(jpi,jpj), STAT=sbc_blk_cswl_alloc ) 148 !PRINT *, '*** LB: done!' 149 CALL mpp_sum ( 'sbcblk', sbc_blk_cswl_alloc ) 150 IF( sbc_blk_cswl_alloc /= 0 ) CALL ctl_stop( 'STOP', 'sbc_blk_cswl_alloc: failed to allocate arrays' ) 151 END FUNCTION sbc_blk_cswl_alloc 152 !LB. 152 153 153 154 … … 173 174 & ln_NCAR, ln_COARE_3p0, ln_COARE_3p5, ln_ECMWF, & ! bulk algorithm 174 175 & cn_dir , ln_taudif, rn_zqt, rn_zu, & 175 & rn_pfac, rn_efac, rn_vfac, ln_Cd_L12, ln_Cd_L15 176 & rn_pfac, rn_efac, rn_vfac, ln_Cd_L12, ln_Cd_L15, & 177 & ln_skin ! cool-skin / warm-layer param. !LB 176 178 !!--------------------------------------------------------------------- 177 179 ! … … 198 200 IF( ln_ECMWF ) THEN ; nblk = np_ECMWF ; ioptio = ioptio + 1 ; ENDIF 199 201 ! 202 !LB: 203 ! !** initialization of the cool-skin / warm-layer parametrization 204 IF( ln_skin ) THEN 205 IF ( ln_NCAR ) CALL ctl_stop( 'sbc_blk_init: Cool-skin/warm-layer param. not compatible with NCAR algorithm!' ) 206 ! ! allocate array(s) for cool-skin/warm-layer param. 207 IF( sbc_blk_cswl_alloc() /= 0 ) CALL ctl_stop( 'STOP', 'sbc_blk : unable to allocate standard arrays' ) 208 END IF 209 !LB. 210 200 211 IF( ioptio /= 1 ) CALL ctl_stop( 'sbc_blk_init: Choose one and only one bulk algorithm' ) 201 212 ! … … 279 290 END SELECT 280 291 ! 292 WRITE(numout,*) 293 WRITE(numout,*) ' use cool-skin/warm-layer parameterization (SSST) ln_skin = ', ln_skin !LB 294 ! 281 295 ENDIF 282 296 ! … … 413 427 414 428 ! ----------------------------------------------------------------------------- ! 415 ! I Radiative FLUXES!429 ! I Solar FLUX ! 416 430 ! ----------------------------------------------------------------------------- ! 417 431 … … 422 436 ENDIF 423 437 424 zqlw(:,:) = ( sf(jp_qlw)%fnow(:,:,1) - Stef * zst(:,:)*zst(:,:)*zst(:,:)*zst(:,:) ) * tmask(:,:,1) ! Long Wave425 438 426 439 ! ----------------------------------------------------------------------------- ! … … 429 442 430 443 ! ... specific humidity at SST and IST tmask( 431 zsq(:,:) = 0.98* q_sat( zst(:,:), sf(jp_slp)%fnow(:,:,1) )444 zsq(:,:) = rdct_qsat_salt * q_sat( zst(:,:), sf(jp_slp)%fnow(:,:,1) ) 432 445 !! 433 446 !! Estimate of potential temperature at z=rn_zqt, based on adiabatic lapse-rate … … 444 457 CASE( np_COARE_3p5 ) ; CALL turb_coare3p5( rn_zqt, rn_zu, zst, ztpot, zsq, sf(jp_humi)%fnow, wndm, & ! COARE v3.5 445 458 & Cd_atm, Ch_atm, Ce_atm, t_zu, q_zu, zU_zu, cdn_oce, chn_oce, cen_oce ) 446 CASE( np_ECMWF ) ; CALL turb_ecmwf ( rn_zqt, rn_zu, zst, ztpot, zsq, sf(jp_humi)%fnow, wndm, & ! ECMWF 447 & Cd_atm, Ch_atm, Ce_atm, t_zu, q_zu, zU_zu, cdn_oce, chn_oce, cen_oce ) 459 460 !LB: Skin!!! 461 CASE( np_ECMWF ) 462 IF ( ln_skin ) THEN 463 !IF (lwp) PRINT *, ' *** LB(sbcblk.F90) => calling "turb_ecmwf" WITH CSWL optional arrays!!!' 464 !IF (lwp) PRINT *, ' *** LB(sbcblk.F90) => BEFORE ZST(40:50,30) =', zst(40:50,30) 465 CALL turb_ecmwf ( rn_zqt, rn_zu, zst, ztpot, zsq, sf(jp_humi)%fnow, wndm, & ! ECMWF 466 & Cd_atm, Ch_atm, Ce_atm, t_zu, q_zu, zU_zu, cdn_oce, chn_oce, cen_oce, & 467 & Qsw=qsr(:,:), rad_lw=sf(jp_qlw)%fnow(:,:,1), slp=sf(jp_slp)%fnow(:,:,1)) 468 !LB: "zst" and "zsq" have been updated with skin temp. !!! (from bulk SST)... 469 zst(:,:) = zst(:,:)*tmask(:,:,1) 470 zsq(:,:) = zsq(:,:)*tmask(:,:,1) 471 !IF (lwp) PRINT *, ' *** LB(sbcblk.F90) => AFTER ZST(40:50,30) =', zst(40:50,30) 472 ELSE 473 IF (lwp) PRINT *, ' *** LB(sbcblk.F90) => calling "turb_ecmwf" WITHOUT CSWL optional arrays!!!' 474 CALL turb_ecmwf ( rn_zqt, rn_zu, zst, ztpot, zsq, sf(jp_humi)%fnow, wndm, & ! ECMWF 475 & Cd_atm, Ch_atm, Ce_atm, t_zu, q_zu, zU_zu, cdn_oce, chn_oce, cen_oce ) 476 END IF 477 !LB. 478 448 479 CASE DEFAULT 449 480 CALL ctl_stop( 'STOP', 'sbc_oce: non-existing bulk formula selected' ) 450 481 END SELECT 482 483 484 !LB: cool-skin/warm-layer: 485 IF ( ln_skin ) tsk(:,:) = zst(:,:) 486 451 487 452 488 ! ! Compute true air density : … … 507 543 508 544 545 ! ----------------------------------------------------------------------------- ! 546 ! III Net longwave radiative FLUX ! 547 ! ----------------------------------------------------------------------------- ! 548 549 !! LB: now after Turbulent fluxes because must use the skin temperature rather that the SST ! (zst is skin temperature if ln_skin==.TRUE.) 550 zqlw(:,:) = ( sf(jp_qlw)%fnow(:,:,1) - emiss_w * stefan * zst(:,:)*zst(:,:)*zst(:,:)*zst(:,:) ) * tmask(:,:,1) ! Long Wave 551 552 509 553 IF(ln_ctl) THEN 510 554 CALL prt_ctl( tab2d_1=zqla , clinfo1=' blk_oce: zqla : ', tab2d_2=Ce_atm , clinfo2=' Ce_oce : ' ) … … 519 563 520 564 ! ----------------------------------------------------------------------------- ! 521 ! I IITotal FLUXES !565 ! IV Total FLUXES ! 522 566 ! ----------------------------------------------------------------------------- ! 523 567 ! … … 527 571 qns(:,:) = zqlw(:,:) - zqsb(:,:) - zqla(:,:) & ! Downward Non Solar 528 572 & - sf(jp_snow)%fnow(:,:,1) * rn_pfac * rLfus & ! remove latent melting heat for solid precip 529 & - zevap(:,:) * pst(:,:) * rcp & ! remove evap heat content at SST 573 & - zevap(:,:) * pst(:,:) * rcp & ! remove evap heat content at SST !LB??? pst is Celsius !? 530 574 & + ( sf(jp_prec)%fnow(:,:,1) - sf(jp_snow)%fnow(:,:,1) ) * rn_pfac & ! add liquid precip heat content at Tair 531 575 & * ( sf(jp_tair)%fnow(:,:,1) - rt0 ) * rcp & … … 539 583 #endif 540 584 ! 541 IF ( nn_ice == 0 ) THEN 585 !!#LB: NO WHY???? IF ( nn_ice == 0 ) THEN 586 CALL iom_put( "rho_air" , zrhoa ) ! output air density (kg/m^3) !#LB 542 587 CALL iom_put( "qlw_oce" , zqlw ) ! output downward longwave heat over the ocean 543 588 CALL iom_put( "qsb_oce" , - zqsb ) ! output downward sensible heat over the ocean … … 551 596 CALL iom_put( 'snowpre', sprecip ) ! Snow 552 597 CALL iom_put( 'precip' , tprecip ) ! Total precipitation 553 ENDIF 598 IF ( ln_skin ) THEN 599 CALL iom_put( "t_skin" , (zst - rt0) * tmask(:,:,1) ) ! T_skin in Celsius 600 CALL iom_put( "dt_skin" , (zst - pst - rt0) * tmask(:,:,1) ) ! T_skin - SST temperature difference... 601 END IF 602 !!#LB. ENDIF 554 603 ! 555 604 IF(ln_ctl) THEN … … 564 613 565 614 566 567 FUNCTION rho_air( ptak, pqa, pslp )568 !!-------------------------------------------------------------------------------569 !! *** FUNCTION rho_air ***570 !!571 !! ** Purpose : compute density of (moist) air using the eq. of state of the atmosphere572 !!573 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk)574 !!-------------------------------------------------------------------------------575 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: ptak ! air temperature [K]576 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pqa ! air specific humidity [kg/kg]577 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pslp ! pressure in [Pa]578 REAL(wp), DIMENSION(jpi,jpj) :: rho_air ! density of moist air [kg/m^3]579 !!-------------------------------------------------------------------------------580 !581 rho_air = pslp / ( R_dry*ptak * ( 1._wp + rctv0*pqa ) )582 !583 END FUNCTION rho_air584 585 586 FUNCTION cp_air( pqa )587 !!-------------------------------------------------------------------------------588 !! *** FUNCTION cp_air ***589 !!590 !! ** Purpose : Compute specific heat (Cp) of moist air591 !!592 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk)593 !!-------------------------------------------------------------------------------594 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pqa ! air specific humidity [kg/kg]595 REAL(wp), DIMENSION(jpi,jpj) :: cp_air ! specific heat of moist air [J/K/kg]596 !!-------------------------------------------------------------------------------597 !598 Cp_air = Cp_dry + Cp_vap * pqa599 !600 END FUNCTION cp_air601 602 603 FUNCTION q_sat( ptak, pslp )604 !!----------------------------------------------------------------------------------605 !! *** FUNCTION q_sat ***606 !!607 !! ** Purpose : Specific humidity at saturation in [kg/kg]608 !! Based on accurate estimate of "e_sat"609 !! aka saturation water vapor (Goff, 1957)610 !!611 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk)612 !!----------------------------------------------------------------------------------613 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: ptak ! air temperature [K]614 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pslp ! sea level atmospheric pressure [Pa]615 REAL(wp), DIMENSION(jpi,jpj) :: q_sat ! Specific humidity at saturation [kg/kg]616 !617 INTEGER :: ji, jj ! dummy loop indices618 REAL(wp) :: ze_sat, ztmp ! local scalar619 !!----------------------------------------------------------------------------------620 !621 DO jj = 1, jpj622 DO ji = 1, jpi623 !624 ztmp = rt0 / ptak(ji,jj)625 !626 ! Vapour pressure at saturation [hPa] : WMO, (Goff, 1957)627 ze_sat = 10.**( 10.79574*(1. - ztmp) - 5.028*LOG10(ptak(ji,jj)/rt0) &628 & + 1.50475*10.**(-4)*(1. - 10.**(-8.2969*(ptak(ji,jj)/rt0 - 1.)) ) &629 & + 0.42873*10.**(-3)*(10.**(4.76955*(1. - ztmp)) - 1.) + 0.78614 )630 !631 q_sat(ji,jj) = reps0 * ze_sat/( 0.01_wp*pslp(ji,jj) - (1._wp - reps0)*ze_sat ) ! 0.01 because SLP is in [Pa]632 !633 END DO634 END DO635 !636 END FUNCTION q_sat637 638 639 FUNCTION gamma_moist( ptak, pqa )640 !!----------------------------------------------------------------------------------641 !! *** FUNCTION gamma_moist ***642 !!643 !! ** Purpose : Compute the moist adiabatic lapse-rate.644 !! => http://glossary.ametsoc.org/wiki/Moist-adiabatic_lapse_rate645 !! => http://www.geog.ucsb.edu/~joel/g266_s10/lecture_notes/chapt03/oh10_3_01/oh10_3_01.html646 !!647 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk)648 !!----------------------------------------------------------------------------------649 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: ptak ! air temperature [K]650 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pqa ! specific humidity [kg/kg]651 REAL(wp), DIMENSION(jpi,jpj) :: gamma_moist ! moist adiabatic lapse-rate652 !653 INTEGER :: ji, jj ! dummy loop indices654 REAL(wp) :: zrv, ziRT ! local scalar655 !!----------------------------------------------------------------------------------656 !657 DO jj = 1, jpj658 DO ji = 1, jpi659 zrv = pqa(ji,jj) / (1. - pqa(ji,jj))660 ziRT = 1. / (R_dry*ptak(ji,jj)) ! 1/RT661 gamma_moist(ji,jj) = grav * ( 1. + rLevap*zrv*ziRT ) / ( Cp_dry + rLevap*rLevap*zrv*reps0*ziRT/ptak(ji,jj) )662 END DO663 END DO664 !665 END FUNCTION gamma_moist666 667 668 FUNCTION L_vap( psst )669 !!---------------------------------------------------------------------------------670 !! *** FUNCTION L_vap ***671 !!672 !! ** Purpose : Compute the latent heat of vaporization of water from temperature673 !!674 !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://sourceforge.net/p/aerobulk)675 !!----------------------------------------------------------------------------------676 REAL(wp), DIMENSION(jpi,jpj) :: L_vap ! latent heat of vaporization [J/kg]677 REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: psst ! water temperature [K]678 !!----------------------------------------------------------------------------------679 !680 L_vap = ( 2.501 - 0.00237 * ( psst(:,:) - rt0) ) * 1.e6681 !682 END FUNCTION L_vap683 615 684 616 #if defined key_si3 … … 710 642 ! 711 643 ! set transfer coefficients to default sea-ice values 712 Cd_atm(:,:) = Cd_ice713 Ch_atm(:,:) = Cd_ice714 Ce_atm(:,:) = Cd_ice644 Cd_atm(:,:) = rCd_ice 645 Ch_atm(:,:) = rCd_ice 646 Ce_atm(:,:) = rCd_ice 715 647 716 648 wndm_ice(:,:) = 0._wp !!gm brutal.... … … 737 669 ENDIF 738 670 739 !! CALL iom_put( " Cd_ice", Cd_atm) ! output value of pure ice-atm. transfer coef.671 !! CALL iom_put( "rCd_ice", Cd_atm) ! output value of pure ice-atm. transfer coef. 740 672 !! CALL iom_put( "Ch_ice", Ch_atm) ! output value of pure ice-atm. transfer coef. 741 673 … … 792 724 REAL(wp) :: zst3 ! local variable 793 725 REAL(wp) :: zcoef_dqlw, zcoef_dqla ! - - 794 REAL(wp) :: zztmp, z1_rLsub 726 REAL(wp) :: zztmp, z1_rLsub ! - - 795 727 REAL(wp) :: zfr1, zfr2 ! local variables 796 728 REAL(wp), DIMENSION(jpi,jpj,jpl) :: z1_st ! inverse of surface temperature … … 803 735 !!--------------------------------------------------------------------- 804 736 ! 805 zcoef_dqlw = 4. 0 * 0.95 * Stef! local scalars806 zcoef_dqla = - Ls * 11637800. * (-5897.8)737 zcoef_dqlw = 4._wp * 0.95_wp * stefan ! local scalars 738 zcoef_dqla = -rLsub * 11637800._wp * (-5897.8_wp) !LB: BAD! 807 739 ! 808 740 zrhoa(:,:) = rho_air( sf(jp_tair)%fnow(:,:,1), sf(jp_humi)%fnow(:,:,1), sf(jp_slp)%fnow(:,:,1) ) … … 824 756 qsr_ice(ji,jj,jl) = zztmp * ( 1. - palb(ji,jj,jl) ) * qsr(ji,jj) 825 757 ! Long Wave (lw) 826 z_qlw(ji,jj,jl) = 0.95 * ( sf(jp_qlw)%fnow(ji,jj,1) - Stef* ptsu(ji,jj,jl) * zst3 ) * tmask(ji,jj,1)758 z_qlw(ji,jj,jl) = 0.95 * ( sf(jp_qlw)%fnow(ji,jj,1) - stefan * ptsu(ji,jj,jl) * zst3 ) * tmask(ji,jj,1) 827 759 ! lw sensitivity 828 760 z_dqlw(ji,jj,jl) = zcoef_dqlw * zst3 … … 834 766 ! ... turbulent heat fluxes with Ch_atm recalculated in blk_ice_tau 835 767 ! Sensible Heat 836 z_qsb(ji,jj,jl) = zrhoa(ji,jj) * cpa* Ch_atm(ji,jj) * wndm_ice(ji,jj) * (ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1))768 z_qsb(ji,jj,jl) = zrhoa(ji,jj) * rCp_air * Ch_atm(ji,jj) * wndm_ice(ji,jj) * (ptsu(ji,jj,jl) - sf(jp_tair)%fnow(ji,jj,1)) 837 769 ! Latent Heat 838 qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, zrhoa(ji,jj) * Ls* Ch_atm(ji,jj) * wndm_ice(ji,jj) * &770 qla_ice(ji,jj,jl) = rn_efac * MAX( 0.e0, zrhoa(ji,jj) * rLsub * Ch_atm(ji,jj) * wndm_ice(ji,jj) * & 839 771 & ( 11637800. * EXP( -5897.8 * z1_st(ji,jj,jl) ) / zrhoa(ji,jj) - sf(jp_humi)%fnow(ji,jj,1) ) ) 840 772 ! Latent heat sensitivity for ice (Dqla/Dt) … … 847 779 848 780 ! Sensible heat sensitivity (Dqsb_ice/Dtn_ice) 849 z_dqsb(ji,jj,jl) = zrhoa(ji,jj) * cpa* Ch_atm(ji,jj) * wndm_ice(ji,jj)781 z_dqsb(ji,jj,jl) = zrhoa(ji,jj) * rCp_air * Ch_atm(ji,jj) * wndm_ice(ji,jj) 850 782 851 783 ! ----------------------------! … … 1064 996 ! generic drag over a cell partly covered by ice 1065 997 !!Cd(:,:) = Cd_oce(:,:) * ( 1._wp - at_i_b(:,:) ) + & ! pure ocean drag 1066 !! & Cd_ice * at_i_b(:,:) + & ! pure ice drag998 !! & rCd_ice * at_i_b(:,:) + & ! pure ice drag 1067 999 !! & zCe * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**zmu ! change due to sea-ice morphology 1068 1000 1069 1001 ! ice-atm drag 1070 Cd(:,:) = Cd_ice + & ! pure ice drag1002 Cd(:,:) = rCd_ice + & ! pure ice drag 1071 1003 & zCe * ( 1._wp - at_i_b(:,:) )**zcoef * at_i_b(:,:)**(zmu-1._wp) ! change due to sea-ice morphology 1072 1004 … … 1141 1073 ! Atmospheric and Surface Variables 1142 1074 zst(:,:) = sst_m(:,:) + rt0 ! convert SST from Celcius to Kelvin 1143 zqo_sat(:,:) = 0.98_wp* q_sat( zst(:,:) , sf(jp_slp)%fnow(:,:,1) ) ! saturation humidity over ocean [kg/kg]1144 zqi_sat(:,:) = 0.98_wp * q_sat( ztm_su(:,:), sf(jp_slp)%fnow(:,:,1) ) ! saturation humidity over ice [kg/kg]1075 zqo_sat(:,:) = rdct_qsat_salt * q_sat( zst(:,:) , sf(jp_slp)%fnow(:,:,1) ) ! saturation humidity over ocean [kg/kg] 1076 zqi_sat(:,:) = q_sat( ztm_su(:,:), sf(jp_slp)%fnow(:,:,1) ) ! saturation humidity over ice [kg/kg] !LB: no 0.98 !!(rdct_qsat_salt) 1145 1077 ! 1146 1078 DO jj = 2, jpjm1 ! reduced loop is necessary for reproducibility
Note: See TracChangeset
for help on using the changeset viewer.