UPP (develop)
Loading...
Searching...
No Matches
UPP_PHYSICS.f File Reference

UPP_PHYSICS is a collection of UPP subroutines for physics variables calculation. More...

Go to the source code of this file.

Functions/Subroutines

subroutine, public upp_physics::calcape (itype, dpbnd, p1d, t1d, q1d, l1d, cape, cins, pparc, zeql, thund)
 Computes CAPE and CINS.
 
subroutine, public upp_physics::calcape2 (itype, dpbnd, p1d, t1d, q1d, l1d, cape, cins, lfc, esrhl, esrhh, dcape, dgld, esp)
 Computes CAPE and CINS.
 
subroutine, public upp_physics::caldiv (uwnd, vwnd, div)
 Computes divergence.
 
subroutine, public upp_physics::calgradps (ps, psx, psy)
 Computes gradients of a scalar field PS or LNPS.
 
subroutine, public upp_physics::calrh (p1, t1, q1, rh)
 Computes relative humidity.
 
subroutine, public upp_physics::calrh_gfs (p1, t1, q1, rh)
 Computes relative humidity.
 
subroutine, public upp_physics::calrh_gsd (p1, t1, q1, rhb)
 Compute RH with the NOAA GSL (formerly NOAA GSD) algorithm used for RUC and Rapid Refresh.
 
subroutine, public upp_physics::calrh_nam (p1, t1, q1, rh)
 Computes relative humidity.
 
subroutine, public upp_physics::calrh_pw (rhpw)
 Algorithm used at GSL for RUC and Rapid Refresh.
 
subroutine, public upp_physics::calslr_roebber (tprs, rhprs, slr)
 Computes snow solid-liquid-ratio (SLR) using the Roebber algorithm.
 
subroutine, public upp_physics::calslr_uutah (slr)
 Computes snow solid-liquid-ratio slr using the Steenburgh algorithm.
 
subroutine, public upp_physics::calvor (uwnd, vwnd, absv)
 Computes absolute vorticity.
 
elemental real function, public upp_physics::fpvsnew (t)
 Computes saturation vapor pressure.
 
elemental real function, public upp_physics::tvirtual (t, q)
 Computes virtual temperature.
 

Detailed Description

UPP_PHYSICS is a collection of UPP subroutines for physics variables calculation.

calcape() computes CAPE/CINS and other storm related variables.

calcape2() computes additional storm related variables.

caldiv() computes divergence.

calgradps() computes gardients of a scalar field PS or LNPS.

calrh(), calrh_nam(), calrh_gfs(), calrh_gsd() compute RH using various algorithms.

The NAM v4.1.18 algorithm (calrh_nam()) is selected as default for NMMB and FV3GFS, FV3GEFS, and FV3R for the UPP 2020 unification.

calrh_pw() algorithm use at GSD for RUC and Rapid Refresh.

calslr_roebber() computes snow solid-liquid-ratio slr using the Roebber algorithm.

calslr_uutah() computes snow solid-liquid-ratio slr using the UUtah Steenburgh algorithm.

calvor() computes absolute vorticity.

fpvsnew() computes saturation vapor pressure.

tvirtual() computes virtual temperature.

Program history log:

Date Programmer Comments
2020-05-20 Jesse Meng Initial
2022-07-11 Jesse Meng CALSLR_ROEBBER
2023-02-14 Jesse Meng CALSLR_UUTAH
2023-03-22 Sam Trahan Fix out-of-bounds access by not calling BOUND
Author
Jesse Meng
Date
2020-05-20

Definition in file UPP_PHYSICS.f.

Function/Subroutine Documentation

◆ calcape()

subroutine, public upp_physics::calcape ( integer, intent(in)  itype,
real, intent(in)  dpbnd,
real, dimension(ista:iend,jsta:jend), intent(in)  p1d,
real, dimension(ista:iend,jsta:jend), intent(in)  t1d,
real, dimension(ista:iend,jsta:jend), intent(inout)  q1d,
integer, dimension(ista:iend,jsta:jend), intent(in)  l1d,
real, dimension(ista:iend,jsta:jend), intent(inout)  cape,
real, dimension(ista:iend,jsta:jend), intent(inout)  cins,
real, dimension(ista:iend,jsta:jend), intent(inout)  pparc,
real, dimension(ista:iend,jsta:jend), intent(inout)  zeql,
real, dimension(ista:iend,jsta:jend)  thund 
)

