Ignore:
Timestamp:
2019-07-18T11:35:28+02:00 (14 months ago)
Author:
laurent
Message:

LB: gooodbye "coare3p5", hello "coare3p6" + cleaner air humidity handling + support for Relative humidity.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r11085_ASINTER-05_Brodeau_Advanced_Bulk/src/OCE/SBC/sbcblk_phy.F90

    r11266 r11284  
    1717   !!   Ri_bulk       : bulk Richardson number aka BRN 
    1818   !!   q_sat         : saturation humidity as a function of SLP and temperature 
    19  
     19   !!   q_air_rh      : specific humidity as a function of RH, t_air and SLP 
     20    
    2021   USE dom_oce        ! ocean space and time domain 
    2122   USE phycst         ! physical constants 
     
    3738   PUBLIC Ri_bulk 
    3839   PUBLIC q_sat 
     40   PUBLIC q_air_rh 
    3941 
    4042   !!---------------------------------------------------------------------- 
     
    115117      !!                           ***  FUNCTION L_vap  *** 
    116118      !! 
    117       !! ** Purpose : Compute the latent heat of vaporization of water from temperature 
     119      !! ** Purpose : Compute the latent heat of vaporization of water out of temperature 
    118120      !! 
    119121      !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://github.com/brodeau/aerobulk/) 
     
    273275 
    274276 
     277   FUNCTION e_sat_sclr( ptak, pslp ) 
     278      !!---------------------------------------------------------------------------------- 
     279      !!                   ***  FUNCTION e_sat_sclr  *** 
     280      !!                  < SCALAR argument version > 
     281      !! ** Purpose : water vapor at saturation in [Pa] 
     282      !!              Based on accurate estimate by Goff, 1957 
     283      !! 
     284      !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://github.com/brodeau/aerobulk/) 
     285      !! 
     286      !!    Note: what rt0 should be here, is 273.16 (triple point of water) and not 273.15 like here 
     287      !!---------------------------------------------------------------------------------- 
     288      REAL(wp), INTENT(in) ::   ptak    ! air temperature                  [K] 
     289      REAL(wp), INTENT(in) ::   pslp    ! sea level atmospheric pressure   [Pa] 
     290      REAL(wp)             ::   e_sat_sclr   ! water vapor at saturation   [kg/kg] 
     291      ! 
     292      REAL(wp) ::   zta, ztmp   ! local scalar 
     293      !!---------------------------------------------------------------------------------- 
     294      ! 
     295      zta = MAX( ptak , 180._wp )   ! air temp., prevents fpe0 errors dute to unrealistically low values over masked regions... 
     296      ztmp = rt0 / zta 
     297      ! 
     298      ! Vapour pressure at saturation [Pa] : WMO, (Goff, 1957) 
     299      e_sat_sclr = 100.*( 10.**( 10.79574*(1. - ztmp) - 5.028*LOG10(zta/rt0)        & 
     300         &    + 1.50475*10.**(-4)*(1. - 10.**(-8.2969*(zta/rt0 - 1.)) )  & 
     301         &    + 0.42873*10.**(-3)*(10.**(4.76955*(1. - ztmp)) - 1.) + 0.78614) ) 
     302      ! 
     303   END FUNCTION e_sat_sclr 
     304 
     305 
    275306   FUNCTION q_sat( ptak, pslp ) 
    276307      !!---------------------------------------------------------------------------------- 
     
    288319      ! 
    289320      INTEGER  ::   ji, jj         ! dummy loop indices 
    290       REAL(wp) ::   zta, ze_sat, ztmp   ! local scalar 
    291       !!---------------------------------------------------------------------------------- 
    292       ! 
    293       DO jj = 1, jpj 
    294          DO ji = 1, jpi 
    295             ! 
    296             zta = MAX( ptak(ji,jj) , 180._wp )   ! air temp., prevents fpe0 errors dute to unrealistically low values over masked regions... 
    297             ztmp = rt0 / zta 
    298             ! 
    299             ! Vapour pressure at saturation [hPa] : WMO, (Goff, 1957) 
    300             ze_sat = 10.**( 10.79574*(1. - ztmp) - 5.028*LOG10(zta/rt0)        & 
    301                &    + 1.50475*10.**(-4)*(1. - 10.**(-8.2969*(zta/rt0 - 1.)) )  & 
    302                &    + 0.42873*10.**(-3)*(10.**(4.76955*(1. - ztmp)) - 1.) + 0.78614  ) 
    303             ! 
    304             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] 
     321      REAL(wp) ::   ze_sat   ! local scalar 
     322      !!---------------------------------------------------------------------------------- 
     323      ! 
     324      DO jj = 1, jpj 
     325         DO ji = 1, jpi 
     326            ! 
     327            ze_sat =  e_sat_sclr( ptak(ji,jj), pslp(ji,jj) ) 
     328            ! 
     329            q_sat(ji,jj) = reps0 * ze_sat/( pslp(ji,jj) - (1._wp - reps0)*ze_sat ) 
    305330            ! 
    306331         END DO 
     
    308333      ! 
    309334   END FUNCTION q_sat 
     335 
     336   FUNCTION q_air_rh(prha, ptak, pslp) 
     337      !!---------------------------------------------------------------------------------- 
     338      !! Specific humidity of air out of Relative Humidity 
     339      !! 
     340      !! ** Author: L. Brodeau, june 2016 / AeroBulk (https://github.com/brodeau/aerobulk/) 
     341      !!---------------------------------------------------------------------------------- 
     342      REAL(wp), DIMENSION(jpi,jpj)             :: q_air_rh 
     343      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: prha        !: relative humidity      [fraction, not %!!!] 
     344      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: ptak        !: air temperature        [K] 
     345      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) :: pslp        !: atmospheric pressure   [Pa] 
     346      ! 
     347      INTEGER  ::   ji, jj      ! dummy loop indices 
     348      REAL(wp) ::   ze      ! local scalar 
     349      !!---------------------------------------------------------------------------------- 
     350      ! 
     351      DO jj = 1, jpj 
     352         DO ji = 1, jpi 
     353            ze = prha(ji,jj)*e_sat_sclr(ptak(ji,jj), pslp(ji,jj)) 
     354            q_air_rh(ji,jj) = ze*reps0/(pslp(ji,jj) - (1. - reps0)*ze) 
     355         END DO 
     356      END DO 
     357      ! 
     358   END FUNCTION q_air_rh 
    310359 
    311360 
Note: See TracChangeset for help on using the changeset viewer.