Water Vapor - Other Formulas

My programs implement a number of formulas which claim to compute saturated water pressure at a given temperature. However, from time to time I run across new formulas - and occasionally, flat out wrong formulas. Since I don't want to be continuously updating my program, I will simply add these "new" formulas here.

The Engineering ToolBox

The Engineering ToolBox is a source of fairly reliable data. As a result, I was shocked (and disappointed) to find that their Saturation Pressure of Water Vapor formula is totally different than all the other formulas I've collected.

p_ws = e^{(77.3450 + 0.0057 T - 7235 / T)} / T^8.2

T in kelvins

In general, the results agree fairly well with GoffGratch up to about 0°C, then it begins to strongly diverge. For instance, (based on this table)

 T (C)  GoffGratch (mb)   EngToolBox (mb)
 -15        1.91              1.9 
 - 9        3.09              3.0
   2        7.05              6.9
  25       31.652            31.30
  46      100.84            101.4

Only GISS-AR4 (NASA) is worse.

Unfortunately, they do not reference their source for this equation - it is merely presented as "a statement of fact".

A Google search for that equation returned numerous hits. Unfortunately, all the examples I checked are obviously just quoting The Engineering ToolBox - either an explicit reference or the identical text and formatting.

In addition, on their Dry Bulb, Wet Bulb and Dew Point Temperatures page they say

The SI unit is Kelvin (K). Zero Kelvin equals to -273°C.

which is, of course, wrong! It should be -273.15°C. (Yes - it does make a difference! Particularly in a reference book.)

NASA Tech Note D-8401

Equations for the Determination of Humidity from Dewpoint and Psychrometric Data - NASA Tech Note D-8401, Parish and Putnam, January 1977 - provides yet another set of equations with derivations. They claim that the following is the Magnus formula.

log (e_SW) = -4.9283 Log T - 2937.4/T + 22.5518

which has the general form

e_SW = 10^{(c + b/T)} T^a   where e_SW is the saturation water vapor pressure

None of the other formulas I have seen to date have that form. The MagnusTetens formula is

log(P) = A + B*T/(C+T)

To be clear, several formulas contain T^a, but a is usually positive and there are always a lot more terms.

Their derivations of the constants - a, b, c - are provided in Appendix A. They always use 273 to convert between K and C (should be 273.15) and they use 2 vapor pressure calibration points - both are at 0°C (32°F).

Over water	Over Ice
6.11 mb	6.107 mb

Appendix A, Example 2 computes relative humidity from a dew point below freezing. Their value is 17.3%, my program, using MurphyKoop, computes 19.2% over water and 17.1% over ice. Example 3 uses 2 temperatures above freezing - they compute 84.6%, my program finds the same for most of the provided algorithms.

In addition, their definition of relative humidity (based on that giving in the Smithsonian Meteorological Tables, 1971) is a bit strange.

It is the ratio of the actual partial pressure of water vapor over ice or water (depending on whether the air temperature is above or below freezing) divided by the saturation pressure over water.

I have not seen any other papers specify that the denominator always be with respect to the saturation pressure over liquid water.

New Magnus Coefficients

Improved Magnus Form Approximation of Saturation Vapor Pressure, Alduchov (Ru) and Eskridge (US), April 1996, Journal of Applied Meteorology

This paper argues that the standard Magnus coefficients should be replaced so that the Magnus equation can be used instead of the much more complex equations of Wexler (1976), Sonntag (1990), and Goff and Gratch (1946).

As an aside, these dates do not match those provided by Vömel - Wexler (1983), and Sonntag (1994) - and the MurphyKoop (2005) formulation was published after this paper.

Part of the justification of this approach is that atmospheric temperature is reported to only 0.1°C accuracy world wide and only 0.5°C in the US. That, plus the general inaccuracies of the reported humidity proxies (such as dew point) means that the true accuracy of any computation is low enough that several formulations produce calculated saturation pressures that are "the same" within some margin of error.

This is the general form of the Magnus formula.

P = C e^AT/(B+T)  where T is in °C

Alduchov provides 6 different sets of Magnus coefficients by others, and one "better formulation" suggested by himself. The following table only shows 3 - MagnusTetens from Vömel, one of the six "other" formulations, and the new "suggested" Alduchov formulation.

Constant	MagnusTetens 1967		Sonntag 1990		Alduchov 1996
Over	Water	Ice	Water	Ice	Water	Ice
A	17.269388	21.8745584	17.62	22.46	17.625	22.587
B	237.3	265.5	243.12	272.62	243.04	273.86
C	6.1078	6.1078	6.112	6.112	6.1094	6.1121

MagnusTetens (1967) via Vömel

	e_w = 6.1078 e^{17.269388 *}^{(T-273.16) / (T - 35.86)}		Over water
	e_i = 6.1078 e^{21.8745584 * (T-273.16) / (T - 7.66)}		Over ice
	with T in [K] and e_w & e_i in [hPa]