Computes CAPE and CINS.

This routine computes CAPE and CINS given temperature, pressure, and specific humidty. In "storm and cloud dynamics" (1989, academic press) cotton and anthes define CAPE (equation 9.16, p501) as

el
cape = sum g * ln(thetap/thetaa) dz
lcl
Where,
el = equilibrium level,
lcl = lifting condenstation level,
g = gravitational acceleration,
thetap = lifted parcel potential temperature,
thetaa = ambient potential temperature.
Note that the integrand ln(THETAP/THETAA) approximately
equals (THETAP-THETAA)/THETAA.  This ratio is often used
in the definition of CAPE/CINS.

Two types of CAPE/CINS can be computed by this routine.  The
summation process is the same For both cases.  What differs
is the definition of the parcel to lift.  FOR ITYPE=1 the
parcel with the warmest THETA-E in A DPBND pascal layer above
the model surface is lifted.  the arrays P1D, T1D, and Q1D
are not used.  For itype=2 the arrays P1D, T1D, and Q1D
define the parcel to lift in each column.  Both types of
CAPE/CINS may be computed in a single execution of the post
processor.

This algorithm proceeds as follows.
For each column, 
   (1)  Initialize running CAPE and CINS SUM TO 0.0
   (2)  Compute temperature and pressure at the LCL using
        look up table (PTBL).  Use either parcel that gives
        max THETAE in lowest DPBND above ground (ITYPE=1)
        or given parcel from t1D,Q1D,...(ITYPE=2).
   (3)  Compute the temp of a parcel lifted from the LCL.
        We know that the parcel's
        equivalent potential temperature (THESP) remains
        constant through this process.  we can
        compute tpar using this knowledge using look
        up table (subroutine TTBLEX).
   (4)  Find the equilibrium level.  This is defined as the
        highest positively buoyant layer.
        (If there is no positively buoyant layer, CAPE/CINS
         will be zero)
   (5)  Compute CAPE/CINS.  
        (A) Compute THETAP.  We know TPAR and P.
        (B) Compute THETAA.  We know T and P.  
   (6)  Add G*(THETAP-THETAA)*DZ to the running CAPE or CINS sum.
        (A) If THETAP > THETAA, add to the CAPE sum.
        (B) If THETAP < THETAA, add to the CINS sum.
   (7)  Are we at equilibrium level? 
        (A) If yes, stop the summation.
        (b) if no, contiunue the summation.
   (8)  Enforce limits on CAPE and CINS (i.e. no negative CAPE)
Parameters
[in]ITYPEINTEGER Flag specifying how parcel to lift is identified. See comments above.
[in]DPBNDDepth over which one searches for most unstable parcel.
[in]P1DArray of pressure of parcels to lift.
[in]T1DArray of temperature of parcels to lift.
[in]Q1DArray of specific humidity of parcels to lift.
[in]L1DArray of model level of parcels to lift.
[out]CAPEConvective available potential energy (J/kg).
[out]CINSConvective inhibition (J/kg).
[out]PPARCPressure level of parcel lifted when one searches over a particular depth to compute CAPE/CIN.
[in,out]ZEQLEquivalent level height.
THUNDThunder parameter.

Program history log:

Date Programmer Comments
1993-02-10 Russ Treadon Initial
1993-06-19 Russ Treadon Generalized routine to allow for type 2 CAPE/CINS calculations
1994-09-23 Mike Baldwin Modified to use look up tables instead of complicated equations
1994-10-13 Mike Baldwin Modified to continue CAPE/CINS calc up to at highest buoyant layer
1998-06-12 T Black Conversion from 1-D TO 2-D
1998-08-18 T Black Compute APE internally
2000-01-04 Jim Tuccillo MPI Version
2002-01-15 Mike Baldwin WRF Version
2003-08-24 G Manikin Added level of parcel being lifted as output from the routine and added the depth over which one searches for the most unstable parcel as input
2010-09-09 G Manikin Changed computation to use virtual temp added eq lvl hght and thunder parameter
2015-??-?? S Moorthi Optimization and threading
2021-07-28 W Meng Restrict computation from undefined grids
2021-09-01 E Colon Equivalent level height index for RTMA
Author
Russ Treadon W/NP2
Date
1993-02-10

Definition at line 561 of file UPP_PHYSICS.f.

◆ calcape2()

