LCOV - code coverage report
Current view: | top level - physics/cam - wv_sat_methods.F90 (source / functions) |  |  | Hit | Total | Coverage
Test: | coverage.info |  | Lines: | 103 | 169 | 60.9 %
Date: | 2024-12-17 17:57:11 |  | Functions: | 18 | 28 | 64.3 %

              Line data    Source code
           1             : module wv_sat_methods
           2             :
           3             : ! This portable module contains all CAM methods for estimating
           4             : ! the saturation vapor pressure of water.
           5             : !
           6             : ! wv_saturation provides CAM-specific interfaces and utilities
           7             : ! based on these formulae.
           8             : !
           9             : ! Typical usage of this module:
          10             : !
          11             : ! Init:
          12             : ! call wv_sat_methods_init(r8, <constants>, errstring)
          13             : !
          14             : ! Get scheme index from a name string:
          15             : ! scheme_idx = wv_sat_get_scheme_idx(scheme_name)
          16             : ! if (.not. wv_sat_valid_idx(scheme_idx)) <throw some error>
          17             : !
          18             : ! Get pressures:
          19             : ! es = wv_sat_svp_water(t, scheme_idx)
          20             : ! es = wv_sat_svp_ice(t, scheme_idx)
          21             : !
          22             : ! Use ice/water transition range:
          23             : ! es = wv_sat_svp_trice(t, ttrice, scheme_idx)
          24             : !
          25             : ! Note that elemental functions cannot be pointed to, nor passed
          26             : ! as arguments. If you need to do either, it is recommended to
          27             : ! wrap the function so that it can be given an explicit (non-
          28             : ! elemental) interface.
          29             :
          30             : implicit none
          31             : private
          32             : save
          33             :
          34             : integer, parameter :: r8 = selected_real_kind(12) ! 8 byte real
          35             :
          36             : integer, parameter :: VLENS = 128 ! vector length for a GPU kernel
          37             :
          38             : real(r8) :: tmelt   ! Melting point of water at 1 atm (K)
          39             : real(r8) :: h2otrip ! Triple point temperature of water (K)
          40             : real(r8) :: tboil   ! Boiling point of water at 1 atm (K)
          41             :
          42             : real(r8) :: ttrice  ! Ice-water transition range
          43             :
          44             : real(r8) :: epsilo  ! Ice-water transition range
          45             : real(r8) :: omeps   ! 1._r8 - epsilo
          46             :
          47             : ! Indices representing individual schemes
          48             : integer, parameter :: Invalid_idx = -1
          49             : integer, parameter :: GoffGratch_idx = 1
          50             : integer, parameter :: MurphyKoop_idx = 2
          51             : integer, parameter :: Bolton_idx = 3
          52             :
          53             : ! Index representing the current default scheme.
          54             : integer, parameter :: initial_default_idx = GoffGratch_idx
          55             : integer :: default_idx = initial_default_idx
          56             :
          57             : !$acc declare create (epsilo, tmelt, tboil, omeps, h2otrip, ttrice)
          58             :
          59             : public wv_sat_methods_init
          60             : public wv_sat_get_scheme_idx
          61             : public wv_sat_valid_idx
          62             :
          63             : public wv_sat_set_default
          64             : public wv_sat_reset_default
          65             :
          66             : public wv_sat_qsat_water, wv_sat_qsat_water_vect
          67             : public wv_sat_qsat_ice, wv_sat_qsat_ice_vect
          68             :
          69             : public wv_sat_svp_trans, wv_sat_svp_trans_vect
          70             :
          71             : ! pressure -> humidity conversion
          72             : public wv_sat_svp_to_qsat, wv_sat_svp_to_qsat_vect
          73             :
          74             : ! Combined qsat operations
          75             : public wv_sat_qsat_trans
          76             :
          77             : public wv_sat_svp_water, wv_sat_svp_water_vect
          78             : public wv_sat_svp_ice, wv_sat_svp_ice_vect
          79             :
          80             : contains
          81             :
          82             : !---------------------------------------------------------------------
          83             : ! ADMINISTRATIVE FUNCTIONS
          84             : !---------------------------------------------------------------------
          85             :
          86             : ! Get physical constants
          87        1536 : subroutine wv_sat_methods_init(kind, tmelt_in, h2otrip_in, tboil_in, &
          88             :      ttrice_in, epsilo_in, errstring)
          89             :   integer, intent(in) :: kind
          90             :   real(r8), intent(in) :: tmelt_in
          91             :   real(r8), intent(in) :: h2otrip_in
          92             :   real(r8), intent(in) :: tboil_in
          93             :   real(r8), intent(in) :: ttrice_in
          94             :   real(r8), intent(in) :: epsilo_in
          95             :   character(len=*), intent(out)  :: errstring
          96             :
          97        1536 :   errstring = ' '
          98             :
          99        1536 :   if (kind /= r8) then
         100           0 :      write(errstring,*) 'wv_sat_methods_init: ERROR: ', &
         101           0 :           kind,' was input kind but ',r8,' is internal kind.'
         102           0 :      return
         103             :   end if
         104             :
         105        1536 :   if (ttrice_in < 0._r8) then
         106           0 :      write(errstring,*) 'wv_sat_methods_init: ERROR: ', &
         107           0 :           ttrice_in,' was input for ttrice, but negative range is invalid.'
         108           0 :      return
         109             :   end if
         110             :
         111        1536 :   tmelt = tmelt_in
         112        1536 :   h2otrip = h2otrip_in
         113        1536 :   tboil = tboil_in
         114        1536 :   ttrice = ttrice_in
         115        1536 :   epsilo = epsilo_in
         116             :
         117        1536 :   omeps = 1._r8 - epsilo
         118             :
         119             :   !$acc update device (tmelt,h2otrip,tboil,ttrice,epsilo,omeps)
         120             :
         121        1536 : end subroutine wv_sat_methods_init
         122             :
         123             : ! Look up index by name.
         124        1536 : pure function wv_sat_get_scheme_idx(name) result(idx)
         125             :   character(len=*), intent(in) :: name
         126             :   integer :: idx
         127             :
         128             :   select case (name)
         129             :   case("GoffGratch")
         130           0 :      idx = GoffGratch_idx
         131             :   case("MurphyKoop")
         132           0 :      idx = MurphyKoop_idx
         133             :   case("Bolton")
         134           0 :      idx = Bolton_idx
         135             :   case default
         136        1536 :      idx = Invalid_idx
         137             :   end select
         138             :
         139        1536 : end function wv_sat_get_scheme_idx
         140             :
         141             : ! Check validity of an index from the above routine.
         142        1536 : pure function wv_sat_valid_idx(idx) result(status)
         143             :   integer, intent(in) :: idx
         144             :   logical :: status
         145             :
         146        1536 :   status = (idx /= Invalid_idx)
         147             :
         148        1536 : end function wv_sat_valid_idx
         149             :
         150             : ! Set default scheme (otherwise, Goff & Gratch is default)
         151             : ! Returns a logical representing success (.true.) or
         152             : ! failure (.false.).
         153        1536 : function wv_sat_set_default(name) result(status)
         154             :   character(len=*), intent(in) :: name
         155             :   logical :: status
         156             :
         157             :   ! Don't want to overwrite valid default with invalid,
         158             :   ! so assign to temporary and check it first.
         159             :   integer :: tmp_idx
         160             :
         161        1536 :   tmp_idx = wv_sat_get_scheme_idx(name)
         162             :
         163        1536 :   status = wv_sat_valid_idx(tmp_idx)
         164             :
         165        1536 :   if (status) default_idx = tmp_idx
         166             :
         167        1536 : end function wv_sat_set_default
         168             :
         169             : ! Reset default scheme to initial value.
         170             : ! The same thing can be accomplished with wv_sat_set_default;
         171             : ! the real reason to provide this routine is to reset the
         172             : ! module for testing purposes.
         173           0 : subroutine wv_sat_reset_default()
         174             :
         175           0 :   default_idx = initial_default_idx
         176             :
         177           0 : end subroutine wv_sat_reset_default
         178             :
         179             : !---------------------------------------------------------------------
         180             : ! UTILITIES
         181             : !---------------------------------------------------------------------
         182             :
         183             : ! Get saturation specific humidity given pressure and SVP.
         184             : ! Specific humidity is limited to range 0-1.
         185 29532646494 : function wv_sat_svp_to_qsat(es, p) result(qs)
         186             :   real(r8), intent(in) :: es  ! SVP
         187             :   real(r8), intent(in) :: p   ! Current pressure.
         188             :   real(r8)             :: qs
         189             :
         190             :   ! If pressure is less than SVP, set qs to maximum of 1.
         191 29532646494 :   if ( (p - es) <= 0._r8 ) then
         192             :      qs = 1.0_r8
         193             :   else
         194 29499338570 :      qs = epsilo*es / (p - omeps*es)
         195             :   end if
         196             :
         197 29532646494 : end function wv_sat_svp_to_qsat
         198             :
         199             : ! Get saturation specific humidity given pressure and SVP.
         200             : ! Specific humidity is limited to range 0-1.
         201  1195182936 : subroutine  wv_sat_svp_to_qsat_vect(es, p, qs, vlen)
         202             :
         203             :   integer,  intent(in) :: vlen
         204             :   real(r8), intent(in)  :: es(vlen)  ! SVP
         205             :   real(r8), intent(in)  :: p(vlen)   ! Current pressure.
         206             :   real(r8), intent(out) :: qs(vlen)
         207             :   integer :: i
         208             :
         209             :   ! If pressure is less than SVP, set qs to maximum of 1.
         210             :
         211             :   !$acc data present (es,p,qs)
         212             :
         213             :   !$acc parallel vector_length(VLENS) default(present)
         214             :   !$acc loop gang vector
         215 54881415936 :   do i=1,vlen
         216 54881415936 :      if ( (p(i) - es(i)) <= 0._r8 ) then
         217  1351491030 :         qs(i) = 1.0_r8
         218             :      else
         219 52334741970 :         qs(i) = epsilo*es(i) / (p(i) - omeps*es(i))
         220             :      end if
         221             :   end do
         222             :   !$acc end parallel
         223             :
         224             :   !$acc end data
         225  1195182936 : end subroutine wv_sat_svp_to_qsat_vect
         226             :
         227 29020232765 : subroutine wv_sat_qsat_water(t, p, es, qs, idx)
         228             :   !------------------------------------------------------------------!
         229             :   ! Purpose:                                                         !
         230             :   !   Calculate SVP over water at a given temperature, and then      !
         231             :   !   calculate and return saturation specific humidity.             !
         232             :   !------------------------------------------------------------------!
         233             :
         234             :   ! Inputs
         235             :   real(r8), intent(in) :: t    ! Temperature
         236             :   real(r8), intent(in) :: p    ! Pressure
         237             :   ! Outputs
         238             :   real(r8), intent(out) :: es  ! Saturation vapor pressure
         239             :   real(r8), intent(out) :: qs  ! Saturation specific humidity
         240             :
         241             :   integer,  intent(in), optional :: idx ! Scheme index
         242             :
         243 29020232765 :   es = wv_sat_svp_water(t, idx)
         244             :
         245 29020232765 :   qs = wv_sat_svp_to_qsat(es, p)
         246             :
         247             :   ! Ensures returned es is consistent with limiters on qs.
         248 29020232765 :   es = min(es, p)
         249             :
         250 29020232765 : end subroutine wv_sat_qsat_water
         251             :
         252   768414816 : subroutine wv_sat_qsat_water_vect(t, p, es, qs, vlen, idx)
         253             :   !------------------------------------------------------------------!
         254             :   ! Purpose:                                                         !
         255             :   !   Calculate SVP over water at a given temperature, and then      !
         256             :   !   calculate and return saturation specific humidity.             !
         257             :   !------------------------------------------------------------------!
         258             :   ! Inputs
         259             :
         260             :   integer,  intent(in) :: vlen
         261             :   real(r8), intent(in) :: t(vlen)    ! Temperature
         262             :   real(r8), intent(in) :: p(vlen)    ! Pressure
         263             :   ! Outputs
         264             :   real(r8), intent(out) :: es(vlen)  ! Saturation vapor pressure
         265             :   real(r8), intent(out) :: qs(vlen)  ! Saturation specific humidity
         266             :
         267             :   integer,  intent(in), optional :: idx ! Scheme index
         268             :   integer :: i
         269             :
         270             :   !$acc data present (t,p,es,qs)
         271             :
         272   768414816 :   call wv_sat_svp_water_vect(t, es, vlen, idx)
         273   768414816 :   call wv_sat_svp_to_qsat_vect(es, p, qs, vlen)
         274             :
         275             :   !$acc parallel vector_length(VLENS) default(present)
         276             :   !$acc loop gang vector
         277 36113834016 :   do i=1,vlen
         278             :      ! Ensures returned es is consistent with limiters on qs.
         279 36113834016 :      es(i) = min(es(i), p(i))
         280             :   enddo
         281             :   !$acc end parallel
         282             :
         283             :   !$acc end data
         284   768414816 : end subroutine wv_sat_qsat_water_vect
         285             :
         286           0 : subroutine wv_sat_qsat_ice(t, p, es, qs, idx)
         287             :   !------------------------------------------------------------------!
         288             :   ! Purpose:                                                         !
         289             :   !   Calculate SVP over ice at a given temperature, and then        !
         290             :   !   calculate and return saturation specific humidity.             !
         291             :   !------------------------------------------------------------------!
         292             :
         293             :   ! Inputs
         294             :   real(r8), intent(in) :: t    ! Temperature
         295             :   real(r8), intent(in) :: p    ! Pressure
         296             :   ! Outputs
         297             :   real(r8), intent(out) :: es  ! Saturation vapor pressure
         298             :   real(r8), intent(out) :: qs  ! Saturation specific humidity
         299             :
         300             :   integer,  intent(in), optional :: idx ! Scheme index
         301             :
         302           0 :   es = wv_sat_svp_ice(t, idx)
         303             :
         304           0 :   qs = wv_sat_svp_to_qsat(es, p)
         305             :
         306             :   ! Ensures returned es is consistent with limiters on qs.
         307           0 :   es = min(es, p)
         308             :
         309           0 : end subroutine wv_sat_qsat_ice
         310             :
         311     8935056 : subroutine wv_sat_qsat_ice_vect(t, p, es, qs, vlen, idx)
         312             :   !------------------------------------------------------------------!
         313             :   ! Purpose:                                                         !
         314             :   !   Calculate SVP over ice at a given temperature, and then        !
         315             :   !   calculate and return saturation specific humidity.             !
         316             :   !------------------------------------------------------------------!
         317             :   ! Inputs
         318             :
         319             :   integer,  intent(in) :: vlen
         320             :   real(r8), intent(in) :: t(vlen)    ! Temperature
         321             :   real(r8), intent(in) :: p(vlen)    ! Pressure
         322             :   ! Outputs
         323             :   real(r8), intent(out) :: es(vlen)  ! Saturation vapor pressure
         324             :   real(r8), intent(out) :: qs(vlen)  ! Saturation specific humidity
         325             :
         326             :   integer,  intent(in), optional :: idx ! Scheme index
         327             :   integer :: i
         328             :
         329             :   !$acc data present (t,p,es,qs)
         330             :
         331     8935056 :   call wv_sat_svp_ice_vect(t, es, vlen, idx)
         332     8935056 :   call wv_sat_svp_to_qsat_vect(es, p, qs, vlen)
         333             :
         334             :   !$acc parallel vector_length(VLENS) default(present)
         335             :   !$acc loop gang vector
         336 11790741456 :   do i=1,vlen
         337             :      ! Ensures returned es is consistent with limiters on qs.
         338 11790741456 :      es(i) = min(es(i), p(i))
         339             :   enddo
         340             :   !$acc end parallel
         341             :
         342             :   !$acc end data
         343     8935056 : end subroutine wv_sat_qsat_ice_vect
         344             :
         345           0 : subroutine wv_sat_qsat_trans(t, p, es, qs, idx)
         346             :   !------------------------------------------------------------------!
         347             :   ! Purpose:                                                         !
         348             :   !   Calculate SVP over ice at a given temperature, and then        !
         349             :   !   calculate and return saturation specific humidity.             !
         350             :   !------------------------------------------------------------------!
         351             :
         352             :   ! Inputs
         353             :   real(r8), intent(in) :: t    ! Temperature
         354             :   real(r8), intent(in) :: p    ! Pressure
         355             :   ! Outputs
         356             :   real(r8), intent(out) :: es  ! Saturation vapor pressure
         357             :   real(r8), intent(out) :: qs  ! Saturation specific humidity
         358             :
         359             :   integer,  intent(in), optional :: idx ! Scheme index
         360             :
         361           0 :   es = wv_sat_svp_trans(t, idx)
         362             :
         363           0 :   qs = wv_sat_svp_to_qsat(es, p)
         364             :
         365             :   ! Ensures returned es is consistent with limiters on qs.
         366           0 :   es = min(es, p)
         367             :
         368           0 : end subroutine wv_sat_qsat_trans
         369             :
         370             : !---------------------------------------------------------------------
         371             : ! SVP INTERFACE FUNCTIONS
         372             : !---------------------------------------------------------------------
         373             :
         374 33670967524 : function wv_sat_svp_water(t, idx) result(es)
         375             :   real(r8), intent(in) :: t
         376             :   integer,  intent(in), optional :: idx
         377             :   real(r8)             :: es
         378             :
         379             :   integer :: use_idx
         380             :
         381 33670967524 :   if (present(idx)) then
         382           0 :      use_idx = idx
         383             :   else
         384 33670967524 :      use_idx = default_idx
         385             :   end if
         386             :
         387 33670967524 :   select case (use_idx)
         388             :   case(GoffGratch_idx)
         389 33670967524 :      es = GoffGratch_svp_water(t)
         390             :   case(MurphyKoop_idx)
         391           0 :      es = MurphyKoop_svp_water(t)
         392             :   case(Bolton_idx)
         393 33670967524 :      es = Bolton_svp_water(t)
         394             :   end select
         395             :
         396 33670967524 : end function wv_sat_svp_water
         397             :
         398  1143687168 : subroutine  wv_sat_svp_water_vect(t, es, vlen, idx)
         399             :   integer,  intent(in) :: vlen
         400             :   real(r8), intent(in) :: t(vlen)
         401             :   integer,  intent(in), optional :: idx
         402             :   real(r8), intent(out) :: es(vlen)
         403             :   integer :: i
         404             :   integer :: use_idx
         405             :
         406             :   !$acc data present (t,es)
         407             :
         408  1143687168 :   if (present(idx)) then
         409           0 :      use_idx = idx
         410             :   else
         411  1143687168 :      use_idx = default_idx
         412             :   end if
         413             :
         414  1143687168 :   select case (use_idx)
         415             :   case(GoffGratch_idx)
         416  1143687168 :      call GoffGratch_svp_water_vect(t,es,vlen)
         417             :   case(MurphyKoop_idx)
         418           0 :      call MurphyKoop_svp_water_vect(t,es,vlen)
         419             :   case(Bolton_idx)
         420  1143687168 :      call Bolton_svp_water_vect(t,es,vlen)
         421             :   end select
         422             :
         423             :   !$acc end data
         424  1143687168 : end subroutine wv_sat_svp_water_vect
         425             :
         426  4650767015 : function wv_sat_svp_ice(t, idx) result(es)
         427             :   real(r8), intent(in)  :: t
         428             :   integer,  intent(in), optional :: idx
         429             :   real(r8)              :: es
         430             :
         431             :   integer :: use_idx
         432             :
         433  4650767015 :   if (present(idx)) then
         434           0 :      use_idx = idx
         435             :   else
         436  4650767015 :      use_idx = default_idx
         437             :   end if
         438             :
         439  4650767015 :   select case (use_idx)
         440             :   case(GoffGratch_idx)
         441  4650767015 :      es = GoffGratch_svp_ice(t)
         442             :   case(MurphyKoop_idx)
         443           0 :      es = MurphyKoop_svp_ice(t)
         444             :   case(Bolton_idx)
         445  4650767015 :      es = Bolton_svp_water(t)
         446             :   end select
         447             :
         448  4650767015 : end function wv_sat_svp_ice
         449             :
         450   384207408 : subroutine wv_sat_svp_ice_vect(t, es, vlen, idx)
         451             :   integer,  intent(in) :: vlen
         452             :   real(r8), intent(in) :: t(vlen)
         453             :   integer,  intent(in), optional :: idx
         454             :   real(r8), intent(out) :: es(vlen)
         455             :   integer :: i
         456             :
         457             :   integer :: use_idx
         458             :
         459             :   !$acc data present (t,es)
         460             :
         461   384207408 :   if (present(idx)) then
         462           0 :      use_idx = idx
         463             :   else
         464   384207408 :      use_idx = default_idx
         465             :   end if
         466             :
         467   384207408 :   select case (use_idx)
         468             :   case(GoffGratch_idx)
         469   384207408 :      call GoffGratch_svp_ice_vect(t,es,vlen)
         470             :   case(MurphyKoop_idx)
         471           0 :      call MurphyKoop_svp_ice_vect(t,es,vlen)
         472             :   case(Bolton_idx)
         473   384207408 :      call Bolton_svp_water_vect(t,es,vlen)
         474             :   end select
         475             :
         476             :   !$acc end data
         477   384207408 : end subroutine wv_sat_svp_ice_vect
         478             :
         479      385536 : function wv_sat_svp_trans(t, idx) result(es)
         480             :
         481             :   real(r8), intent(in) :: t
         482             :   integer,  intent(in), optional :: idx
         483             :   real(r8) :: es
         484             :
         485             :   real(r8) :: esice      ! Saturation vapor pressure over ice
         486             :   real(r8) :: weight     ! Intermediate scratch variable for es transition
         487             :
         488             : !
         489             : ! Water
         490             : !
         491      385536 :   if (t >= (tmelt - ttrice)) then
         492      192000 :      es = wv_sat_svp_water(t,idx)
         493             :   else
         494             :      es = 0.0_r8
         495             :   end if
         496             :
         497             : !
         498             : ! Ice
         499             : !
         500      385536 :   if (t < tmelt) then
         501             :
         502      224256 :      esice = wv_sat_svp_ice(t,idx)
         503             :
         504      224256 :      if ( (tmelt - t) > ttrice ) then
         505             :         weight = 1.0_r8
         506             :      else
         507       30720 :         weight = (tmelt - t)/ttrice
         508             :      end if
         509             :
         510      224256 :      es = weight*esice + (1.0_r8 - weight)*es
         511             :   end if
         512             :
         513      385536 : end function wv_sat_svp_trans
         514             :
         515           0 : subroutine wv_sat_svp_trans_vect(t, es, vlen, idx)
         516             :
         517             :   integer,  intent(in)  :: vlen
         518             :   real(r8), intent(in)  :: t(vlen)
         519             :   integer,  intent(in), optional :: idx
         520             :   real(r8), intent(out) :: es(vlen)
         521             :
         522           0 :   real(r8) :: esice(vlen)      ! Saturation vapor pressure over ice
         523             :   real(r8) :: weight           ! Intermediate scratch variable for es transition
         524             :   integer  :: i
         525             :
         526             :   !$acc data present (t,es) &
         527             :   !$acc      create  (esice)
         528             :
         529             : !
         530             : ! Water
         531             : !
         532           0 :   call wv_sat_svp_water_vect(t,es,vlen,idx)
         533             :   !$acc parallel vector_length(VLENS) default(present)
         534             :   !$acc loop gang vector
         535           0 :   do i = 1, vlen
         536           0 :      if (t(i) < (tmelt - ttrice)) then
         537           0 :         es(i) = 0.0_r8
         538             :      end if
         539             :   end do
         540             :   !$acc end parallel
         541             : !
         542             : ! Ice
         543             : !
         544           0 :   call wv_sat_svp_ice_vect(t,esice,vlen,idx)
         545             :   !$acc parallel vector_length(VLENS) default(present)
         546             :   !$acc loop gang vector
         547           0 :   do i = 1, vlen
         548           0 :      if (t(i) < tmelt) then
         549           0 :         if ( (tmelt - t(i)) > ttrice ) then
         550             :            weight = 1.0_r8
         551             :         else
         552           0 :            weight = (tmelt - t(i))/ttrice
         553             :         end if
         554             :
         555           0 :         es(i) = weight*esice(i) + (1.0_r8 - weight)*es(i)
         556             :      end if
         557             :   end do
         558             :   !$acc end parallel
         559             :
         560             :   !$acc end data
         561           0 : end subroutine wv_sat_svp_trans_vect
         562             :
         563             : !---------------------------------------------------------------------
         564             : ! SVP METHODS
         565             : !---------------------------------------------------------------------
         566             :
         567             : ! Goff & Gratch (1946)
         568             :
         569 33670967524 : function GoffGratch_svp_water(t) result(es)
         570             :   real(r8), intent(in) :: t  ! Temperature in Kelvin
         571             :   real(r8) :: es             ! SVP in Pa
         572             :
         573             :   ! uncertain below -70 C
         574             :   es = 10._r8**(-7.90298_r8*(tboil/t-1._r8)+ &
         575             :        5.02808_r8*log10(tboil/t)- &
         576             :        1.3816e-7_r8*(10._r8**(11.344_r8*(1._r8-t/tboil))-1._r8)+ &
         577             :        8.1328e-3_r8*(10._r8**(-3.49149_r8*(tboil/t-1._r8))-1._r8)+ &
         578 33670967524 :        log10(1013.246_r8))*100._r8
         579             :
         580 33670967524 : end function GoffGratch_svp_water
         581             :
         582  1143687168 : subroutine GoffGratch_svp_water_vect(t, es, vlen)
         583             :   integer, intent(in) :: vlen
         584             :   real(r8), intent(in)  :: t(vlen)  ! Temperature in Kelvin
         585             :   real(r8), intent(out) :: es(vlen) ! SVP in Pa
         586             :   real(r8) :: log_tboil
         587             :   integer :: i
         588             :
         589             :   !$acc data present (t,es)
         590             :
         591             :   ! Goff, J. A., and S. Gratch. “Low-Pressure Properties of Water from -160F
         592             :   ! to 212F.” Trans. Am. Soc. Heat. Vent. Eng. 52 (1946): 95–121.
         593             :   ! uncertain below -70 C
         594             :
         595  1143687168 :   log_tboil = log10(tboil)
         596             :
         597             :   !$acc parallel vector_length(VLENS) default(present)
         598             :   !$acc loop gang vector
         599 42380009568 :   do i=1,vlen
         600 41236322400 :      es(i) = 10._r8**(-7.90298_r8*(tboil/t(i)-1._r8)+ &
         601             :        5.02808_r8*(log_tboil-log10(t(i)))- &
         602             :        1.3816e-7_r8*(10._r8**(11.344_r8*(1._r8-t(i)/tboil))-1._r8)+ &
         603             :        8.1328e-3_r8*(10._r8**(-3.49149_r8*(tboil/t(i)-1._r8))-1._r8)+ &
         604 83616331968 :        log10(1013.246_r8))*100._r8
         605             :   enddo
         606             :   !$acc end parallel
         607             :
         608             :   !$acc end data
         609  1143687168 : end subroutine GoffGratch_svp_water_vect
         610             :
         611  4650767015 : function GoffGratch_svp_ice(t) result(es)
         612             :   real(r8), intent(in) :: t  ! Temperature in Kelvin
         613             :   real(r8) :: es             ! SVP in Pa
         614             :
         615             :   ! good down to -100 C
         616             :   es = 10._r8**(-9.09718_r8*(h2otrip/t-1._r8)-3.56654_r8* &
         617             :        log10(h2otrip/t)+0.876793_r8*(1._r8-t/h2otrip)+ &
         618  4650767015 :        log10(6.1071_r8))*100._r8
         619             :
         620  4650767015 : end function GoffGratch_svp_ice
         621             :
         622   384207408 : subroutine GoffGratch_svp_ice_vect(t, es, vlen)
         623             :   integer, intent(in)   :: vlen
         624             :   real(r8), intent(in)  :: t(vlen)  ! Temperature in Kelvin
         625             :   real(r8), intent(out) :: es(vlen) ! SVP in Pa
         626             :   real(r8), parameter :: log_param = log10(6.1071_r8)
         627             :   integer :: i
         628             :   ! good down to -100 C
         629             :
         630             :   !$acc data present (t,es)
         631             :
         632             :   !$acc parallel vector_length(VLENS) default(present)
         633             :   !$acc loop gang vector
         634 18056917008 :   do i=1,vlen
         635 17672709600 :      es(i) = 10._r8**(-9.09718_r8*(h2otrip/t(i)-1._r8)-3.56654_r8* &
         636             :           log10(h2otrip/t(i))+0.876793_r8*(1._r8-t(i)/h2otrip)+ &
         637 35729626608 :           log_param)*100._r8
         638             :   enddo
         639             :   !$acc end parallel
         640             :
         641             :   !$acc end data
         642   384207408 : end subroutine GoffGratch_svp_ice_vect
         643             :
         644             : ! Murphy & Koop (2005)
         645             :
         646           0 : function MurphyKoop_svp_water(t) result(es)
         647             :   real(r8), intent(in) :: t  ! Temperature in Kelvin
         648             :   real(r8) :: es             ! SVP in Pa
         649             :
         650             :   ! (good for 123 < T < 332 K)
         651             :   es = exp(54.842763_r8 - (6763.22_r8 / t) - (4.210_r8 * log(t)) + &
         652             :        (0.000367_r8 * t) + (tanh(0.0415_r8 * (t - 218.8_r8)) * &
         653             :        (53.878_r8 - (1331.22_r8 / t) - (9.44523_r8 * log(t)) + &
         654           0 :        0.014025_r8 * t)))
         655             :
         656           0 : end function MurphyKoop_svp_water
         657             :
         658           0 : subroutine MurphyKoop_svp_water_vect(t, es, vlen)
         659             :   integer, intent(in)   :: vlen
         660             :   real(r8), intent(in)  :: t(vlen)  ! Temperature in Kelvin
         661             :   real(r8), intent(out) :: es(vlen)             ! SVP in Pa
         662             :
         663             :   integer :: i
         664             :   ! Murphy, D. M., and T. Koop. “Review of the Vapour Pressure of Ice and
         665             :   ! Supercooled Water for Atmospheric Applications.” Q. J. R. Meteorol.
         666             :   ! Soc. 131, no. 608 (2005): 1539–65. 10.1256/qj.04.94
         667             :   ! (good for 123 < T < 332 K)
         668             :
         669             :   !$acc data present (t,es)
         670             :
         671             :   !$acc parallel vector_length(VLENS) default(present)
         672             :   !$acc loop gang vector
         673           0 :   do i = 1, vlen
         674           0 :      es(i) = exp(54.842763_r8 - (6763.22_r8 / t(i)) - (4.210_r8 * log(t(i))) + &
         675             :           (0.000367_r8 * t(i)) + (tanh(0.0415_r8 * (t(i) - 218.8_r8)) * &
         676             :           (53.878_r8 - (1331.22_r8 / t(i)) - (9.44523_r8 * log(t(i))) + &
         677           0 :           0.014025_r8 * t(i))))
         678             :   end do
         679             :   !$acc end parallel
         680             :
         681             :   !$acc end data
         682           0 : end subroutine MurphyKoop_svp_water_vect
         683             :
         684           0 : function MurphyKoop_svp_ice(t) result(es)
         685             :   real(r8), intent(in) :: t  ! Temperature in Kelvin
         686             :   real(r8) :: es             ! SVP in Pa
         687             :
         688             :   ! (good down to 110 K)
         689             :   es = exp(9.550426_r8 - (5723.265_r8 / t) + (3.53068_r8 * log(t)) &
         690           0 :        - (0.00728332_r8 * t))
         691             :
         692           0 : end function MurphyKoop_svp_ice
         693             :
         694           0 : subroutine MurphyKoop_svp_ice_vect(t, es, vlen)
         695             :   integer, intent(in)   :: vlen
         696             :   real(r8), intent(in) :: t(vlen)  ! Temperature in Kelvin
         697             :   real(r8), intent(out) :: es(vlen)             ! SVP in Pa
         698             :
         699             :   integer :: i
         700             :   ! (good down to 110 K)
         701             :
         702             :   !$acc data present (t,es)
         703             :
         704             :   !$acc parallel vector_length(VLENS) default(present)
         705             :   !$acc loop gang vector
         706           0 :   do i = 1, vlen
         707           0 :      es(i) = exp(9.550426_r8 - (5723.265_r8 / t(i)) + (3.53068_r8 * log(t(i))) &
         708           0 :              - (0.00728332_r8 * t(i)))
         709             :   end do
         710             :   !$acc end parallel
         711             :
         712             :   !$acc end data
         713           0 : end subroutine MurphyKoop_svp_ice_vect
         714             :
         715             : ! Bolton (1980)
         716             : ! zm_conv deep convection scheme contained this SVP calculation.
         717             : ! It appears to be from D. Bolton, 1980, Monthly Weather Review.
         718             : ! Unlike the other schemes, no distinct ice formula is associated
         719             : ! with it. (However, a Bolton ice formula exists in CLUBB.)
         720             :
         721             : ! The original formula used degrees C, but this function
         722             : ! takes Kelvin and internally converts.
         723             :
         724           0 : function Bolton_svp_water(t) result(es)
         725             :   real(r8),parameter :: c1 = 611.2_r8
         726             :   real(r8),parameter :: c2 = 17.67_r8
         727             :   real(r8),parameter :: c3 = 243.5_r8
         728             :
         729             :   real(r8), intent(in) :: t  ! Temperature in Kelvin
         730             :   real(r8) :: es             ! SVP in Pa
         731             :
         732           0 :   es = c1*exp( (c2*(t - tmelt))/((t - tmelt)+c3) )
         733             :
         734           0 : end function Bolton_svp_water
         735             :
         736           0 : subroutine Bolton_svp_water_vect(t, es,vlen)
         737             :   real(r8),parameter :: c1 = 611.2_r8
         738             :   real(r8),parameter :: c2 = 17.67_r8
         739             :   real(r8),parameter :: c3 = 243.5_r8
         740             :
         741             :   integer, intent(in)   :: vlen
         742             :   real(r8), intent(in)  :: t(vlen)  ! Temperature in Kelvin
         743             :   real(r8), intent(out) :: es(vlen)             ! SVP in Pa
         744             :
         745             :   integer :: i
         746             :
         747             :   !$acc data present (t,es)
         748             :
         749             :   !$acc parallel vector_length(VLENS) default(present)
         750             :   !$acc loop gang vector
         751           0 :   do i = 1, vlen
         752           0 :      es(i) = c1*exp( (c2*(t(i) - tmelt))/((t(i) - tmelt)+c3) )
         753             :   end do
         754             :   !$acc end parallel
         755             :
         756             :   !$acc end data
         757           0 : end subroutine Bolton_svp_water_vect
         758             :
         759             : end module wv_sat_methods

Generated by: LCOV version 1.14