Note that MagnusTetens (1967) subtracts 273.16 meaning that the temperature is not actually Celsius and that, therefore, a direct comparison of the coefficients in a table like this is not possible.

Full equations

The Alduchov paper discusses 4 "complex" equations by 3 authors - its purpose was to replace these "difficult to compute" formulations with a simpler Magnus form. (With modern computers, this goal is questionable.)

The full Sonntag equation is identical in both sources (and different from Sonntag's Magnus approximation used in the table above).

The Goff Gratch equation differs between Vömel and this paper in the fourth term.

Vömel - Goff Gratch 1946, 1984
Log10 ew =  -7.90298 (373.16/T-1) 
           + 5.02808 Log10(373.16/T) 
           - 1.3816 10^-7 (10^{11.344 (1-T/373.16)}  -1) 
           + 8.1328 10^-3 (10^{-3.49149 (373.16/T-1)}  -1) 
           + Log10(1013.246) 

Alduchov - Goff Gratch 1945
Log10 ew =   7.90298 (1 - Ts/T) 
           + 5.02808 Log10(Ts/T) 
           + 1.3816 10^-7 (1 - 10^{11.344 (1-T/Ts)}) 
           - 8.1328 10^-3 (1 - 10^{-3.49149 (1-Ts/T)}) 
           + Log10(1013.246) 

Where Ts = 373.16 K  the steam point at 1 atm

This is the problem in the 4th term.

10^{-3.49149 (Ts/T - 1)}  vs 10^{-3.49149 (1 - Ts/T)}

There is no way to know if this is just a typo or if it invalidates the paper's conclusions.

The 1946 Goff Gratch equation is identical to the equation identified by Vömel as WMO (Goff, 1957; WMO, 2012). Note that the 1946 equation is based on the steam point (373.16 K), and that the 1957 equation is based on the triple point (273.16 K).
The Wexler (1976) equation is identical to the equation provided in New Equations for Computing Vapor Pressure and Enhancement Factor, Arden L.Buck, 1981 (see the next section) - except that most for the coefficients were normalized with one digit left of the decimal point.

Wexler

New Equations for Computing Vapor Pressure and Enhancement Factor, Arden L.Buck, 1981

This paper contains a good overview of the progress toward developing equations to compute saturation vapor pressure. It includes a modification to the Magnus equation by Bögel and 2 equations by Wexler.

Wexler

	e_w = 0.01 exp[-2991.2729 T^-2 - 6017.0128 T^-1 + 18.87643854 - 0.028354721 T + 0.17838301 10^-4 T² - 0.84150417 10^-9 T³ + 0.44412543 10^-12 T⁴ + 2.858487 ln(T)]		1976 - Over water
	e_i = 0.01 exp[- 5865.3696 T^-1 + 22.241033 + 0.013749042 T - 0.34031775 10^-4 T² + 0.26967687 10^-7 T³ + 0.6918651 ln(T)]		1977 - Over ice
	with T in [K] and e_w & e_i in [hPa]

washington.edu skewt.js

While developing a skew-T charting application, I found washington.edu skewt.js. This is their saturation vapor pressure subroutine (no source given)

function ESAT(t) {
  // Computes the saturation vapour pressure over water at temperature t.
  // ESAT in mb, t in K. Log to base 10 is needed for this function.
  // Dim a0 As Double, a1 As Double, a2 As Double
  a0 = 23.832241 - 5.02808 * Log10(t);
  a1 = 0.00000013816 * 10 ^ (11.344 - 0.0303998 * t);
  a2 = 0.0081328 * 10 ^ (3.49149 - 1302.8844 / t);
  ESAT = 10 ^ (a0 - a1 + a2 - 2949.076 / t);
  return ESAT;
}

This is almost identical to the Goff Gratch equation (only one constant is different), but with the terms multiplied out (which presumably makes the computation faster) which, unfortunately, obscures the use of the triple point temperature.

Goff Gratch equation
(Smithsonian Tables, 1984, after Goff and Gratch, 1946):

    Log10 ew =  -7.90298 (373.16/T-1)
                       + 5.02808 Log10(373.16/T)
                       - 1.3816 10^-7 (10^[11.344   (1-T/373.16)] -1)
                       + 8.1328 10^-3 (10^[-3.49149 (373.16/T-1)] -1)
                       + Log10(1013.246)
    with T in [K] and ew in [hPa]
voemel

The following computations were performed using a spreadsheet.

 7.90298 * 373.16        = 2949.0760168
 5.02808 * Log10(373.16) = 12.931694236
11.344   / 373.16        = 0.0303998285
 3.49149 * 373.16        = 1302.8844084
Log10(1013.246)          = 3.0057148979

7.90298 + 12.931694236 + 1.3816*10^-7 - 8.1328*10^-3 + 3.0057148979 = 23.8322564721

This is the only significant difference - 23.832241 vs 23.832256

A few test computations show the same results as Goff Gratch to 2 decimal places - 4 significant figures.

Formula	-20 C	0 C	15 C	30 C
Goff Gratch	1.25292	6.10336	17.03281	42.40599
washington.edu	1.252880448	6.1031443925	17.0322055233	42.4044797454
MurphyKoop	1.25504	6.11213	17.0588	42.46814

MODTRAN

MODTRAN is a binned-spectrum transmittance program used by the US Air Force. The code for version 3 (dated April 1985, the current version is 6, but I don't have the code) uses

F(A) is saturated water vapor density at temp T 

A    = 273.15/T

F(A) = EXP(18.9766 - 14.9595*A - 2.43882*A*A) * A

They claim that the formula "is accurate to better than 1 percent from -50 to +50 deg C."

However, since this formula computes density (gm/m³), I need to use the gas law to get the pressure.

PV = nRT

P = ρ/molWt RT

Where ρ is F(A) in the program

Where the molecular weight of water (molWt) used in the program is 18.015 gm/mol.

Via the engineeringtoolbox, the gas constant is

R = 8314.47 J/(kmol K)

And the decimal point needs to be adjusted.

PV = nRT
P = ρ/molWt RT

P N/m² = ρ (gm/m³)/[molWt (gm/mol)] * (R J/(kmol K)     * T K) * (kmol/1000 mol)
       = ρ (gm/m³)/[molWt (gm/mol)] * (R (N m)/(kmol K) * T K) * (kmol/1000 mol)
       = 10^-3 N/m²

1 bar/1000 mbar * 10^5 pascal/bar * 1 N/m²/1 pascal = 10² N/m²/mbar

P mbar = 10^-3 N/m² / (10² N/m²/mbar)
       = ρ/molWt RT * 10^-5 mbar

Assuming T is in kelvin (it works), this table compares the adjusted MODTRAN formula with the MurphyKoop formulation as computed using my calculator.

Temperature		MurphyKoop mbar	MODTRAN
K	C	MurphyKoop mbar	g/m3	mbar
223.15	-50	0.06317	0.0616191052	0.0634618288
273.15	0	6.11213	4.8466124646	6.109980215
300	26.85	35.36764	25.5441184444	35.3681609463
323.15	50	123.51983	83.1908557132	124.0737921266

So, their formula is pretty good, but it does not exactly match any of the other formulations (models) in the calculator.

They do not provide a formula for finding the relative humidity over ice.

As a part of my normal research, I performed a google search for

EXP(18.9766 - 14.9595*A - 2.43882*A*A)

There were a few hits.

EVALUATION OF ATMOSPHERIC EFFECTS FOR OPERATIONAL TACTICAL DECISION AID 1988 (pdf): An Air Force report - Refers to LOTRAN-6 and CTRAN
The formula is "an empirical function for saturated water vapor density (g/m³) at temperature T"
This was very useful since it clearly stated that the MODTRAN formula computed density in g/m³ - the source code comments just said "density" and gave no units.
Calculating and converting between common water vapour measures in meteorological data and their use in support of earth observation validation 2003 (pdf): CSIRO Earth Observation Centre (EOC)
This presents the formula, but without providing a reference
Physics of the Greenhouse Effect and Convection in Warm Oceans 1994 (pdf, not searchable): Equation 1 has the same coefficients as the formula above. It is adapted from Kneizys et al. (1988) - Users Guide to LOTRAN 7
With respect to relative humidity recorded by radio sondes, "errors are seen to range from about 6% near the surface to nearly 13% at about the 200-mb level."

I checked a number of additional references and they all reference that equation to LOTRAN.

Vaisala

In Humidity Conversion Formulas, Vaisala suggests using the IAPWS formula over liquid water and the following over ice

Ps = 6.11657 exp(-13.928169 [1 - (273.16/T)^1.5 ] + 34.707823[1 - (273.16/T)^1.25 ] ) // over ice -100 to 0.01°C

However, the voemel reference says

The Hyland and Wexler formulation is used by Vaisala and is very similar to the formula by Sonntag.

Lowe

Paul R. Lowe wrote at least 2 papers on attempting to find faster algorithms for computing the water vapor saturation pressure. Basically, the accepted standard at the time was Goff-Gratch 1946 which requires multiple logs and powers. Several other formulations also existed then, but they all required at least one exponential. Even today, all the other algorithms I have seen require at least one exponential.

In The Computation of Saturation Vapor Pressure (March 1974), Lowe argues for finding an algorithm with acceptable accuracy and a fast computation speed. The paper presents two 6th order polynomials for computing saturated water vapor over ice and liquid water that, he claims, run faster and are more accurate than Tetens.

In An Approximating Polynomial for the Computation of Saturated Vapor Pressure (1976), he repeats the polynomials and provides additional data supporting his claims.

These are just the coefficients for temperature in C and over water. The papers also contain coefficients for the temperature in kelvin over water and for over ice.

E_s = a0 + T(a1 + T(a2 + T(a3 + T(a4 + T(a5 + T(a6) )))))
a0 = 6.107799961
a1 = 4.436518521 E-1
a2 = 1.428945805 E-2
a3 = 2.650648471 E-4
a4 = 3.031240396 E-6
a5 = 2.034080948 E-8
a6 = 6.136820929 E-11

Range of validity: -50°C to +50°C over liquid water

Accuracy

In 1976, to determine the accuracy of his polynomial, Lowe selected Goff-Gratch 1946 as the gold standard. According to New Equations for Computing Vapor Pressure and Enhancement Factor (Buck 1981), Wexler 1976 is a more accurate formulation.

In Table 1, Lowe 1976 claims that Goff-Gratch computes the saturation pressure as 6.1078 mb at 0°C, but, according to my calculator, that is the pressure computed at 0.01°C (the triple point) and the Goff-Gratch value at 0.00°C is actually 6.10336 mb. Via a0 in the coefficients above, Lowe's formula sets the saturation pressure to 6.107799961 mb at 0.00°C. (The current standard is 6.11657 ± 0.0001 mbar at 273.16 K (0.01°C) - the triple point.)

I checked a number of Goff-Gratch and Tetens values in his table with those computed by my calculator - they do not match. One of us has an obvious math error.

Compared to voemel, there is only a slight difference in one of the Lowe 1976 Tetens formula coefficients, but a fairly significant difference in the values shown in Table 1.

Psat = 6.1078 * EXP(17.2693882*(Temp-273.16) /(Temp-35.86)) Temp in K Murray 1967

Psat = 6.1078 * EXP(17.269388 * (T_K-273.16) / (T_K-35.86))  T_K in K Voemel

However, the computed and reported values differ just a bit too much

Values for Tetens formula
Temp	Voemel	Lowe
-25	0.79856	0.79928
0	6.10336	6.1078
10	12.27069	12.2789
20	23.36647	23.3810
30	42.40202	42.4254

This is bad! Using 20°C as an example

exp( [ ln(23.3664700) / 17.2693880 ] * 17.2693882 ) = 23.3664709

In my opinion, the table values and those computed via my calculator should have been identical ... and they are not !!

Whether or not these errors are significant depends on your application. But, in a peer-reviewed paper about algorithm accuracy, I expect better. On the other hand, this is an old paper and some of the standards and/or accepted calibration points have changed since then. As a result, his coefficients may need to be tweaked anyway.

At any rate, the Lowe paper is important since it is the main reference on the wikipedia Vapour pressure of water page and since these are the only formulations I've seen that use simple polynomials.

At some point, I want to compare the computation speeds and accuracies of all the formulations I have access to, in a manner similar to Lowe, but using MurphyKoop or HylandWexler as the accuracy standard.

Additional Formulas

Lowe 1976 mentions 3 other formulas I have not seen previously

2 versions of Tabata (1973)

log10(E_s) = 9.28603523 - 2.32237885(1000/T)

log10(E_s) = 8.42926609 - 1.82717843(1000/T) - 0.071208271(1000/T)²

T is in kelvins

Richards formula (from Wigley 1974) (Richards April 1971 article - paywalled)

E_s = 1013.25 exp(13.3185 S - 1.9760 S² - 0.6445 S³ - 0.1299 S⁴)

S = 1 - (373.16/T)
T is in kelvins
373.16 K is identified as the steam point

A number of formulations, including Goff-Gratch, identify 373.16 K as the steam point. However, before 1990 the steam point was actually 373.15 K, after 1990, it became 373.124 K.

The Richards 1971 formula is identical to the one used in Atmospheric Chemistry and Physics by Seinfeld and Pandis with the exception that they used 373.15 K and not 373.16 K.

CRC

The CRC Handbook of Chemistry and Physics (85th ed. 2004) has tables of vapor pressures over ice from 0.01°C to -80°C (page 6-9) and over water from 0°C to 370°C (pages 6-10 & 6-11). (The number of displayed pages is not predictable - if it is hidden, scroll up a few pages and then you can scroll down and see page 6-11.)

There is no formula for finding the values - just tables.

The values in the "over water" table are from the NBS/NRC Steam Tables, 1984 and do not agree with any of the formulas used in my main program. I have not checked the values against those produced by formulas on this page.

The "over ice" values closely track MurphyKoop and HylandWexler from 0.01°C to -40°C and then they begin to diverge and follow Sonntag from -50°C to -80°C.

The tables do not cover "over water" below freezing.

Author: Robert Clemenzi