subroutine, public upp_physics::calcape2 ( integer, intent(in)  itype,
real, intent(in)  dpbnd,
real, dimension(ista:iend,jsta:jend), intent(in)  p1d,
real, dimension(ista:iend,jsta:jend), intent(in)  t1d,
real, dimension(ista:iend,jsta:jend), intent(inout)  q1d,
integer, dimension(ista:iend,jsta:jend), intent(in)  l1d,
real, dimension(ista:iend,jsta:jend), intent(inout)  cape,
real, dimension(ista:iend,jsta:jend), intent(inout)  cins,
real, dimension(ista:iend,jsta:jend), intent(inout)  lfc,
real, dimension(ista:iend,jsta:jend), intent(inout)  esrhl,
real, dimension(ista:iend,jsta:jend), intent(inout)  esrhh,
real, dimension(ista:iend,jsta:jend), intent(inout)  dcape,
real, dimension(ista:iend,jsta:jend), intent(inout)  dgld,
real, dimension(ista:iend,jsta:jend), intent(inout)  esp 
)

Computes CAPE and CINS.

This routine computes CAPE and CINS given temperature, pressure, and specific humidty. In "storm and cloud dynamics" (1989, academic press) cotton and anthes define CAPE (equation 9.16, p501) as

el
cape = sum g * ln(thetap/thetaa) dz
lcl
Where,
el = equilibrium level,
lcl = lifting condenstation level,
g = gravitational acceleration,
thetap = lifted parcel potential temperature,
thetaa = ambient potential temperature.
Note that the integrand ln(THETAP/THETAA) approximately
equals (THETAP-THETAA)/THETAA.  This ratio is often used
in the definition of CAPE/CINS.

Two types of CAPE/CINS can be computed by this routine.  The
summation process is the same For both cases.  What differs
is the definition of the parcel to lift.  FOR ITYPE=1 the
parcel with the warmest THETA-E in A DPBND pascal layer above
the model surface is lifted.  the arrays P1D, T1D, and Q1D
are not used.  For itype=2 the arrays P1D, T1D, and Q1D
define the parcel to lift in each column.  Both types of
CAPE/CINS may be computed in a single execution of the post
processor.

This algorithm proceeds as follows.
For each column, 
   (1)  Initialize running CAPE and CINS SUM TO 0.0
   (2)  Compute temperature and pressure at the LCL using
        look up table (PTBL).  Use either parcel that gives
        max THETAE in lowest DPBND above ground (ITYPE=1)
        or given parcel from t1D,Q1D,...(ITYPE=2).
   (3)  Compute the temp of a parcel lifted from the LCL.
        We know that the parcel's
        equivalent potential temperature (THESP) remains
        constant through this process.  we can
        compute tpar using this knowledge using look
        up table (subroutine TTBLEX).
   (4)  Find the equilibrium level.  This is defined as the
        highest positively buoyant layer.
        (If there is no positively buoyant layer, CAPE/CINS
         will be zero)
   (5)  Compute CAPE/CINS.  
        (A) Compute THETAP.  We know TPAR and P.
        (B) Compute THETAA.  We know T and P.  
   (6)  Add G*(THETAP-THETAA)*DZ to the running CAPE or CINS sum.
        (A) If THETAP > THETAA, add to the CAPE sum.
        (B) If THETAP < THETAA, add to the CINS sum.
   (7)  Are we at equilibrium level? 
        (A) If yes, stop the summation.
        (b) if no, contiunue the summation.
   (8)  Enforce limits on CAPE and CINS (i.e. no negative CAPE)
Parameters
[in]ITYPEINTEGER Flag specifying how parcel to lift is identified. See comments above.
[in]DPBNDDepth over which one searches for most unstable parcel.
[in]P1DArray of pressure of parcels to lift.
[in]T1DArray of temperature of parcels to lift.
[in]Q1DArray of specific humidity of parcels to lift.
[in]L1DArray of model level of parcels to lift.
[out]CAPEConvective available potential energy (J/kg).
[out]CINSConvective inhibition (J/kg).
[out]LFClevel of free convection (m).
[out]ESRHLLower bound to account for effective helicity calculation.
[out]ESRHHUpper bound to account for effective helicity calculation.
[out]DCAPEdowndraft CAPE (J/KG).
[out]DGLDDendritic growth layer depth (m).
[out]ESPEnhanced stretching potential.

Program history log:

