! Version 1.0 released by David Romps on April 15, 2021.
! 
! When using this code, please cite:
! 
! @article{20dewpoint,
!   Title   = {Accurate expressions for the dew point and frost point derived from the {Rankine-Kirchhoff} approximations},
!   Author  = {David M. Romps},
!   Journal = {Journal of the Atmospheric Sciences},
!   Year    = {2021},
!   Month   = jul,
!   Number  = {7},
!   Pages   = {2113--2116},
!   Volume  = {78}
! }
!
! This dew-point function returns the dewpoint (Td) in K.
! The inputs are:
! - p in Pascals
! - T in Kelvins
! - Exactly one of rh, rhl, and rhs (dimensionless, from 0 to 1):
!    * The value of rh is interpreted to be the relative humidity with
!      respect to liquid water if T >= 273.15 K and with respect to ice if
!      T < 273.15 K. 
!    * The value of rhl is interpreted to be the relative humidity with
!      respect to liquid water
!    * The value of rhs is interpreted to be the relative humidity with
!      respect to ice
! - return_fp is an optional logical flag.  If true, the frost point (Tf)
!   is returned instead of the dew point (Td). 
! - return_max_dp_fp is an optional logical flag.  If true, the maximum of the
!   dew point (Td) and frost point (Tf) is returned.

program test

   use dewpoint_mod

   implicit none

   double precision :: p, rh, T1, T2

   p = 1.e5
   rh = 0.5
   T1 = 300.
   T2 = 200.

   if ( &
      abs(dewpoint(T1,rhl=rh,return_fp=.false.)/288.7153070587-1.) < 1.e-10 .and. &
      abs(dewpoint(T1,rhs=rh,return_fp=.false.)/292.8006747771-1.) < 1.e-10 .and. &
      abs(dewpoint(T2,rhl=rh,return_fp=.false.)/195.3232095158-1.) < 1.e-10 .and. &
      abs(dewpoint(T2,rhs=rh,return_fp=.false.)/190.8053512419-1.) < 1.e-10 .and. &
      abs(dewpoint(T1,rhl=rh,return_fp=.true. )/286.6935602684-1.) < 1.e-10 .and. &
      abs(dewpoint(T1,rhs=rh,return_fp=.true. )/290.1843594643-1.) < 1.e-10 .and. &
      abs(dewpoint(T2,rhl=rh,return_fp=.true. )/200.0743709369-1.) < 1.e-10 .and. &
      abs(dewpoint(T2,rhs=rh,return_fp=.true. )/195.5867236563-1.) < 1.e-10 ) then
      print *,'Success'
   else
      print *,'Failure'
   end if

end program test