Date Programmer Comments
1993-02-10 Russ Treadon Initial
1993-06-19 Russ Treadon Generalized routine to allow for type 2 CAPE/CINS calculations
1994-09-23 Mike Baldwin Modified to use look up tables instead of complicated equations
1994-10-13 Mike Baldwin Modified to continue CAPE/CINS calc up to at highest buoyant layer
1998-06-12 T Black Conversion from 1-D TO 2-D
1998-08-18 T Black Compute APE internally
2000-01-04 Jim Tuccillo MPI Version
2002-01-15 Mike Baldwin WRF Version
2003-08-24 G Manikin Added level of parcel being lifted as output from the routine and added the depth over which one searches for the most unstable parcel as input
2010-09-09 G Manikin Changed computation to use virtual temp added eq lvl hght and thunder parameter
2015-??-?? S Moorthi Optimization and threading
2021-09-03 J Meng Modified to add 0-3km CAPE/CINS, LFC, effective helicity, downdraft CAPE, dendritic growth layer depth, ESP
2021-09-01 E Colon Equivalent level height index for RTMA
2022-08-27 S Trahan Fixed bug in CALCAPE2 where extreme atmospheric conditions cause an out-of-bounds access
2022-09-01 S Trahan Fixed another bug in CALCAPE2 where extreme atmospheric conditions cause an out-of-bounds access
Author
Russ Treadon W/NP2
Date
1993-02-10

Definition at line 1039 of file UPP_PHYSICS.f.

References exch().

◆ caldiv()

subroutine, public upp_physics::caldiv ( real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u,lm), intent(in)  uwnd,
real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u,lm), intent(in)  vwnd,
real, dimension(ista:iend,jsta:jend,lm), intent(inout)  div 
)

Computes divergence.

For GFS, this routine copmutes the horizontal divergence using 2nd-order centered scheme on a lat-lon grid

Parameters
[in]UWNDU wind (m/s) mass-points.
[in]VWNDV wind (m/s) mass-points.
[out]DIVdivergence (1/s) mass-points.

Program history log:

Date Programmer Comments
2016-05-05 Sajal Kar Modified CALVORT to compute divergence from wind components
2016-07-22 S Moorthi Modified polar divergence calculation
Author
Sajal Kar W/NP2
Date
2016-05-05

Definition at line 2164 of file UPP_PHYSICS.f.

References exch(), and fullpole().

◆ calgradps()

subroutine, public upp_physics::calgradps ( real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u), intent(in)  ps,
real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u), intent(inout)  psx,
real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u), intent(inout)  psy 
)

Computes gradients of a scalar field PS or LNPS.

For GFS, this routine computes horizontal gradients of PS or LNPS. Using 2nd-order centered scheme on a lat-lon grid.

Parameters
[in]PSSurface pressure (Pa) mass-points.
[out]PSXZonal gradient of PS at mass-points.
[out]PSYMeridional gradient of PS at mass-points.

Program history log:

Date Programmer Comments
2016-05-05 Sajal Kar Reduced from CALVORT to zonal and meridional gradients of given surface pressure PS, or LNPS
Author
Sajal Kar W/NP2
Date
2016-05-05

Definition at line 2442 of file UPP_PHYSICS.f.

References exch().

◆ calrh()

subroutine, public upp_physics::calrh ( real, dimension(ista:iend,jsta:jend), intent(in)  p1,
real, dimension(ista:iend,jsta:jend), intent(in)  t1,
real, dimension(ista:iend,jsta:jend), intent(inout)  q1,
real, dimension(ista:iend,jsta:jend), intent(out)  rh 
)

Computes relative humidity.

Parameters
[in]P1real Pressure (Pa)
[in]T1real Temperature (K)
[in,out]Q1real Specific humidity (kg/kg)
Note
P1/T1/Q1 refer to pressure/temperature/specific humidity at the selected /requested level. For example, if the user wants relative humidity at 500mb, then T/P/Q at 500mb is input into the call to output RH@500mb
Parameters
[out]RHreal Relative humidity (decimal form)

Definition at line 70 of file UPP_PHYSICS.f.

◆ calrh_gfs()

subroutine, public upp_physics::calrh_gfs ( real, dimension(ista:iend,jsta:jend), intent(in)  p1,
real, dimension(ista:iend,jsta:jend), intent(in)  t1,
real, dimension(ista:iend,jsta:jend), intent(inout)  q1,
real, dimension(ista:iend,jsta:jend), intent(inout)  rh 
)

Computes relative humidity.

This routine computes relative humidity given pressure, temperature, specific humidity. an upper and lower bound of 100 and 1 percent relative humidity is enforced. When these bounds are applied the passed specific humidity array is adjusted as necessary to produce the set relative humidity.

Parameters
[in]P1Pressure (pa)
[in]T1Temperature (K)
[in,out]Q1Specific humidity (kg/kg)
[out]RHRelative humidity (decimal form)

Program History Log

Date Programmer Comments
????-??-?? DENNIS DEAVEN Initial
1992-12-22 Russ Treadon Modified as described above
1998-06-08 T Black Conversion from 1-D to 2-D
1998-08-18 Mike Baldwin Modify to compute RH over ice as in model
1998-12-16 Geoff Manikin undo RH computation over ice
2000-01-04 Jim Tuccillo MPI Version
2002-06-11 Mike Baldwin WRF Version
2013-08-13 S. Moorthi Threading
2006-03-19 Wen Meng Modify top pressure to 1 pa
Author
Russ Treadon W/NP2
Date
1992-12-22

Definition at line 195 of file UPP_PHYSICS.f.

◆ calrh_gsd()

subroutine, public upp_physics::calrh_gsd ( real, dimension(ista:iend,jsta:jend)  p1,
real, dimension(ista:iend,jsta:jend)  t1,
real, dimension(ista:iend,jsta:jend)  q1,
real, dimension(ista:iend,jsta:jend)  rhb 
)

Compute RH with the NOAA GSL (formerly NOAA GSD) algorithm used for RUC and Rapid Refresh.

Parameters
P1real Pressure (Pa)
T1real Temperature (K)
Q1real Specific humidity (kg/kg)
Note
P1/T1/Q1 refer to pressure/temperature/specific humidity at the selected /requested level. For example, if the user wants relative humidity at 500mb, then T/P/Q at 500mb is input into the call to output RH@500mb
Parameters
RHBreal Relative humidity (decimal form)

Definition at line 268 of file UPP_PHYSICS.f.

◆ calrh_nam()

subroutine, public upp_physics::calrh_nam ( real, dimension(ista:iend,jsta:jend), intent(in)  p1,
real, dimension(ista:iend,jsta:jend), intent(in)  t1,
real, dimension(ista:iend,jsta:jend), intent(inout)  q1,
real, dimension(ista:iend,jsta:jend), intent(out)  rh 
)

Computes relative humidity.

This routine computes relative humidity given pressure, temperature, specific humidity. an upper and lower bound of 100 and 1 percent relative humidity is enforced. When these bounds are applied the passed specific humidity array is adjusted as necessary to produce the set relative humidity.

Parameters
[in]P1Pressure (pa)
[in]T1Temperature (K)
[in,out]Q1Specific humidity (kg/kg)
[out]RHRelative humidity (decimal form)

Program History Log

Date Programmer Comments
????-??-?? DENNIS DEAVEN Initial
1992-12-22 Russ Treadon Modified as described above
1998-06-08 T Black Conversion from 1-D to 2-D
1998-08-18 Mike Baldwin Modify to compute RH over ice as in model
1998-12-16 Geoff Manikin undo RH computation over ice
2000-01-04 Jim Tuccillo MPI Version
2002-06-11 Mike Baldwin WRF Version
2006-03-19 Wen Meng Modify top pressure to 1 pa
Author
Russ Treadon W/NP2
Date
1992-12-22

Definition at line 116 of file UPP_PHYSICS.f.

◆ calrh_pw()

subroutine, public upp_physics::calrh_pw ( real, dimension(ista:iend,jsta:jend)  rhpw)

Algorithm used at GSL for RUC and Rapid Refresh.

Parameters
RHPWreal Relative humidity with respect to precipitable water (entire atmosphere)

Definition at line 310 of file UPP_PHYSICS.f.

◆ calslr_roebber()

subroutine, public upp_physics::calslr_roebber ( real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u,lsm), intent(in)  tprs,
real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u,lsm), intent(in)  rhprs,
real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u), intent(out)  slr 
)

Computes snow solid-liquid-ratio (SLR) using the Roebber algorithm.

Obtained the code and data from WPC. WPC's SLR products include SLR computed from GFS and NAM, SLR climotology, and averaged SLR. UPP computes SLR for GFS and RRFS. SLR climatology is not used in UPP calculation but the data is saved in fix directory for reference. Breadboard coefficients are included in this module to enhance the performance. Original Breadboard coefficients files are also saved in fix directory.

Parameters
[in]tprsreal Temperature on pressure levels.
[in]rhprsreal Relative humidity on pressure levels.
[out]slrreal Solid snow to liquid ratio.

Program history log:

Date Programmer Comments
2022-07-11 Jesse Meng Initial

2023-01-06 | Jesse Meng ! Import Breadboard coefficients into module

Author
Jesse Meng
Date
2022-07-11

Definition at line 2678 of file UPP_PHYSICS.f.

◆ calslr_uutah()

subroutine, public upp_physics::calslr_uutah ( real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u), intent(out)  slr)

Computes snow solid-liquid-ratio slr using the Steenburgh algorithm.

Obtained the code and data from U of Utah Jim Steenburgh and Peter Veals. SLR = m1X1 + m2X2 + m3X3 + m4X4 + m5X5 + m6X6 + b.

 X1 = wind speed at at 1km above ground level (AGL) in m/s
 m1 = -0.174848

 X2 = temperature at 2km AGL in Kelvin
 m2 = -0.52644

 X3 = wind speed at 2 km AGL in m/s
 m3 = 0.034911

 X4 = wind speed at 500 m AGL in m/s
 m4 = -0.270473

 X5 = temperature at 1 km AGL in Kelvin
 m5 = 0.028299

 X6 = temperature at 500 m AGL in m/s
 m6 = 0.096273

 b =  118.35844
Parameters
[out]SLRreal Solid snow to liquid ratio

Program history log:

Date Programmer Comments
2023-01-23 Jesse Meng Initial
Author
Jesse Meng
Date
2023-01-03

Definition at line 4347 of file UPP_PHYSICS.f.

◆ calvor()

subroutine, public upp_physics::calvor ( real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u), intent(in)  uwnd,
real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u), intent(in)  vwnd,
real, dimension(ista_2l:iend_2u,jsta_2l:jend_2u), intent(inout)  absv 
)

Computes absolute vorticity.

Parameters
[in]UWNDU wind (m/s) mass-points.
[in]VWNDV wind (m/s) mass-points.
[out]ABSVabsolute vorticity (1/s) mass-points.

Program history log:

Date Programmer Comments
1992-12-22 Russ Treadon Initial
1998-06-08 T Black Convesion from 1-D to 2-D
2000-01-04 Jim Tuccillo MPI Version
2002-01-15 Mike Baldwin WRF Version C-grid
2005-03-01 H Chuang Add NMM E grid
2005-05-17 H Chuang Add Potential vorticity calculation
2005-07-07 B Zhou Add RSM in computing DVDX, DUDY and UAVG
2013-08-09 S Moorthi Optimize the vorticity loop including threading
2016-08-05 S Moorthi add zonal filetering
2019-10-17 Y Mao Skip calculation when U/V is SPVAL
2020-11-06 J Meng Use UPP_MATH Module
2022-05-26 H Chuang Use GSL approach for FV3R
Author
Russ Treadon W/NP2
Date
1992-12-22

Definition at line 1743 of file UPP_PHYSICS.f.

References exch(), and fullpole().

◆ fpvsnew()

elemental real function, public upp_physics::fpvsnew ( real, intent(in)  t)

Computes saturation vapor pressure.

Compute saturation vapor pressure from the temperature. A linear interpolation is done between values in a lookup table computed in gpvs. See documentation for fpvsx for details. Input values outside table range are reset to table extrema. The interpolation accuracy is almost 6 decimal places. On the Cray, fpvs is about 4 times faster than exact calculation. This function should be expanded inline in the calling routine.

Parameters
[in]tReal(krealfp) Temperature in Kelvin.
Returns
fpvsnew Real(krealfp) Saturation vapor pressure in Pascals.

Program history log:

Date Programmer Comments
1991-05-07 Iredell Initial. Made into inlinable function
1994-12-30 Iredell Expand table
1999-03-01 Iredell F90 module
2001-02-26 Iredell Ice phase
Author
N Phillips w/NMC2X2
Date
1982-12-30

Definition at line 398 of file UPP_PHYSICS.f.

◆ tvirtual()

elemental real function, public upp_physics::tvirtual ( real, intent(in)  t,
real, intent(in)  q 
)

Computes virtual temperature.

Returns
Parameters
[in]Treal Temperature
[in]Qreal Specific humidity
Returns
virtual temperature

Definition at line 1701 of file UPP_PHYSICS.f.