Reappraising the appropriate calculation of a common meteorological quantity: potential temperature

The potential temperature is a widely used quantity in atmospheric science since it is conserved for dry air's adiabatic changes of state. Its definition involves the specific heat capacity of dry air, which is traditionally assumed as constant. However, the literature provides different values of this allegedly constant parameter, which are reviewed and discussed in this study. Furthermore, we derive the potential temperature for a temperature-dependent parameterisation of the specific heat capacity of dry air, thus providing a new reference potential temperature with a more rigorous basis. This new reference shows different values and vertical gradients, in particular in the stratosphere and above, compared to the potential temperature that assumes constant heat capacity. The application of the new reference potential temperature is discussed for computations of the Brunt–Väisälä frequency, Ertel's potential vorticity, diabatic heating rates, and for the vertical sorting of observational data.

1 Introduction

According to the book Thermodynamics of the Atmosphere by Alfred Wegener (1911), the first published use of the expression potential temperature in meteorology is credited to Wladimir Köppen (1888)¹ and Wilhelm von Bezold (1888), both following the conclusions of Hermann von Helmholtz (1888) (see also Kutzbach, 2016). Even prior to the introduction of the entropy, Poisson (1833) and Thomson (1862) used the “adiabatic equation”, the basis of what is understood today as “potential temperature”², to describe adiabatic processes, e.g. the coincident variation of temperature and pressure on the movement of air, which is “independent of the effects produced by the radiation or conduction of heat” (Thomson, 1862)³. Approximately 26 years later, von Helmholtz perceived that within the atmosphere the heat exchange between air masses of different temperatures, which are relatively moved, is insufficiently explained by heat transfer due only to radiation and convection. He argued that wind phenomena (e.g. the trade winds), storm events, and the atmospheric circulation were more intense, of larger extent, and more persistent than observed if the air's heat exchange within the discontinuity region (the friction surface of the different air masses) was not mainly due to eddy-driven mixing. On his way to analytically describe the heat exchange of different air masses within the atmosphere, in May 1880, von Helmholtz introduced the air's immanent heat while its absolute temperature changes with changing pressure (von Helmholtz, 1888). In essence, von Helmholtz concluded that the temperature gained by a volume of dry air due to its adiabatic descent from a certain initial pressure level (p) to ground pressure (p₀) corresponds to the air's immanent heat. In November of the same year, in agreement with von Helmholtz and probably inspired by a presentation that was given in June by Köppen (1888), this property was renamed and reintroduced as the air's potential temperature (θ in the following) by von Bezold (1888) with the following definition for strictly adiabatic changes of state:

$$\begin{array}{}\text{(1)}& \mathit{\theta }=T{\left({\displaystyle \frac{p_{\mathrm{0}}}{p}}\right)}^{\frac{\mathit{\gamma }-\mathrm{1}}{\mathit{\gamma }}},\end{array}$$

where T and p are the absolute temperature and pressure, respectively, of an air parcel at a certain initial (pressure) altitude level. The quantities θ and p₀ are corresponding values of the same air parcel's absolute temperature and pressure if the air was exposed to conditions at ground level. The dimensionless coefficient γ, nowadays called the isentropic exponent, was specified as 1.41 (von Bezold, 1888).

Moreover, in the same publication, von Bezold concluded that for moist air's adiabatic changes of state, its potential temperature remains unchanged as long as the change of state occurs within dry-adiabatic limits; and further, if there is condensation and precipitation, the potential temperature changes by a magnitude that is determined by the amount of water that falls out of the air parcel. From a modern perspective, it is clear that the air parcel is an isolated thermodynamic system, and adiabatic processes correspond to processes with conserved entropy (i.e. isentropic processes). The description of the immanent heat is then equivalent to the thermodynamic state function entropy, which corresponds to potential temperature of dry air in a one-to-one relationship.

In general, the potential temperature has the benefit of providing a practicable vertical coordinate (equivalent to the pressure level or the altitude above, e.g. sea level) to visualise and analyse the vertical distribution and variability of (measured) data related to any type of atmospheric parameter. Admittedly, the use of the potential temperature as a vertical coordinate is initially less intuitive than applying altitude or pressure coordinates. Indeed, the potential temperature bears a certain abstractness to describe an air parcel's state at a certain altitude level by its imaginary dry-adiabatic descent to ground conditions. However, one major advantage of using the potential temperature as a vertical coordinate is that the (measured) data are sortable with respect to the entropy state at which the atmospheric samples were taken. Thus, comparing repeated measurements of an atmospheric parameter on an isentropic surface or layer excludes any diabatic change in the probed air mass.

Apart from characterising the isentropes, the vertical profiles of the potential temperature (θ as a function of height z) are used as the reference for evaluating the atmosphere's actual vertical temperature gradient, which allows characterising its static stability. Notably, von Bezold (1888) already proposed the potential temperature as an atmospheric stability criterion. In its basic formulation, the potential temperature exclusively refers to the state of dry air, and thus the potential temperature characterises the atmosphere's static stability with respect to vertical displacements of a dry air parcel. In meteorology, the static stability parameter is expressed in terms of the (squared) Brunt–Väisälä frequency N, often written in the form

$$\begin{array}{}\text{(2)}& N^{\mathrm{2}}={\displaystyle \frac{g}{\mathit{\theta }}}{\displaystyle \frac{\partial \mathit{\theta }}{\partial z}},\end{array}$$

where $g=\mathrm{9.81}\mathrm{m}{\mathrm{s}}^{-\mathrm{2}}$ is the gravitational acceleration. The potential temperature twice enters the formulation of the stability parameter, as the denominator (θ⁻¹) and as the vertical gradient $\frac{\partial \mathit{\theta }}{\partial z}$. In the research field of dynamical meteorology, the potential vorticity (PV) is often used (Ertel, 1942; Hoskins et al., 1985; Schubert et al., 2004). The PV is proportional to the scalar product of the atmosphere's vorticity (the air's local spinning motion) and its stratification (the air's tendency to spread in layers of diminished exchange). More concretely, the PV is the scalar product of the absolute vorticity vector and the three-dimensional gradient of θ, i.e. not only the potential temperature's vertical gradient but also its partial derivatives on the horizontal plane add to the resulting PV, although, particularly at stratospheric altitudes, the vertical gradient constitutes the dominant contribution. For the analytical description of a fluid's motion within a rotational system, as is the atmosphere, the PV provides a quantity that varies exclusively due to diabatic processes. Frequently, the PV is used to define the tropopause height (usually at 2 PV units; see e.g. Gettelman et al., 2011) or the edge of a large-scale cyclone such as the polar winter vortex on specific θ levels (cf. Curtius et al., 2005).

While for a dry atmosphere (i.e. with little or no water vapour) the potential temperature is the correct conserved quantity (corresponding to entropy) for reversible processes, for an atmosphere containing water in two or more phases (vapour, liquid, and/or solid phases) energy transfers due to phase changes play a major role. Thus, the formulation of the potential temperature has to be extended since entropy is still the right quantity for reversible processes, including phase changes. Starting from the equation for the moist specific entropy, as derived from the first law of thermodynamics and the Gibbs equation, further extensions of the dry-air potential temperature have been developed (Hauf and Höller, 1987; Emanuel, 1994; Marquet, 2011; Marquet and Geleyn, 2015) to account for phase changes and deviations from thermodynamic equilibrium, e.g. by irreversible processes. By assuming only reversible processes (i.e. conserved entropy), approximate formulas can be derived (e.g. Emanuel, 1994). However, in the case of large hydrometeors, liquid or solid particles are removed due to gravitational acceleration, leading to an irreversible process; hence the formulas based on the assumption of a reversible process are no longer applicable. Sometimes for this situation a so-called pseudo-adiabatic potential temperature is defined, assuming instantaneous removal of hydrometeors from the air parcel; usually, meaningful approximations to this quantity are given, since generally it cannot be derived from first principles. Equivalent potential temperature including phase changes for vapour and liquid water is often used for the determination of convective instabilities. The general formulation can be easily adapted for an ice equivalent potential temperature, i.e. for reversible processes in pure ice clouds (see e.g. Spichtinger, 2014). Although the latent heat of sublimation is larger than the latent heat of vaporisation, the absolute mass content of water vapour decreases exponentially with decreasing temperature, leading to only small corrections due to phase changes in pure ice clouds.

At altitudes above the clouds' top, within the upper troposphere and across the tropopause, the air is substantially dried out compared to tropospheric in-cloud conditions. Therefore, above clouds and further aloft, e.g. within the stratosphere, the conventional dry-air potential temperature may suffice to provide a meaningful vertical coordinate. Moreover, the potential temperature or the virtual potential temperature, which includes water vapour, are commonly used as prognostic variables in numerical models for the formulations of the energy equation (e.g. Skamarock et al., 2005; Skamarock and Klemp, 2008; Zängl et al., 2015; Borchert et al., 2019). Thereby, very often both variants, the potential temperature as well as the equivalent potential temperature, are involved to account for dry-air situations and cloud conditions.

In any case, the use of the potential temperature requires the following preconditions to be fulfilled:

1. θ should be based on a rigorous derivation to ensure its validity as a function of atmospheric altitude in order not to corrupt its character as a vertical coordinate that allows for appropriately comparing (measured) atmospheric parameters, and

2. θ should approximate to the greatest possible extent the true entropy state of a probed air mass and should preferably account for the implied dependencies on atmospheric variables, even under the assumption that air behaves as an ideal gas,

with the aim that the potential temperature behaves as a rational physical variable. Thus, still abiding by the ideal-gas assumption, a re-assessment of the fundamental atmospheric quantity θ is suggested, which is based on the state of knowledge of air's thermodynamic properties, and this re-assessed θ is comprehensively examined concerning its ability to hold also for atmospheric conditions above the troposphere.

In principle, the concept of the potential temperature is transferable to all systems of thermally stratified fluids such as a planetary gas atmosphere or an ocean, to investigate heat fluxes (advection or diffusion) or the static stability of the fluid. In astrophysics, the potential temperature is used almost identically as in atmospheric sciences to describe dynamic processes and thermodynamic properties (e.g. static stability or vorticity) in the atmosphere of planets other than the Earth. Here, the same value p₀=1000 hPa, as applied to the Earth's atmosphere, is frequently used as a reference pressure for the atmosphere of other planets (Catling, 2015, Table 4), whereby the formulations of the specific heat capacity require adaptations to account for the individual gas composition of the respective planetary atmosphere. In order to simulate the weather in the atmosphere of other planets, the Weather Research and Forecasting Model (WRF) was extended to planetWRF (Richardson et al., 2007), and the governing equations considered within the WRF model include a prognostic equation for the potential temperature (Skamarock et al., 2005; Skamarock and Klemp, 2008). However, the temperature dependency of the isobaric heat capacity c_p is not generally negligible, especially when taking “deep atmospheres, such as on Venus” (Catling, 2015, p. 436) into account or the temperature lapse rates on other planets (Li et al., 2018). The atmosphere of Saturn's moon Titan, the only known moon with a substantial atmosphere, was comprehensively studied with frequent application of the potential temperature based on profile measurement of temperature and pressure in Titan's atmosphere by the Huygens probe (Müller-Wodarg et al., 2014).

Moreover, the potential temperature is a frequently used quantity in oceanography (e.g. McDougall et al., 2003; Feistel, 2008), while here the consideration of seawater's salinity and its impact on the specific heat capacity of seawater implies additional complexity. In particular, McDougall et al. (2003) suggests a re-assessment of the potential temperature as applied in oceanography to approximate the adiabatic lapse rate; thus this study bears certain parallels to the present investigation aiming at the reappraisal of the potential temperature for atmosphere-related purposes. These studies from other disciplines motivate the need for a re-assessment of the potential temperature for the atmospheric sciences. Thus, the approach provided herein proposes a modified calculation of the widely used quantity of the potential temperature by additionally accounting for the current state of knowledge concerning air's properties.

The study is organised as follows. The derivation of the potential temperature for an ideal gas with constant specific heat capacity c_p is recalled in Sect. 2. In Sect. 3 the assumption of a constant c_p is discussed together with a synopsis of various c_p values as provided in the literature. The temperature dependency of c_p is examined in Sect. 4, and a parameterisation is given. Section 5 is devoted to the definition and computation of a new reference potential temperature θ_ref based on the temperature-dependent specific heat capacity, while Sect. 6 focuses on the influence of real-gas effects on the resulting potential temperature. Section 7 presents some implications of the use of θ_ref, and concluding remarks are given in Sect. 8.

2 Derivation of the potential temperature for an ideal gas

The Gibbs equation (see e.g. Kondepudi and Prigogine, 1998) is a general thermodynamic relation to describe the state of a system with m components and reads as

$$\begin{array}{}\text{(3)}& T\mathrm{d}S=\mathrm{d}H-V\mathrm{d}p-\sum _{k=\mathrm{1}}^m{\mathit{\mu }}_k\mathrm{d}{M}_k,\end{array}$$

where T denotes the absolute temperature in kelvin (K), S the entropy in joules per kelvin (J K⁻¹), H the enthalpy in joules (J), V the volume in cubic metres (m³), μ_k the chemical potential of component k in joules per kilogram (J kg⁻¹), M_k the mass of component k in kilograms (kg), and p the static pressure in pascals (Pa). Assuming no phase conversion or chemical reaction within the system, the mass of each component does not change; hence dM_k=0 for each component k.

In the following, dry air is assumed to be the single component in the system. Expressing the Gibbs equation in its specific form (i.e. division by the total mass M_a of dry air; note that lowercase letters indicate specific variables – e.g. $h=H/M_{\mathrm{a}}$) leads to

$$\begin{array}{}\text{(4)}& T\mathrm{d}s=\mathrm{d}h-{\displaystyle \frac{V}{M_{\mathrm{a}}}}\mathrm{d}p\iff \mathrm{d}s={\displaystyle \frac{\mathrm{1}}{T}}\mathrm{d}h-{\displaystyle \frac{V}{M_{\mathrm{a}}T}}\mathrm{d}p.\end{array}$$

Furthermore, approximating dry air as an ideal gas leads to the following simplifications.

• The ideal-gas law

  $$\begin{array}{}\text{(5)}& pV=M_{\mathrm{a}}{R}_{\mathrm{a}}T\end{array}$$

  can be applied with the specific gas constant R_a of dry air, which is

  $$\begin{array}{}\text{(6)}& \begin{array}{rl}{R}_{\mathrm{a}}& ={\displaystyle \frac{R}{M_{\text{mol,a}}}}\\ & ={\displaystyle \frac{\mathrm{8.31446261815324}\mathrm{J}{\mathrm{mol}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}}{\mathrm{0.0289586}\mathrm{kg}{\mathrm{mol}}^{-\mathrm{1}}\pm \mathrm{0.0000002}\mathrm{kg}{\mathrm{mol}}^{-\mathrm{1}}}}\\ & \in \left[\mathrm{287.11350}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}},\mathrm{287.11748}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}\right],\end{array}\end{array}$$

  with the molar gas constant R in $\mathrm{J}{\mathrm{mol}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ (Tiesinga et al., 2020; Newell et al., 2018) and M_mol,a the molar mass of dry air (Lemmon et al., 2000), composed of nitrogen N₂, oxygen O₂, and argon Ar.

• The specific enthalpy is given by

  $$\begin{array}{}\text{(7)}& \mathrm{d}h=c_p\mathrm{d}T,\end{array}$$

  with c_p the specific heat capacity of dry air.

Based on these assumptions, the change in the specific entropy (within the fluid dry air) is given by

$$\begin{array}{}\text{(8)}& \mathrm{d}s={\displaystyle \frac{c_p}{T}}\mathrm{d}T-R_{\mathrm{a}}{\displaystyle \frac{\mathrm{d}p}{p}}.\end{array}$$

For isentropic changes of state, i.e. ds=0, Eq. (8) reduces to

$$\begin{array}{}\text{(9)}& {\displaystyle \frac{c_p}{T}}\mathrm{d}T=R_{\mathrm{a}}{\displaystyle \frac{\mathrm{d}p}{p}}.\end{array}$$

Note that the assumption of dry air being an ideal gas does not imply that in Eq. (9) the specific heat capacity c_p is constant. While statistical mechanics excludes any pressure dependence in the ideal-gas heat capacity, the general derivation (cf. Appendix A) permits a temperature dependence of c_p. However, usually the temperature dependence is neglected in atmospheric physics, and, instead, c_p is assumed as constant (see e.g. Ambaum, 2010, p. 48/49, where vibrational modes of the air molecules are neglected). Immediately below and in Sect. 3, the treatment of c_p as a temperature-independent constant is discussed. The introduction of the temperature dependence then follows in Sect. 4.

Treating c_p as a constant, rearrangement of Eq. (9) leads to

$$\begin{array}{}\text{(10)}& {\displaystyle \frac{\mathrm{d}T}{T}}={\displaystyle \frac{R_{\mathrm{a}}}{c_p}}{\displaystyle \frac{\mathrm{d}p}{p}}.\end{array}$$

Integration of Eq. (10) over the range from ground-level pressure and temperature (p₀, T₀) to the pressure and temperature at a specific height (p, T) yields

$$\begin{array}{}\text{(11)}& \ln \left({\displaystyle \frac{T}{T_{\mathrm{0}}}}\right)=\underset{T_{\mathrm{0}}}{\overset{T}{\int }}{\displaystyle \frac{\mathrm{d}{T}^{\prime }}{T^{\prime }}}={\displaystyle \frac{R_{\mathrm{a}}}{c_p}}\underset{p_{\mathrm{0}}}{\overset{p}{\int }}{\displaystyle \frac{\mathrm{d}{p}^{\prime }}{p^{\prime }}}={\displaystyle \frac{R_{\mathrm{a}}}{c_p}}\ln \left({\displaystyle \frac{p}{p_{\mathrm{0}}}}\right),\end{array}$$

and, after another straightforward conversion, one arrives at

$$\begin{array}{}\text{(12)}& \ln \left({\displaystyle \frac{T_{\mathrm{0}}}{T}}\right)={\displaystyle \frac{R_{\mathrm{a}}}{c_p}}\ln \left({\displaystyle \frac{p_{\mathrm{0}}}{p}}\right).\end{array}$$

With the definition ${\mathit{\theta }}_{c_p}=T_{\mathrm{0}}$, Eq. (12) is transformed into the commonly used expression for determining the potential temperature

$$\begin{array}{}\text{(13)}& {\mathit{\theta }}_{c_p}=T{\left({\displaystyle \frac{p_{\mathrm{0}}}{p}}\right)}^{\frac{R_{\mathrm{a}}}{c_p}},\end{array}$$

for which the ground-level pressure p₀ is arbitrary but usually set to p₀=1000 hPa. This choice coincides with the definition of the World Meteorological Organisation (WMO, 1966) and the standard-state pressure (Tiesinga et al., 2020) but should not be confused with the standard atmosphere 101 325 Pa (Tiesinga et al., 2020). In the following, ${\mathit{\theta }}_{c_p}$ denotes the potential temperature based on a constant c_p, and, when a specific value of c_p is applied, the subscript c_p in the potential temperature's notation is replaced by the corresponding c_p value.

3 Examining the assumption of constant c_p for dry air

The general theory of thermodynamics, assuming dry air as an ideal gas, gives the expression

$$\begin{array}{}\text{(14)}& c_p=\left(\mathrm{1}+{\displaystyle \frac{f}{\mathrm{2}}}\right)R_{\mathrm{a}}\end{array}$$

for the constant specific heat capacity, which is based on the results of statistical mechanics and the equipartition theorem (e.g. Huang, 1987). In Eq. (14), the parameter $f=f_{\mathrm{trans}}+f_{\mathrm{rot}}+f_{\mathrm{vib}}$ is equal to the total number of degrees of freedom of the gas molecules of which dry air consists. The individual contributions to f comprise the degrees of freedom of translation f_trans, rotation f_rot, and vibration f_vib. Assuming further that dry air exclusively consists of the linear molecules N₂ and O₂ (implying f_trans=3 and f_rot=2, while the contribution of Ar remains disregarded) and additionally neglecting the vibrational degrees of freedom (f_vib=0), the general relation Eq. (14) reduces to

$$\begin{array}{}\text{(15)}& c_p=\left(\mathrm{1}+{\displaystyle \frac{\mathrm{3}+\mathrm{2}}{\mathrm{2}}}\right)R_{\mathrm{a}}={\displaystyle \frac{\mathrm{7}}{\mathrm{2}}}{R}_{\mathrm{a}}.\end{array}$$

Although the neglect of vibrational excitation, particularly at very low temperatures, seems plausible and appropriate, errors are already introduced by this assumption for the temperature range relevant in the atmosphere.

In atmospheric sciences, for the majority of computations that require the specific heat capacity of dry air, a constant value of c_p may be appropriate. According to the WMO (1966), the recommended value for c_p of dry air is 1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$, and, furthermore (ibid.), it is defined that $\mathit{\gamma }=\frac{c_p}{c_v}=\frac{\mathrm{7}}{\mathrm{5}}=\mathrm{1.4}$ (cf. Eq. 1). This definition is consistent with the general thermodynamic theory together with all aforementioned additional assumptions and results in Eq. (15) as well.

Even assuming a universally valid constant c_p, a single consistently used value of c_p was not found. Instead, the specified values of c_p vary among different textbooks and other sources. In Table 1, some of the available values of constant specific heat capacity for dry air are compiled, indicating a variability of c_p that ranges from 994 to 1011 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$. However, the extremes in Table 1 are from old references of historical interest only; to reflect recently stated values the narrower range 1000 to 1010 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ is considered.

Wegener and Wegener (1935, converted from units other than SI)Roedel and Wagner (2011, p. 66)(WMO, 1966)Tripoli and Cotton (1981, the Appendix therein)Holton (2004, p. 491)Wallace and Hobbs (2006, p. 75)Schumann (2012)Wendisch and Brenguier (2013, p. 24)Liou (2002, Appendix F)Ambaum (2010, table “Useful Data”)Pruppacher and Klett (2010, converted from units other than SI; p. 489)Curry and Webster (1998, p. 62)WMO (1966)Bohren et al. (1998, p. 384)Houghton (2002, p. 275)Zdunkowski and Bott (2003, p. 705)Brasseur and Solomon (2005, p. 426)Seinfeld and Pandis (2006, p. 1178)Cotton et al. (2011, Table 2.1)Bolton (1980)Emanuel (1994, Appendix 2)Wendisch and Brenguier (2013, p. 69)Stamnes et al. (2017, p. 14)Chang et al. (2006)Tiwary and Williams (2019, beneath Eq. 8.8; possibly a typo, as indicated by inconsistencies on reproducing their conclusions based on this value)Brusseau et al. (2019, p. 59)(WMO, 1966)

Table 1 Synopsis of temperature-independent constant values given mainly in textbooks for the specific heat capacity c_p of dry air from various sources (non-exhaustive). Note that the WMO (1966) indicates a minimum and maximum “range of actual values” together with their recommended value c_p=1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$.

Figure 1 Vertical profiles of (a) atmospheric pressure and (b) temperature as functions of height, corresponding to the US Standard Atmosphere.

Figure 2 Computed vertical course of the potential temperature ${\mathit{\theta }}_{c_p}$ based on the two extremes of constant values for the specific heat capacity c_p provided in the literature including the historical extreme values (a; cf. Table 1), and (b) the absolute differences $\mathrm{\Delta }{\mathit{\theta }}_{c_p}={\mathit{\theta }}_{\mathrm{994}}$–θ₁₀₁₁ and $\mathrm{\Delta }{\mathit{\theta }}_{c_p}={\mathit{\theta }}_{\mathrm{1000}}$–θ₁₀₁₀ between the two resulting curves of ${\mathit{\theta }}_{c_p}$. The absolute difference $\mathrm{\Delta }{\mathit{\theta }}_{c_p}={\mathit{\theta }}_{\mathrm{1003.5}}$–θ_1006.5 is also shown (green curve), corresponding to a more realistic interval of c_p values.

These different values of constant c_p scatter within a small range (below ±1.1 %) around the WMO's recommendation 1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$, which may seem negligible if c_p contributes only as a linear coefficient within an equation (e.g. in the expression of a correction factor; cf. Weigel et al., 2016). However, in the formulation of the potential temperature ${\mathit{\theta }}_{c_p}$, cf. Eq. (13), the specific heat capacity c_p does not contribute linearly but rather as the denominator in the exponent. Thus, the variety of different c_p values, although scattering within a small range, impacts the resulting ${\mathit{\theta }}_{c_p}$ significantly. To illustrate this impact, a computation of ${\mathit{\theta }}_{c_p}$ by using Eq. (13) was based on the values of static pressure (p; cf. Fig. 1a) and absolute temperature (T; cf. Fig. 1b) corresponding to the US Standard Atmosphere (United States Committee on Extension to the Standard Atmosphere, 1976). From the list of the different c_p in Table 1, the extreme values were selected in order to initially illustrate the sensitivity of the resulting ${\mathit{\theta }}_{c_p}$ to variations in c_p in the range of ∼1 %, as seen in the literature. In Fig. 2a, the individual profiles of ${\mathit{\theta }}_{c_p}$ are shown for the extremes of the historic c_p values (Table 1), while Fig. 2b illustrates the absolute differences $\mathrm{\Delta }{\mathit{\theta }}_{c_p}={\mathit{\theta }}_{\mathrm{994}}$–θ₁₀₁₁ (red curve), $\mathrm{\Delta }{\mathit{\theta }}_{c_p}={\mathit{\theta }}_{\mathrm{1000}}$–θ₁₀₁₀ (blue curve), and $\mathrm{\Delta }{\mathit{\theta }}_{c_p}={\mathit{\theta }}_{\mathrm{1003.5}}$–θ_1006.5 (green curve). The absolute error exhibited with the blue curve in Fig. 2b is based on the extremes of most recently referred c_p values in the literature (Table 1). At an altitude of 8.5 km, the difference $\mathrm{\Delta }{\mathit{\theta }}_{c_p}$ already exceeds 1 K (blue curve). The values of $\mathrm{\Delta }{\mathit{\theta }}_{c_p}$ reach approximately 1.2 K at 10 km and rise further, above 4 K, with increasing altitude up to 20 km. At 50 km, approximately where the stratopause is located, which is the chosen upper height limit for this investigation, the computed $\mathrm{\Delta }{\mathit{\theta }}_{c_p}$ reaches 43 K. The green curve corresponds to the more realistic c_p interval $\mathrm{1005}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}\pm \mathrm{1.5}$ $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ as recommended by the WMO; the difference reaches approximately 13 K at the stratopause.

Figure 2 illustrates the possible spread of ${\mathit{\theta }}_{c_p}$ based on a range of c_p values from different literature references; hence, if one uses a different value for c_p from the literature than that defined by WMO (1966), the difference ${\mathit{\theta }}_{\mathrm{1005}}-{\mathit{\theta }}_{c_p}$ might be significant. Since the c_p values provided by some literature references are close to the value c_p=1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ recommended by the WMO (1966), the subsequent comparisons will be made to θ₁₀₀₅. The ${\mathit{\theta }}_{c_p}$ values based on c_p values other than 1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ are only used to illustrate respective deviations. Although the curves in Fig. 2b depict extremes in the deviation of potential temperatures, as they are based on the extremes of c_p values (cf. Table 1), they nevertheless illustrate the sensitive response of ${\mathit{\theta }}_{c_p}$ to even small variations in c_p, on the order of 1 %. Further proof of this sensitivity from the mathematical perspective is provided in Appendix B. The impact of this sensitivity becomes important at altitudes of ∼10 km and above, thus where the use of the potential temperature becomes increasingly meaningful. Here, and in particular above the cloud tops, the small-scale and comparatively fast tropospheric dynamics (causing vertical transport and implying diabatic processes) become diminished, while further above, towards the stratosphere, an increasingly layered vertical structure of the atmosphere is taking over.

As indicated above, the reason for this sensitivity to small variations of air's specific heat capacity is that it affects the exponent of the equation for ${\mathit{\theta }}_{c_p}$. The studies of Ooyama (1990, 2001) document an interesting attempt to formulate for example the energy balance equations for the moist atmosphere, wherein entropy replaces the more common formulation using the potential temperature. This substitution avoids the use of the potential temperature, which “is merely an exponential transform of the entropy expressed in units of temperature” (Ooyama, 2001); thus, within this equation, air’s specific heat capacity is implied exclusively as a linear coefficient. Consequently, a parameterisation for the temperature dependence of the specific heat capacity (c_p(T); cf. Sect. 4) may be easily adopted. However, the crucial drawback of the entropy-based equations is that to gain a numerical model for weather forecast purposes, the parameterisations of most of the physical processes within the atmosphere would require a reformulation.

It should be noted that not only do literature values of air's specific heat capacity c_p vary, but also the values of the gas constant R_a vary slightly due to different historical approximations for the molar gas constant⁴ R and for the composition of dry air. The variation of values for R_a is typically only on the order of 0.1 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$, whereas the variability in c_p is on the order of a few $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ (cf. Table 1). Therefore, within the exponent of the expression (13) for ${\mathit{\theta }}_{c_p}$, the variability of c_p has by far a stronger impact on the resulting ${\mathit{\theta }}_{c_p}$ value than the variability of R_a.

However, accepting for a moment the WMO's definition (15) of c_p (WMO, 1966), the variability of air's c_p should naturally be constrained to certain limits. With the specific gas constant R_a=287.05 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ (WMO, 1966), the WMO's definition leads to c_p=1004.675 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$. In contrast, taking into account the uncertainty introduced in R_a by the molar mass of dry air, cf. Eq. (6), the resulting range for air's specific heat capacity is $\mathrm{1004.897}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}\le c_p\le \mathrm{1004.912}$ $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$. It may be surmised that the rounded value c_p=1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ as recommended by the WMO (1966) had the main goal to simplify certain calculations, which at the time may have been mostly done by hand.

4 Accounting for the temperature dependence of air's specific heat capacity

Next, while retaining the ideal-gas assumption, we consider the dependence of air's c_p on temperature, mainly over the atmospherically relevant range (180 to 300 K). The temperature dependence of c_p is, of course, not a new finding. Experimental approaches for determining the calorimetric properties of air and the temperature dependence of a fluid's specific heat capacity are described by Witkowski (1896), who investigated the change in the mean c_p as a function of temperature intervals between room temperature (as a fixed reference) and various warmer and colder temperatures, for atmospheric pressures and slightly beyond. Despite the potentially high uncertainty of the experimental results from these times, Witkowski (1896) already indicated that with decreasing temperature the experimentally determined c_p values initially decline, then pass a minimum, and subsequently increase again at lower temperatures (T<170 K). The description of refined experiments and ascertainable data of air's c_p(T) for temperatures below 293 K is summarised by Scheel and Heuse (1912), Jakob (1923), and Roebuck (1925, 1930), illustrating in comprehensive detail the experimental effort and providing the resulting data. The review by Awano (1936) compiled and compared the data of c_p(T) of dry air (“air containing neither carbon-dioxide nor steam”, Awano, 1936), and he attested – at that time – the previously mentioned studies to constitute “the most reliable experiments”. During the decades following these experiments, further insights were gained and landmarks were reached, which are summarised in the comprehensive survey by Lemmon et al. (2000) of the progress of modern formulations for the thermodynamic properties of air and about the experiments the previous formulations were based on.

Figure 3 Variety of suggested values for the specific heat capacity of air. Ranges of constant values for c_p (including the historical) together with the recommended value by the WMO (1966) are displayed as given in Table 1 (dashed lines). The parameterisations of air's $c_p^{\mathrm{0}}\left(T\right)$, assuming dry air as an ideal gas, accounting for its temperature dependence by Lemmon et al. (2000, solid magenta curve) and by Dixon (2007, solid cyan curve) are displayed. Discrete measurement and literature data at about 1000 hPa (i.e. as often specified, at “about one atmosphere”) are indicated by dots. In addition, the studies by Awano (1936) and Vasserman et al. (1966) provide data at other atmospheric pressures, as indicated by squares, diamonds, and triangles.

Figure 3 illustrates the range of suggested constant values for the specific heat capacity as indicated in Table 1 (dashed curves) together with the measurements that were made to obtain air's behaviour as a function of temperature and pressure. Note that Fig. 3 includes data at other atmospheric pressures, indicated by squares, diamonds, and triangles. In the same figure, calculated values of c_p(T) of dry air are displayed, resulting from the equation of state which was derived from experimental p, V, and T data by Vasserman et al. (1966), who provided an extensive review of previous experimental and theoretical works and of the state of knowledge at that time. In addition, Fig. 3 exhibits two different parameterisations, by Lemmon et al. (2000) and by Dixon (2007, see p. 376 in his book – the accuracy is “within 0.1 % from 200 to 450 K”), which account for the temperature dependence of the specific heat capacity c_p(T). The parameterisation by Lemmon et al. (2000), to be discussed in detail in Sect. 4.2, is valid for dry air assumed as an ideal gas, whereas this distinction is not explicitly made in Dixon (2007). Moreover, Fig. 3 contains discrete values of dry air's c_p(T) extracted from the database REFPROP (Reference Fluid Thermodynamic and Transport Properties Database by NIST, the National Institute of Standards and Technology, Lemmon et al., 2018), which is based on parameterisations resulting from thermodynamic considerations discussed later.

The measurement data, as well as the parameterisations, clearly indicate a dependence of air's specific heat capacity on the temperature. At temperatures above 300 K, the data points by Jakob (1923) are surprisingly well captured by the parameterisations, while below 270 K the course of the parameterised and measured c_p(T) diverges significantly. Possible reasons for this include the following:

• the measurements of c_p(T) have a precision likely no better than 1 % (in particular the historical measurements), and there could be systematic errors, especially at low temperatures;

• the measured data reflect the true thermodynamic behaviour of the real gas rather than that of an ideal gas.

However, it is immediately obvious from Fig. 3 that a good agreement among (i) the experimentally determined c_p(T) data, (ii) a constant c_p (e.g. 1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$; WMO, 1966), and (iii) the parameterised c_p(T) is found only for a temperature interval ranging from 270 to 300 K. For air temperatures below 270 K, the constant value c_p=1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ is only comparable with the values from Vasserman et al. (1966) but fails to coincide with other parameterised or experimentally determined values of c_p(T).

4.1 The temperature dependence of the ideal-gas specific heat capacity

As already indicated by the data depicted in Fig. 3, the specific heat capacity c_p depends on the gas temperature. With regard to measured values, the lack of constancy may be due to real-gas effects or to a dependence of the ideal-gas heat capacity on temperature. In this section, we focus on the latter effect, denoting the ideal-gas isobaric specific heat capacity by $c_p^{\mathrm{0}}\left(T\right)$, where the superscript 0 indicates the underlying ideal-gas assumption. For an individual gas, there is always a contribution from the three translational degrees of freedom, $c_{p,\mathrm{trans}}^{\mathrm{0}}=\frac{\mathrm{5}}{\mathrm{2}}{R}_i$, where R_i is the specific gas constant of the gas. If the molecule is assumed to be a rigid rotor, there is also a rotational contribution given by

$$\begin{array}{}\text{(16)}& c_{p,\mathrm{rot}}^{\mathrm{0}}=\left\{\begin{array}{ll}{R}_i,& \text{for linear (e.g. diatomic) molecules},\\ \frac{\mathrm{3}}{\mathrm{2}}{R}_i,& \text{for nonlinear molecules}.\end{array}\right.\end{array}$$

As mentioned previously, at finite temperatures molecules also have contributions to $c_p^{\mathrm{0}}\left(T\right)$ from intramolecular vibrations (and, at high temperatures, excited electronic states). To arrive at a temperature-dependent parameterisation for the ideal-gas specific heat capacity of dry air, the compounds' individual contributions, considering all degrees of freedom, need to be parameterised and then combined according to each compound's proportion in the mixture. For the following, dry air is considered a three-component mixture: the diatomic gases nitrogen (N₂) and oxygen (O₂) and the monatomic gas argon (Ar).

To determine the contribution of N₂ to $c_p^{\mathrm{0}}\left(T\right)$, both Bücker et al. (2002) and Lemmon et al. (2000) use the ideal-gas heat capacity from the reference equation of state of Span et al. (2000) that compares well with the findings from other studies within an uncertainty $\mathrm{\Delta }{c}_p^{\mathrm{0}}$ of less than 0.02 %.

For the contribution of O₂, Lemmon et al. (2000) use the formulation given by Schmidt and Wagner (1985). Alternatively, Bücker et al. (2002) provide a slightly different formulation from the International Union of Pure and Applied Chemistry (IUPAC, Wagner and de Reuck, 1987), after refitting it to more recently obtained data, thereby achieving an overall uncertainty $\mathrm{\Delta }{c}_p^{\mathrm{0}}$ of less than ±0.015 % for O₂ (Bücker et al., 2002). However, the difference in the resulting specific heat capacity contribution by O₂ between the two approaches (Lemmon et al., 2000, or Bücker et al., 2002) is comparatively small. The recent work of Furtenbacher et al. (2019) leads to values of $c_p^{\mathrm{0}}$ for O₂ with even smaller uncertainties, but the differences from the values used here are negligible in our context.

For a monatomic gas such as Ar, vibrational and rotational contributions to the heat capacity do not exist, and Bücker et al. (2002) consider that argon's excited electronic states are relevant only at temperatures above 3500 K. Hence, the contribution of Ar to the specific heat capacity of air reduces to $c_p^{\mathrm{0}}=\frac{\mathrm{5}}{\mathrm{2}}{R}_{\mathrm{Ar}}$.

The approach by Bücker et al. (2002) additionally considers the contribution of further constituents of air, such as water, carbon monoxide, carbon dioxide, and sulfur dioxide. These authors provide an analytical expression for specific heat capacity, accounting for this more complex but proportionally invariant air composition which is specified to deviate from the used reference by $\mathrm{\Delta }{c}_p^{\mathrm{0}}\le \pm \mathrm{0.015}\mathrm{\%}$ in the temperature range of $\mathrm{200}\mathrm{K}\le T\le \mathrm{3300}\mathrm{K}$. At atmospheric altitudes above the clouds' top, i.e. on average above ∼11 km, the air is assumed to have lost most of its water and is deemed as dry. Furthermore, for the following, trace gases that contribute to air's composition by molar fractions of less than that of Ar are neglected.

4.2 NIST's parameterisation of $c_p^{\mathrm{0}}\left(T\right)$

Besides a comprehensive survey of the available experimental data for the specific heat capacity of air, Lemmon et al. (2000) also provide state-of-the-art knowledge for other thermodynamic properties (isochoric heat capacity, speed of sound, vapour–liquid equilibrium, etc.). Additionally, they give two approaches to derive air's thermodynamic properties, including the vapour–liquid equilibrium:

1. an empirical model-based equation of state for standard (dry) air considered as a pseudo-pure fluid and

2. assembly of a mixture model from equations of state for each pure fluid.

Each approach allows calculating the thermodynamic properties, e.g. c_p, of gas mixtures such as dry air, and both are real-gas models with the ideal-gas behaviour as a boundary condition. The major difference between the models is that the first approach considers air as a pseudo-pure fluid while the second, more rigorous approach treats air as a mixture composed of N₂, O₂, and Ar, in molar fractions of 0.7812, 0.2096, and 0.0092, respectively, following Lemmon et al. (2000, their Table 3). This fractional composition of dry air is assumed to be constant from ground level up to 80 km height (United States Committee on Extension to the Standard Atmosphere, 1976), and its fractional composition would have to be shifted significantly to cause a serious deviation of the resulting potential temperature. The contribution to the composition by carbon dioxide (CO₂) and of any other trace species was assumed to be negligible. The validity of both approaches is specified for various states of dry air, from its solidification point (59.75 K) up to temperatures of 1000 K, and for pressures up to 100 MPa and even much further beyond the pressure range that is relevant for atmospheric investigations. Both the pseudo-pure fluid model and the mixture model are implemented in NIST's REFPROP database (cf. https://www.nist.gov/srd/refprop, last access: 2 November 2020) for various physical properties of fluids over a wide range of temperatures and pressures.

Both the pseudo-pure fluid model and the mixture model of Lemmon et al. (2000) use the same expression for the ideal-gas heat capacity, which is rigorously given as a sum of the pure-component contributions:

$$\begin{array}{}\text{(17)}& \begin{array}{rl}{\displaystyle \frac{C_p^{\mathrm{0}}\left(T\right)}{R}}& =x_{{\mathrm{N}}_{\mathrm{2}}}{\left({\displaystyle \frac{C_p^{\mathrm{0}}\left(T\right)}{R}}\right)}_{{\mathrm{N}}_{\mathrm{2}}}+x_{\mathrm{Ar}}{\left({\displaystyle \frac{C_p^{\mathrm{0}}\left(T\right)}{R}}\right)}_{\mathrm{Ar}}\\ & +x_{{\mathrm{O}}_{\mathrm{2}}}{\left({\displaystyle \frac{C_p^{\mathrm{0}}\left(T\right)}{R}}\right)}_{{\mathrm{O}}_{\mathrm{2}}},\end{array}\end{array}$$

where x_i denotes the molar fraction of species i, and $C_p^{\mathrm{0}}$ and the molar gas constant R are given in units of $\mathrm{J}{\mathrm{mol}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$.

Like Bücker et al. (2002), Lemmon et al. (2000) use the expression of Span et al. (2000) for the contribution of N₂ to the heat capacity and adopt $C_p^{\mathrm{0}}=\frac{\mathrm{5}}{\mathrm{2}}R$ for Ar. Together with the contribution by O₂ according to the formulation by Schmidt and Wagner (1985), the expression provided by Lemmon et al. (2000, Eq. 18 therein) for the ideal-gas heat capacity of dry air is

$$\begin{array}{}\text{(18)}& \begin{array}{rl}{\displaystyle \frac{C_p^{\mathrm{0}}\left(T\right)}{R}}& =N_{\mathrm{1}}+N_{\mathrm{2}}T+N_{\mathrm{3}}{T}^{\mathrm{2}}+N_{\mathrm{4}}{T}^{\mathrm{3}}+N_{\mathrm{5}}{T}^{-\frac{\mathrm{3}}{\mathrm{2}}}\\ & +N_{\mathrm{6}}{\displaystyle \frac{\frac{N_{\mathrm{9}}^{\mathrm{2}}}{T^{\mathrm{2}}}\exp \left(\frac{N_{\mathrm{9}}}{T}\right)}{{\left(\exp \left(\frac{N_{\mathrm{9}}}{T}\right)-\mathrm{1}\right)}^{\mathrm{2}}}}+N_{\mathrm{7}}{\displaystyle \frac{\frac{N_{\mathrm{10}}^{\mathrm{2}}}{T^{\mathrm{2}}}\exp \left(\frac{N_{\mathrm{10}}}{T}\right)}{{\left(\exp \left(\frac{N_{\mathrm{10}}}{T}\right)-\mathrm{1}\right)}^{\mathrm{2}}}}\\ & +{\displaystyle \frac{\mathrm{2}{N}_{\mathrm{8}}}{\mathrm{3}}}{\displaystyle \frac{\frac{N_{\mathrm{11}}^{\mathrm{2}}}{T^{\mathrm{2}}}\exp \left(-\frac{N_{\mathrm{11}}}{T}\right)}{{\left(\frac{\mathrm{2}}{\mathrm{3}}\exp \left(-\frac{N_{\mathrm{11}}}{T}\right)+\mathrm{1}\right)}^{\mathrm{2}}}},\end{array}\end{array}$$

with the scalar coefficients N_i for dry air (ibid.),

$$\begin{array}{}\text{(19)}& \begin{array}{ll}{N}_{\mathrm{1}}=\mathrm{3.490888032},& N_{\mathrm{2}}=\mathrm{2.395525583}\times {\mathrm{10}}^{-\mathrm{6}},\\ N_{\mathrm{3}}=\mathrm{7.172111248}\times {\mathrm{10}}^{-\mathrm{9}},& N_{\mathrm{4}}=-\mathrm{3.115413101}\times {\mathrm{10}}^{-\mathrm{13}},\\ N_{\mathrm{5}}=\mathrm{0.223806688},& N_{\mathrm{6}}=\mathrm{0.791309509},\\ N_{\mathrm{7}}=\mathrm{0.212236768},& N_{\mathrm{8}}=\mathrm{0.197938904},\\ N_{\mathrm{9}}=\mathrm{3364.011},& N_{\mathrm{10}}=\mathrm{2242.45},\\ N_{\mathrm{11}}=\mathrm{11580.4},& \end{array}\end{array}$$

which is specified as valid for temperatures from 60 to 2000 K. Because the underlying calculations are based on rigorous statistical mechanics and accurate spectroscopic data, $\frac{C_p^{\mathrm{0}}\left(T\right)}{R}$ should be accurate to within 0.01 % throughout this range, as discussed by Span et al. (2000).

The parameterisation (18) provides the isobaric specific heat capacity of dry air, considered as a mixture of ideal gases. This represents a more rigorous and accurate behaviour than assuming it to be a constant.

4.3 The parameterisation of $c_p^{\mathrm{0}}\left(T\right)$ from an engineer's perspective

The parameterisation from Dixon (2007)

$$\begin{array}{}\text{(20)}& c_p\left(T\right)=\mathrm{1002.5}+\mathrm{275}\times {\mathrm{10}}^{-\mathrm{6}}\cdot {\left(T-\mathrm{200}\right)}^{\mathrm{2}}\end{array}$$

for $\mathrm{200}\mathrm{K}\le T\le \mathrm{450}\mathrm{K}$ is not explicitly described to be based on particular assumptions or data sets. The author indicates his suggested parameterisation to hold within 0.1 % for temperatures between 200 and 450 K. For elevated air temperatures, the deviation between the ideal-gas limit $c_p^{\mathrm{0}}\left(T\right)$ (Lemmon et al., 2000) and Dixon's parameterisation substantially increases. This is most likely due to the chosen type of polynomial approximation (Dixon, 2007), which increasingly departs from the reference $c_p^{\mathrm{0}}\left(T\right)$ for gas temperatures exceeding 450 K.

Concerning the thermophysical properties of humid air, the study by Tsilingiris (2008) provides further insight. Its purpose was to evaluate the transport properties as a function of different levels of the relative humidity and as a function of temperature (from 273 to 373 K) for the gas mixture of air with water vapour at a constant pressure (1013 hPa). The atmospherically relevant pressure range below 1013 hPa and temperatures smaller than 273 K were not considered. Although this study focused on providing a comprehensive account of moisture within air, mainly for technical purposes and engineering calculations, the possible usefulness of these findings to atmospheric investigations is also apparent. However, the impact of water vapour on the resulting gas mixture's c_p(T) is significantly larger (cf. Tsilingiris, 2008) than the uncertainty of dry air's c_p(T) that is discussed in the present work. Furthermore, the consideration of water vapour as a component of air requires very individual and case-specific computations of c_p(T) of moist air, as water vapour is among the most variable constituents of the atmosphere.

The effort required to produce an analytical formulation for gas properties which best reflects the true gas behaviour may indicate that for engineering purposes (pneumatic shock absorbers, engines' combustion efficiency, improvements of turbofan/turboprop propulsion, aerodynamics, material sciences, etc.), especially where pressures exceed atmospheric, the assumption of ideal-gas behaviour introduces excessive uncertainty.

5 The ${\mathit{\theta }}_{c_p\left(T\right)}$ from the temperature-dependent specific heat capacity of air

Previously introduced approaches for computing the specific heat capacity of dry air call for a brief discussion on how to use the obtained c_p(T) to derive the potential temperature. In the following, ${\mathit{\theta }}_{c_p\left(T\right)}$ denotes the derived potential temperature that accounts for the temperature dependence of dry air's specific heat capacity. Furthermore, it should be noted that simply substituting any c_p(T) value into the conventionally used equation and defining Eq. (13) for ${\mathit{\theta }}_{c_p}$ (WMO, 1966) may appear tempting but definitely leads to results inconsistent with ${\mathit{\theta }}_{c_p\left(T\right)}$, which is based on the reference parameterisation of dry air's c_p(T). Therefore, the thermodynamically consistent use of c_p(T) in the derivation of θ is described in the following.

5.1 Derivation of ${\mathit{\theta }}_{c_p\left(T\right)}$ based on the temperature-dependent specific heat capacity of dry air

In the derivation of the potential temperature (cf. Sect. 2), we note that, until reaching the expression for isentropic changes of state (9), no assumption was made about the specific heat capacity. As soon as the temperature dependence of the specific heat capacity comes into play, the re-assessment of Eq. (9) leads to

$$\begin{array}{}\text{(21)}& {\displaystyle \frac{c_p\left(T\right)}{T}}\mathrm{d}T=R_{\mathrm{a}}{\displaystyle \frac{\mathrm{d}p}{p}}.\end{array}$$

Integration of Eq. (21) from the basic state $(p_{\mathrm{0}},{\mathit{\theta }}_{c_p\left(T\right)})$ to any other state (p, T) yields

$$\begin{array}{}\text{(22)}& R_{\mathrm{a}}\ln \left({\displaystyle \frac{p}{p_{\mathrm{0}}}}\right)=R_{\mathrm{a}}\underset{p_{\mathrm{0}}}{\overset{p}{\int }}{\displaystyle \frac{\mathrm{d}{p}^{\prime }}{p^{\prime }}}=\underset{{\mathit{\theta }}_{c_p\left(T\right)}}{\overset{T}{\int }}{\displaystyle \frac{c_p\left(T^{\prime }\right)}{T^{\prime }}}\mathrm{d}{T}^{\prime },\end{array}$$

where ${\mathit{\theta }}_{c_p\left(T\right)}$ is the desired potential temperature.

The rearrangement of Eq. (22) makes evident that the desired potential temperature is a zero of the function F(x), given by

$$\begin{array}{}\text{(23)}& F\left(x\right)=\underset{x}{\overset{T}{\int }}{\displaystyle \frac{c_p\left(T^{\prime }\right)}{T^{\prime }}}\mathrm{d}{T}^{\prime }-R_{\mathrm{a}}\ln \left({\displaystyle \frac{p}{p_{\mathrm{0}}}}\right).\end{array}$$

To arrive at the desired potential temperature ${\mathit{\theta }}_{c_p\left(T\right)}$ for any given temperature and pressure, the equation 0=F(x) must be solved for the variable x, which is the desired ${\mathit{\theta }}_{c_p\left(T\right)}$. Equation (23) has at most only one real zero, since its integrand is strictly positive, which means F(x) is strictly monotonic.

In the following, the ideal-gas reference potential temperature θ_ref is introduced, based on the formulation of the ideal-gas limit of dry air's specific heat capacity $c_p^{\mathrm{0}}\left(T\right)$ in accordance with Eq. (18) as formulated by Lemmon et al. (2000). This reference potential temperature θ_ref represents the zero of F(x) in Eq. (23), wherein c_p(T^′) is to be replaced by $c_p^{\mathrm{0}}\left(T^{\prime }\right)$; i.e. for given p, T the reference potential temperature θ_ref solves the equation

$$\begin{array}{}\text{(24)}& \mathrm{0}=F\left({\mathit{\theta }}_{\mathrm{ref}}\right)=\underset{{\mathit{\theta }}_{\mathrm{ref}}}{\overset{T}{\int }}{\displaystyle \frac{c_p^{\mathrm{0}}\left(T^{\prime }\right)}{T^{\prime }}}\mathrm{d}{T}^{\prime }-R_{\mathrm{a}}\ln \left({\displaystyle \frac{p}{p_{\mathrm{0}}}}\right).\end{array}$$

The parameterisation of $c_p^{\mathrm{0}}\left(T^{\prime }\right)$ is stated to give accurate values for temperatures from 60 to 2000 K (cf. Sect. 4.2); thus values of θ_ref should not exceed 2000 K, since otherwise $c_p^{\mathrm{0}}\left(T^{\prime }\right)$ within the integrand in Eq. (24) is evaluated outside of its range of validity. However, due to the division by T^′, the value of the integrand $\frac{c_p^{\mathrm{0}}\left(T^{\prime }\right)}{T^{\prime }}$ may be expected to give nevertheless a good approximation even if the accuracy of $c_p^{\mathrm{0}}\left(T^{\prime }\right)$ is decreased; hence values θ_ref>2000 K should not be discarded.

It may be noted that further variants of a reference potential temperature are derivable by replacing c_p(T^′) in Eq. (23) by any other expression of the specific heat capacity of air which may appear sufficiently accurate. The steps to compute or approximate the zero of the function (23), described in this study, are independent of the chosen heat capacity formulation.

Unfortunately, for a straightforward solution of the integral (23), the suggested parameterisation of c_p is too complex and an analytically insolvable nonlinear equation 0=F(x) could result. Thus, an approximation of the equation's desired zero is required. Newton's method (cf. e.g. Deuflhard, 2011) provides a standard approach to numerically approximate the zero of a nonlinear equation. Proceeding from an initial guess x₀, Newton's method constructs a sequence {x_k}_k∈ℕ defined by the recursion

$$\begin{array}{}\text{(25)}& \begin{array}{rl}{x}_{k+\mathrm{1}}& =x_k-{\displaystyle \frac{F\left(x_k\right)}{F^{\prime }\left(x_k\right)}}=x_k-{\displaystyle \frac{F\left(x_k\right)}{-\frac{c_p\left(x_k\right)}{x_k}}}\\ & ={\displaystyle \frac{x_k}{c_p\left(x_k\right)}}\left[c_p\left(x_k\right)+F\left(x_k\right)\right]\\ & ={\displaystyle \frac{x_k}{c_p\left(x_k\right)}}\left[c_p\left(x_k\right)-R_{\mathrm{a}}\ln \left({\displaystyle \frac{p}{p_{\mathrm{0}}}}\right)+\underset{x_k}{\overset{T}{\int }}{\displaystyle \frac{c_p\left(T^{\prime }\right)}{T^{\prime }}}\mathrm{d}{T}^{\prime }\right].\end{array}\end{array}$$

The constructed sequence {x_k}_k∈ℕ converges to the equation's desired zero. For the computations described here, the iteration is stopped as soon as the absolute difference $\left|x_{k+\mathrm{1}}-x_k\right|$ of two consecutive iterations falls below 10⁻⁸ K.

For the reference of air's specific heat capacity, $c_p^{\mathrm{0}}\left(T\right)$, the integral (23) turns out not to be explicitly solvable. Therefore, with each iteration, the solution of the integral $\underset{x_k}{\overset{T}{\int }}\frac{c_p^{\mathrm{0}}\left(T^{\prime }\right)}{T^{\prime }}\mathrm{d}{T}^{\prime }$ is approximated by subdividing the entire integration range, [x_k, T], into intermediate intervals with respective size of at most 0.1 K and by applying Simpson's rule on each subinterval.

As a first guess x₀ for the Newton iteration, the conventional definition of ${\mathit{\theta }}_{c_p}$ based on a constant specific heat capacity (WMO, 1966) is inserted:

$$\begin{array}{}\text{(26)}& x_{\mathrm{0}}=T{\left({\displaystyle \frac{p_{\mathrm{0}}}{p}}\right)}^{\frac{R_{\mathrm{a}}}{\mathrm{1005}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}}}={\mathit{\theta }}_{\mathrm{1005}}.\end{array}$$

In the course of Newton's method, the sequence {x_k}_k∈ℕ will converge to the unique zero for any initial guess x₀ due to the monotonicity of F(x). However, the right choice of the initial guess x₀ substantially decreases the error of the first iteration x₁, speeding up convergence to the desired zero of the function F(x). Therefore, it seems wise to use the conventional definition of ${\mathit{\theta }}_{c_p}$ as the first guess for the Newton iteration (25).

Figure 4 (a) Reference potential temperature θ_ref together with the potential temperatures ${\mathit{\theta }}_{\mathrm{994}},{\mathit{\theta }}_{\mathrm{1000}},{\mathit{\theta }}_{\mathrm{1010}}$ and θ₁₀₁₁ relying on constant c_p values (the dashed lines depict the historical extremes for c_p; cf. Table 1). (b) Relative differences $\left({\mathit{\theta }}_{c_p}-{\mathit{\theta }}_{\mathrm{ref}}\right)/{\mathit{\theta }}_{\mathrm{ref}}$ for the same choices $c_p\in \left\{\mathrm{994}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}},\mathrm{1000}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}},\mathrm{1010}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}},\mathrm{1011}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}\right\}$ as in panel (a) between the reference potential temperature and the potential temperatures relying on constant c_p values. For comparison, the relative difference $\left({\mathit{\theta }}_{\mathrm{1005}}-{\mathit{\theta }}_{\mathrm{ref}}\right)/{\mathit{\theta }}_{\mathrm{ref}}$ is displayed, for which c_p=1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ corresponds to the WMO recommendation. In addition, also comparisons with ${\mathit{\theta }}_{\mathrm{1003.5}},{\mathit{\theta }}_{\mathrm{1004}},{\mathit{\theta }}_{\mathrm{1006.5}}$ are included. All profiles are based on the values for temperature and pressure according to the US Standard Atmosphere. Note the linear axis scaling inside and the logarithmic scaling outside of the grey-shaded area in panel (b).

Solving the previously described root-finding problem by Newton's method over the comprehensive range of iteration steps (until the set requirement, i.e. $\left|x_{k+\mathrm{1}}-x_k\right|<{\mathrm{10}}^{-\mathrm{8}}$ K, is fulfilled) finally leads to the reference potential temperature θ_ref. This θ_ref is based on the ideal-gas limit of dry air's specific heat capacity $c_p^{\mathrm{0}}\left(T\right)$, which refers to the current thermodynamic state of knowledge, and, thus, we use θ_ref as our reference for the potential temperature in the following. For evaluating the results, the air temperature and pressure from the US Standard Atmosphere are used once more to set up the vertical profiles of the potential temperature. Figure 4a exhibits the resulting reference profile, i.e. θ_ref (red curve). Additionally, for comparison with the reference, further potential temperature profiles ${\mathit{\theta }}_{c_p}$ are shown based on the two (historical) extremes c_p=994 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ and c_p=1011 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ (dashed curves) and based on the range limits of more recent values c_p=1000 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ and c_p=1010 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ (solid green and magenta curves) of given constant values of air's specific heat capacity (cf. Table 1). Clearly, in particular at elevated altitudes, the courses of θ₁₀₀₀ and θ₁₀₁₀ significantly deviate from the reference. To quantitatively evaluate the match between the different profiles, the relative difference of the profiles based on a constant c_p, with respect to the reference, i.e. $\mathrm{\Delta }\mathit{\theta }/{\mathit{\theta }}_{\mathrm{ref}}=\left({\mathit{\theta }}_{c_p}-{\mathit{\theta }}_{\mathrm{ref}}\right)/{\mathit{\theta }}_{\mathrm{ref}}$, is depicted in Fig. 4b. The comparison demonstrates that the ${\mathit{\theta }}_{c_p}$ profiles significantly depart from the reference by about 300 K at 50 km altitude, corresponding to a relative difference of about 16 %. With both extremes of the recent constant values $c_p\in \left\{\mathrm{1000}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}},\mathrm{1010}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}\right\}$, the relative error level of 0.1 % is exceeded at altitudes about 5 km. While θ₁₀₀₀ continues to increasingly deviate from the reference, θ₁₀₁₀ re-enters and crosses the 0.1 % relative error interval (grey-shaded area) at altitudes between ∼19 and 21 km, before it reaches similar errors to the other ${\mathit{\theta }}_{c_p}$ profiles that are based on a constant c_p. Although the extreme values $c_p\in \left\{\mathrm{1000}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}},\mathrm{1010}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}\right\}$ appear in recent literature, these values may be considered unrealistic. For this reason, Fig. 4b also shows the relative deviations for the values $c_p\in \left\{\mathrm{1003.5}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}},\mathrm{1004}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}},\mathrm{1005}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}\right.$, $\left.\mathrm{1006.5}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}\right\}$, which include the recommended value of the WMO (1966) and a more realistic range; i.e. $c_p=\mathrm{1005}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}\pm \mathrm{1.5}$ $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$. Notably, up to an altitude of 15 km, the reference potential temperature is comparably well matched by both the recommended θ₁₀₀₅ and θ₁₀₀₄ (based on the frequently used alternative c_p=1004 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$; cf. Table 1). Until 15 km altitude, both constant c_p values lead to errors of calculated ${\mathit{\theta }}_{c_p}$ which remain comparatively small within the 0.1 % relative error interval. However, above ∼17.5 km, both θ₁₀₀₄ and θ₁₀₀₅ exceed the 0.1 % relative error interval, and further aloft their relative error with respect to the reference θ_ref increases rapidly.

In the context of numerical models of the atmosphere, the energy balance equation is occasionally formulated based on the potential temperature θ; thus θ constitutes a prognostic model variable. In such a case, the temperature T needs to be calculated from a given pair of values of pressure p and potential temperature θ. Using once more the defining Eq. (22), for given θ a zero of the function

$$\begin{array}{}\text{(27)}& \mathrm{0}=-\underset{T}{\overset{\mathit{\theta }}{\int }}{\displaystyle \frac{c_p\left(T^{\prime }\right)}{T^{\prime }}}\mathrm{d}{T}^{\prime }-R_{\mathrm{a}}\ln \left({\displaystyle \frac{p}{p_{\mathrm{0}}}}\right)\end{array}$$

is to be computed. Since Eq. (27) corresponds to the function F defined in Eq. (23) with the exception of a negative sign, the identical approximation procedure as outlined above in this section for the calculation of $\left(T,p\right)\mapsto \mathit{\theta }$ may be applied mutatis mutandis to calculate the transformation $\left(\mathit{\theta },p\right)\mapsto T$.

In any case, a certain effort is required to implement the new formulation of the potential temperature in an atmospheric model, as this equation should be based on the implicit definition Eq. (22), and such a goal may be the subject of future endeavours.

5.2 Approximations of the reference potential temperature

Of course, the previously described procedure to compute the potential temperature may appear to be anything but practical. Indeed, due to the complications inherent with

• the requirement to numerically solve the integral in the function F(x) and

• the need to use Newton's method for an iteration sequence to approach the zero of F(x),

a convenient approach to re-assess the conventional definition of the potential temperature is not provided at all. This motivates the development of a more practical approximation of the reference potential temperature. To arrive at a practicable approximation procedure, the two principal steps in the suggested procedure are briefly outlined in the following, whereas the comprehensive details and intermediate derivation steps are found in Appendix C.

Proceeding from the definition (23) of the function F(x), the computation of the integral $\underset{x}{\overset{T}{\int }}\frac{c_p^{\mathrm{0}}\left(T^{\prime }\right)}{T^{\prime }}\mathrm{d}{T}^{\prime }$ becomes the first obstacle to a practical approximation. Therefore, a plausible initial step is to replace the integral by an expression that is easier to treat. This expression may be proposed as f(T)−f(x), where the function f is defined as $f\left(x\right)=b_{\mathrm{0}}+b_{\mathrm{1}}\ln \left(x-b_{\mathrm{2}}\right)+b_{\mathrm{3}}x+b_{\mathrm{4}}{x}^{\mathrm{2}}$ and which is recognisable as an approximated primitive of $\frac{c_p^{\mathrm{0}}\left(T^{\prime }\right)}{T^{\prime }}$; see Appendix C1. The choice of the functional form of f is motivated by the exact primitive of the integral in the case of a constant c_p.

As previously discussed (cf. Sect. 5.1), the formulation of a new expression for the potential temperature based on the temperature-dependent specific heat capacity c_p(T) requires finding the zero of the equation 0=F(x), where the function F(x) is defined in Eq. (23). Replacing the exact integral in Eq. (23) by the difference f(T)−f(x) means that F(x) is substituted by the function

$$\begin{array}{}\text{(28)}& \widehat{F}\left(x\right)=f\left(T\right)-f\left(x\right)-R_{\mathrm{a}}\ln \left({\displaystyle \frac{p}{p_{\mathrm{0}}}}\right).\end{array}$$

Consequently, the resulting approximated reference potential temperature, i.e. the respective zero of the function $\widehat{F}\left(x\right)$, is denoted as ${\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}$.

The difference between the approximation result and the reference, i.e.

$$\begin{array}{}\text{(29)}& {\mathit{\theta }}_{\mathrm{ref}}-{\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}},\end{array}$$

is then referred to as the basic error of the approximation. Note that the replacement of the function F by $\widehat{F}$ only circumvents the integration in F; the root-finding problem $\mathrm{0}=\widehat{F}\left(x\right)$ for the approximated reference potential temperature ${\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}$ remains analytically not solvable.

Therefore, the second move towards a practical approximation procedure is to construct approximations θ^(k) to the zero of $\widehat{F}\left(x\right)$ by using Newton's method; see Appendix C2. Newton's method is an iterative procedure; the notation θ^(k) refers to the kth computed iterate. Hence, θ^(k) constitutes an approximation to ${\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}$, and, in the limit k→∞, the approximation error

$$\begin{array}{}\text{(30)}& {\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}-{\mathit{\theta }}^{\left(k\right)}\end{array}$$

vanishes. Two formulations of Newton's method are distinguished in Appendix C2, i.e. the principal application of Newton's method, and its further derivative, called Householder's method. Both formulations require the stipulation of one of the iterates θ^(k) as sufficient to obtain a result of appropriate accuracy. The higher the number of iterations, of course, the smaller is the error (30), whereas the basic error (29) remains unaffected by the number of iterations. Hence, in any case, the basic error (29) is to be accepted as at least implied in the final approximation, even though a well-chosen θ^(k) could result in an approximation error θ_ref−θ^(k) that is smaller than the basic error.

Figure 5 Absolute basic error $\mathrm{\Delta }\mathit{\theta }={\mathit{\theta }}_{\mathrm{ref}}-{\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}$, cf. Eq. (29), from approximating the reference potential temperature along the US Standard Atmosphere (a) and for the extended pressure range 1000 to 0.5 hPa and temperature range 180 to 300 K (b). For orientation, the white solid line indicates the p–T profile from the US Standard Atmosphere. The relative basic error |Δθ|∕θ_ref is shown in panel (c) for the extended pressure and temperature range.

The various errors implied in the proposed approximation procedure combining for the approximation's total error, as well as accompanying details, are discussed in Appendix D. In brief, Fig. 5a illustrates the basic error (29) based on the pressure and temperature profiles of the US Standard Atmosphere, as these provide atmospherically meaningful averages of realistic temperature–pressure data pairs. Based on the parameters of the US Standard Atmosphere, the basic error inherent with the approximation remains below 1.25 K up to altitudes of 50 km. Thus, regarding the subsequent iteration process, a substantial improvement of the error compared to ∼1.5 K is not to be expected for the total error of approximating the reference potential temperature.

An error analysis exclusively based on the US Standard Atmosphere is constrained to specific combinations of the air's pressure and temperature, potentially suppressing latent errors that may emerge if certain fluctuations of the real atmosphere's temperature and pressure profiles are considered. Thus, the error analysis is extended to an atmospheric pressure (p) and temperature (T) range, from 1000 to 0.5 hPa and from 180 to 300 K, such that the conditions within the entire troposphere and stratosphere, including the stratopause, are covered. Figure 5b illustrates the absolute basic error (29) for the extended ranges of pressure and temperature while Fig. 5c illustrates the relative basic error $\left|{\mathit{\theta }}_{\mathrm{ref}}-{\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}\right|/{\mathit{\theta }}_{\mathrm{ref}}$. The contours in Fig. 5b and c mainly highlight two regions: at ∼100 hPa where Δθ never rises above 0.75 K, which corresponds to a maximum relative basic error of 0.15 %, and in a pressure range from ∼5 to 1 hPa where a Δθ of 1.25 K is never exceeded, corresponding to relative errors of at most 0.1 %. Note that the entire Δθ scale ranges up to 3 K, which may only be reached at pressures below 0.8 hPa combined with temperatures above 280 K.

Figure 6 (a) Relative error $\mathrm{\Delta }\mathit{\theta }/{\mathit{\theta }}_{\mathrm{ref}}=\left|{\mathit{\theta }}^{\left(\mathrm{2}\right)}-{\mathit{\theta }}_{\mathrm{ref}}\right|/{\mathit{\theta }}_{\mathrm{ref}}$ of the second iterate θ⁽²⁾, obtained with Newton's method for the ranges of pressure and temperature from 1000 to 0.5 hPa and from 180 to 300 K, respectively. Panels (b) and (c) exhibit the difference $\mathrm{\Delta }\mathit{\theta }=\left|{\mathit{\theta }}_{\mathrm{1005}}-{\mathit{\theta }}_{\mathrm{ref}}\right|$ and relative difference Δθ∕θ_ref, respectively, on a logarithmic scale between the reference potential temperature θ_ref and the potential temperature θ₁₀₀₅ based on a constant specific heat capacity (c_p=1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$). For orientation, the white solid line indicates the p–T profile from the US Standard Atmosphere.

As previously discussed, the basic error is unavoidable and is to be accepted when applying the suggested substitution for the integral in the definition of the function F(x) in Eq. (23). However, as outlined in Appendix C2, the second iterate θ⁽²⁾ of Newton's method (principal application) may thoroughly suffice for the final approximation to the reference potential temperature θ_ref, as this iteration level already features an approximation error (30) which is negligibly small. Figure 6a illustrates the total relative error of the suggested approximation θ⁽²⁾ with respect to the ultimate reference θ_ref for the extended ranges of pressure and temperature. Indeed, the contour pattern in Fig. 6a and the basic relative approximation error shown in Fig. 5c are remarkably similar. Thus, the iteration process itself imparts only a minor contribution to the total error compared to the basic approximation error.

The total approximation error, which is

$$\begin{array}{}\text{(31)}& {\mathit{\theta }}_{\mathrm{ref}}-{\mathit{\theta }}^{\left(\mathrm{2}\right)}=\left({\mathit{\theta }}_{\mathrm{ref}}-{\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}\right)+\left({\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}-{\mathit{\theta }}^{\left(\mathrm{2}\right)}\right),\end{array}$$

is dominated by the unavoidable basic error (first bracket) and augmented by a negligible error inherent to the iteration (second bracket), also supporting the conclusion that the second iterate of Newton's method is an appropriate approximation procedure. Figure 7 presents stepwise instructions for the computation of the second iterate approximation to the reference potential temperature and may serve as a guide to follow the numerous equations and intermediate analytical steps described throughout the derivations in Appendix C.

Figure 7 Flowchart guiding through the process of computing the approximation θ⁽²⁾ by using Newton's formulation (C5) until its second iteration, wherein T (in K) and p (in hPa) are the atmospheric air conditions in terms of temperature and pressure, respectively, and p₀ is set to 1000 hPa (WMO, 1966). Table C1 collects values of θ_ref and the approximation θ⁽²⁾ together with intermediate results for selected pairs of temperature and pressure to verify a computation according to this instruction.

For completeness, Fig. 6b and c exhibit a final comparison by means of the logarithmic difference and the logarithmic relative difference between the reference potential temperature θ_ref and the conventional definition ${\mathit{\theta }}_{c_p}$ (WMO, 1966) based on a constant specific heat capacity c_p=1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$. Notably, over a wide altitude range within the troposphere (i.e. for atmospheric pressures greater than ∼100 hPa), the absolute error $\mathrm{\Delta }\mathit{\theta }=\left|{\mathit{\theta }}_{\mathrm{1005}}-{\mathit{\theta }}_{\mathrm{ref}}\right|$ remains below 1 K, cf. Fig. 6b, corresponding to a relative error Δθ∕θ_ref of at most 0.1 %. However, in the pressure range below ∼100 hPa, deviations of the real atmospheric conditions from those of the US Standard Atmosphere could increase the absolute error Δθ from a few K to up to 10 K, corresponding to an increase in the relative error to 1 %. Further critical pressure levels are at ∼20 and ∼5 hPa, where the error's magnitude increases to several tens and several hundreds of K, respectively. At a pressure of 0.5 hPa, an absolute error Δθ of up to 500 K is reached, which corresponds to a relative error of 10 % or even more.

5.3 Implementation aspects

The use of the new reference potential temperature θ_ref in a numerical model requires additional computational effort to perform corresponding calculations. Hereafter, two aspects are briefly discussed: (i) the formulation of the model equations, which include θ_ref, and (ii) the calculation of θ_ref.

Although it is beyond the scope of the present study to provide a general derivation of an appropriate energy equation based on θ_ref for atmospheric models, a formulation of the total derivative of θ_ref is given by

$$\begin{array}{}\text{(32)}& c_p\left({\mathit{\theta }}_{\mathrm{ref}}\right){\displaystyle \frac{\mathrm{d}{\mathit{\theta }}_{\mathrm{ref}}}{{\mathit{\theta }}_{\mathrm{ref}}}}=c_p\left(T\right){\displaystyle \frac{\mathrm{d}T}{T}}-R_{\mathrm{a}}{\displaystyle \frac{\mathrm{d}p}{p}},\end{array}$$

where the details of its derivation are given in Appendix E. The total derivative of θ_ref may be useful, since the governing equations are commonly formulated as differential equations.

The calculation of both the reference potential temperature θ_ref and its approximation ${\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}$ on the basis of given values of pressure p and temperature T requires an iterative procedure. The additional computational effort inherent with these calculations depends on the number of iterations. If, however, the second iteration θ⁽²⁾ already represents an appropriate approximation of θ_ref (cf. Sect. 5.2), then the flowchart in Fig. 7 immediately conveys the additional computational effort to be expected. The calculation of the starting value x₀ is identical to computing θ₁₀₀₅. An additional effort results from the evaluation of the functions f (three times) and f^′ (two times), respectively, and the combination (two times) of obtained values to determine x₁ and x₂. Since each of these evaluations causes additional numerical steps, the computational effort to obtain θ⁽²⁾ is in total about 7 times more than the calculation of the conventional θ₁₀₀₅, while the algorithmic complexity is constant.

6 The potential temperature for air as a real gas

To account for real-gas effects (that cause a behaviour other than that of an ideal gas; cf. Sect. 4) on the potential temperature, we use the model embedded in REFPROP (Lemmon et al., 2018), a standard reference database from NIST. This model treats air as a mixture and employs state-of-the-art reference equations of state for pure nitrogen (Span et al., 2000), oxygen (Schmidt and Wagner, 1985), and argon (Tegeler et al., 1999). The mixing rule and binary interaction parameters are taken from the GERG-2008 model (Kunz and Wagner, 2012). From its definition in terms of an isentropic process, the potential temperature θ_real(T, p) is defined implicitly by

$$\begin{array}{}\text{(33)}& s({\mathit{\theta }}_{\mathrm{real}},p_{\mathrm{0}})=s(T,p),\end{array}$$

where s is the specific entropy. Calculating θ_real(T, p) is a two-step process. First, the specific entropy s is computed at temperature T and pressure p. Then, the temperature θ_real that gives the same entropy s at the ground pressure p₀ is found. This is an iterative calculation, but it is accomplished automatically within the REFPROP software (Lemmon et al., 2018).

One caveat should be mentioned regarding the computed potential temperatures. The range of validity of the equations of state for the air components (Span et al., 2000; Schmidt and Wagner, 1985; Tegeler et al., 1999) extends only up to 2000 K. At very high altitudes, computed values of θ_real exceed this limit. While all the equations extrapolate in a physically realistic way, their quantitative accuracy is less certain above 2000 K. This caveat also applies to the ideal-gas calculations; the correlations for $c_p^{\mathrm{0}}\left(T\right)$ for N₂ and O₂ are extrapolations beyond 2000 K. However, since the same ideal-gas values are used in the real-gas calculations, any inaccuracy in $c_p^{\mathrm{0}}\left(T\right)$ will cancel when evaluating the difference between ideal-gas and real-gas values of θ.

Figure 8 Difference θ_real−θ_ref reflecting the deviation of the potential temperature θ_real, based on the properties of air behaving as a real gas under variable temperature and pressure, from the herein derived potential temperature expression θ_ref for the ideal-gas limit of the air's specific heat capacity $c_p^{\mathrm{0}}\left(T\right)$. (a) Difference along the profile of the US Standard Atmosphere. (b) Relative difference in p–T coordinates covering any combination of atmospherically relevant temperatures and pressures.

Figure 8 illustrates the comparison between the real-gas potential temperature θ_real and the ideal-gas reference potential temperature θ_ref. Figure 8a shows the difference θ_real−θ_ref along the p–T profile of the US Standard Atmosphere and Fig. 8b accounts for any p–T combination of extended range but shows the relative difference instead. The difference between θ_real and θ_ref never exceeds 0.1 K for the absolute difference or $\mathrm{30}\times {\mathrm{10}}^{-\mathrm{5}}=\mathrm{0.03}\mathrm{\%}$ for the relative difference. As may be anticipated from the deviation of $c_p^{\mathrm{0}}$ shown in Fig. 3 at low temperatures both from the experimentally determined values (which may be inaccurate) and from the REFPROP data, the real-gas effect on the specific heat capacity of dry air tends to increase towards the coldest gas temperatures. However, the difference between the real- and ideal-gas approaches results in essentially no substantial difference between the resulting θ's, neither at ground conditions (for any temperature at ∼1000 hPa) nor at very high altitudes (at pressures below ∼1 hPa). While the negligible difference between θ_real and θ_ref near ground levels is less surprising, the diminished difference at higher altitudes reflects that in this region the potential temperature reaches such high values that the difference between the real-gas and the ideal-gas specific heat capacity becomes insignificant. Within the intermediate (stratospheric) region, the low pressures (and thus the low air densities) cause the ideal-gas assumption to be an accurate approximation even at low temperatures. In general, the degree to which a gas can be treated as ideal is primarily a function of the (molar) density. For an ideal gas, the density is proportional to the quotient $\frac{p}{T}$; this is almost true also for real air. Hence, declining pressures together with rising temperatures both make the air's behaviour increasingly close to ideal.

7 Implications on the use of the potential temperature

As previously shown, the newly defined reference potential temperature θ_ref deviates most from the WMO-defined potential temperature θ₁₀₀₅ at stratospheric altitudes and above (cf. Fig. 6). More particularly, not only do the values from both θ definitions differ, but also their vertical derivatives, i.e. $\frac{\partial {\mathit{\theta }}_{\mathrm{ref}}}{\partial z}$ and $\frac{\partial {\mathit{\theta }}_{\mathrm{1005}}}{\partial z}$. Whether such deviations have a significant effect on an application is very case dependent and requires detailed examination and specific appraisal. Below, four typical applications of the potential temperature were selected and are examined regarding the quantitative effect on the results due to deviations of the introduced reference potential temperature compared to the conventional and commonly used θ₁₀₀₅. The purpose of this examination is to document the magnitude of errors to allow a well-founded, individual decision for each application of the potential temperature as to whether it is worth applying the more rigorous calculation in the particular context.

7.1 The Brunt–Väisälä frequency

The formula for the (squared) Brunt–Väisälä frequency N² is often given in the form of Eq. (2), i.e. a formula involving the potential temperature θ. The substitution of θ in Eq. (2) by the new reference potential temperature θ_ref may be tempting, but it is erroneous, and the resulting quantity is denoted as $N_{\mathrm{false}}^{\mathrm{2}}$. The Brunt–Väisälä frequency is not defined by Eq. (2), since this formula results from various simplifications in its derivation, e.g. by assuming hydrostatic conditions and a constant specific heat capacity. Consequently, the substitution of θ_ref in Eq. (2) leads to a wrong formula for the Brunt–Väisälä frequency that does not correctly consider the temperature dependence of dry air's specific heat capacity.

The Brunt–Väisälä frequency is the oscillation frequency of an air parcel due to a local density perturbation (see e.g. Durran and Klemp, 1982; Marquet and Geleyn, 2013; Wallace and Hobbs, 2006; Ambaum, 2010). Retaining the assumption of hydrostatic conditions, the defining formula yields

$$\begin{array}{}\text{(34)}& N^{\mathrm{2}}={\displaystyle \frac{g}{T}}\left({\displaystyle \frac{\partial T}{\partial z}}+{\displaystyle \frac{g}{c_p\left(T\right)}}\right),\end{array}$$

where the temperature-dependent specific heat capacity c_p(T) was implied, and which quantifies the balance between the actual temperature stratification $\frac{\partial T}{\partial z}$ and the dry-adiabatic lapse rate $-\frac{g}{c_p\left(T\right)}$ (e.g. Holton, 2004).

Figure 9 Vertical profiles of (a) the temperature up to 85 km altitude as typical for mid-latitudes in June (red curve) and December (blue curve). (b) Resulting wrong Brunt–Väisälä frequency $N_{\mathrm{false}}^{\mathrm{2}}$ (dashed lines) and the true Brunt–Väisälä frequency N² (solid lines) for the two temperature profiles from panel (a).

To illustrate the deviation of $N_{\mathrm{false}}^{\mathrm{2}}$ from N², vertical profiles of both variables were calculated based on the temperature profiles shown in Fig. 9a. The temperature data are taken from the Upper Atmosphere Research Satellite Reference Atmosphere Project (URAP; see Swinbank and Ortland, 2003) data and are assumed as typical at mid-latitudes during June and December. The temperature profiles extend up to altitudes of 85 km and thus cover the entire stratosphere and most of the mesosphere. The hydrostatic assumption allowed for computing pressure profiles along the URAP values for the vertical temperature distribution. Subsequently, the reference potential temperature θ_ref and its vertical derivative were calculable. The resulting vertical profiles for $N_{\mathrm{false}}^{\mathrm{2}}$ and the true Brunt–Väisälä frequency N² are shown in Fig. 9b. Evidently, the values of $N_{\mathrm{false}}^{\mathrm{2}}$ (dashed lines) deviate significantly from N² (solid lines) and increasingly so towards higher altitudes above 15 km. However, the absolute deviation $|N^{\mathrm{2}}-N_{\mathrm{1005}}^{\mathrm{2}}|$, using $N_{\mathrm{1005}}^{\mathrm{2}}$ as calculated with θ₁₀₀₅ in accordance with Eq. (2), does not exceed $\mathrm{1.6}\times {\mathrm{10}}^{-\mathrm{6}}{\mathrm{s}}^{-\mathrm{2}}$ (not shown), indicating that $N_{\mathrm{1005}}^{\mathrm{2}}$ is a good representation of N² along these temperature profiles.

For equations involving the potential temperature, however, it should be emphasised that the substitution of θ by θ_ref rarely succeeds and that instead the entire derivation of the equations requires careful consideration of the assumptions, such as the constancy of c_p, to avoid aberrations and erroneous conclusions.

7.2 The potential vorticity

Ertel's potential vorticity (e.g. Ertel, 1942; Hoskins et al., 1985; Schubert et al., 2004; Holton, 2004) may be defined as the potential vorticity of the dry-air potential temperature by

$$\begin{array}{}\text{(35)}& \mathrm{PV}\left(\mathit{\theta }\right)={\displaystyle \frac{\mathrm{1}}{\mathit{\rho }}}\left(\mathrm{2}\mathbf{\Omega }+\mathrm{\nabla }\times \mathit{u}\right)\cdot \mathrm{\nabla }\mathit{\theta }.\end{array}$$

In this definition, $\mathrm{2}\mathbf{\Omega }+\mathrm{\nabla }\times \mathit{u}$ is the absolute vorticity, Ω denotes Earth's angular velocity, u the three-dimensional wind vector, and ρ the air density (see e.g. Hoskins et al., 1985; Cotton et al., 2011; Marquet, 2014). Since Eq. (35) represents the defining equation for Ertel's potential vorticity, the two potential vorticities

$$\begin{array}{}\text{(36)}& \begin{array}{rl}{\mathrm{PV}}_{\mathrm{ref}}& =\mathrm{PV}\left({\mathit{\theta }}_{\mathrm{ref}}\right),\\ {\mathrm{PV}}_{\mathrm{1005}}& =\mathrm{PV}\left({\mathit{\theta }}_{\mathrm{1005}}\right)\end{array}\end{array}$$

based on the new reference potential temperature θ_ref and θ₁₀₀₅, respectively, are considered. To provide a first comparison of these potential vorticities, u=0 is assumed, i.e. an atmosphere at rest. Additionally, the potential temperature is assumed as horizontally constant. Consequently, Eq. (35) reduces to

$$\begin{array}{}\text{(37)}& \mathrm{PV}\left(\mathit{\theta }\right)={\displaystyle \frac{\mathrm{2}\sin \left(\mathit{\varphi }\right)}{\mathit{\rho }}}{\displaystyle \frac{\mathrm{2}\mathit{\pi }}{t_{\mathrm{E}}}}{\displaystyle \frac{\partial \mathit{\theta }}{\partial z}}\end{array}$$

for a position on Earth with geographical latitude ϕ and t_E=24 h, the duration of one rotation of the Earth.

Figure 10 (a) Vertical profiles of the potential vorticity PV_ref computed with θ_ref (solid lines), and PV₁₀₀₅ computed with θ₁₀₀₅ (dashed lines), for an atmosphere at rest along the temperature profiles from Fig. 9a for June (red lines) and December (blue lines). Since the temperature profiles are representative for mid-latitudes on the Northern Hemisphere, the geographical latitude in Eq. (37) was set to 52^∘ N. (b) Relative deviation $|{\mathrm{PV}}_{\mathrm{ref}}-{\mathrm{PV}}_{\mathrm{1005}}|/{\mathrm{PV}}_{\mathrm{ref}}$ of the potential vorticity profiles from panel (a).

Using the temperature profiles from Fig. 9a together with the values of the potential temperatures θ_ref and θ₁₀₀₅, the evaluation of the two potential vorticities (36) by (37) yields the potential vorticity profiles shown in Fig. 10a while their relative deviations are shown in Fig. 10b. Since the temperature profiles are representative for the northern-hemispheric mid-latitudes, the geographical latitude ϕ in Eq. (37) was set to 52^∘ N. At tropospheric altitudes, the relative deviation between θ_ref and θ₁₀₀₅ is small and never exceeds ∼1 %, while it continuously increases towards higher altitudes. According to these profiles, the relative deviation exceeds 10 % at 30 km and reaches 100 % at the highest altitudes.

It is noteworthy, however, that the computations of both N² (cf. Sect. 7.1 and Fig. 9b) and PV (Fig. 10b) are based on the specific temperature profiles from URAP (cf. Sect. 7.1 and Fig. 9a) and thus are not of general validity. The selection of these temperature profiles was entirely arbitrary and exclusively aimed at illustrating possible implications of the use of the developed reference potential temperature. The resulting and indicated deviations are ultimately subject to individual assessment on applying θ_ref.

7.3 Vertical sorting of data

For atmospheric investigations, e.g. in the region of the upper troposphere and lower stratosphere (UT/LS), it is common practice to set vertical profiles of atmospheric parameters in relation to the potential temperature as a vertical coordinate. This way, the increasingly isentropic stratification of the atmosphere above the UT is taken into account. The transport of an air mass along isentropic surfaces, i.e. surfaces of constant potential temperature and entropy, is to be regarded as adiabatic. Hence, the air's composition and properties within the same isentrope interval, regardless of the observation location, are better comparable than they would be if based on other isopleths (i.e. height or pressure coordinates). Investigations of air mass compositions over time and from different regions at the same θ level largely exclude that, during its transport history, the air had experienced vertical displacement and/or diabatic processes (radiative heating, condensation/evaporation) which would result in energy conversion. The tropopause height is often used as a reference height in the θ coordinate system in connection with the vertical sorting of observational data, whereby the assignment of tropospheric and stratospheric processes is made or exchanges across the tropopause are investigated (Holton et al., 1995; Stohl et al., 2003). Consequently, the tropopause height is also determined by the potential vorticity (e.g. Gettelman et al., 2011, and cf. Sect. 7.2) if the conventional tropopause definitions (cold point or lapse rate, WMO, 1957) do not allow for clearly determining the tropopause height, e.g. in the Asian monsoon anticyclone (cf. Höpfner et al., 2019) or in the polar winter vortex (Wilson et al., 1989; Weigel et al., 2014). The conventional definition of θ implies a systematic error in the vertical sorting of observational data in the θ coordinate system, independent of the measurement platform. Investigations with high-altitude research aircraft such as the G-550 HALO (e.g. Wendisch et al., 2016; Voigt et al., 2017), the NASA WB-57 or ER-2 (e.g. Murphy et al., 2007; Dessler, 2002), the M-55 Geophysica (Curtius et al., 2005; Borrmann et al., 2010; Frey, 2011), balloon-borne platforms (Lary et al., 1995; Vernier et al., 2018), or satellite-based vertical profiles (e.g. Davies et al., 2006; Spang et al., 2005) require consideration of the systematic error in θ if calculated as ${\mathit{\theta }}_{c_p}$ in compliance with the definition by the WMO (1966). The possibly inconsistent use of a constant c_p value of 1004 or 1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ (or any other) in different and compared data sets, which could be due to different literature references for this value (cf. Table 1), will not be explored here. At altitudes between 15 and 20 km (ceiling of high-altitude research aircraft), an overestimation by about 0.1 %–0.5 % is to be expected for the potential temperature according to the conventional definition (cf. Fig. 4b). At altitudes of 30–35 km, an overestimation by up to 2 %–5 % results. Whether this error is significant or small compared to the uncertainty of ambient temperature and pressure measurement aboard the respective aircraft is left to individual judgement in the course of data processing. In the case of spacecraft-bound vertical soundings (e.g. from ASTRO-SPAS, SCIAMACHY, or ENVISAT), the error in the potential temperature determined by ${\mathit{\theta }}_{c_p}$ exceeds 10 % at altitudes above 40 km, as shown in Fig. 4b. Finally, we note that the specified errors apply exclusively along the vertical profile of the US standard atmosphere and that deviations of the actual temperature profile from the US standard atmosphere, e.g. warmer temperatures, could lead to larger errors (cf. Fig. 6).

7.4 Diabatic heating rates

Diabatic heating rates refer to the rate of energy $\frac{\mathrm{d}q}{\mathrm{d}t}$ supplied to a given air parcel, e.g. by radiative heating, and are given in units of $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{s}}^{-\mathrm{1}}$. This energy supply causes a temperature change in an air parcel at a rate which hereafter is referred to as the absolute heating rate,

$$\begin{array}{}\text{(38)}& \begin{array}{rl}{\mathrm{AHR}}_{\mathrm{ref}}\left({\displaystyle \frac{\mathrm{d}q}{\mathrm{d}t}}\right)& ={\displaystyle \frac{\mathrm{d}T}{\mathrm{d}t}}={\displaystyle \frac{\mathrm{1}}{c_p^{\mathrm{0}}\left(T\right)}}{\displaystyle \frac{\mathrm{d}q}{\mathrm{d}t}},\\ {\mathrm{AHR}}_{\mathrm{1005}}\left({\displaystyle \frac{\mathrm{d}q}{\mathrm{d}t}}\right)& ={\displaystyle \frac{\mathrm{d}T}{\mathrm{d}t}}={\displaystyle \frac{\mathrm{1}}{\mathrm{1005}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}}}{\displaystyle \frac{\mathrm{d}q}{\mathrm{d}t}}.\end{array}\end{array}$$

Again, the distinction was made between the temperature-dependent $c_p^{\mathrm{0}}\left(T\right)$ and the constant c_p=1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$ specific heat capacity. From the defining Eq. (38), the relative difference between these absolute heating rates, where x designates an arbitrary diabatic heating rate, is

$$\begin{array}{}\text{(39)}& {\displaystyle \frac{{\mathrm{AHR}}_{\mathrm{1005}}\left(x\right)-{\mathrm{AHR}}_{\mathrm{ref}}\left(x\right)}{{\mathrm{AHR}}_{\mathrm{ref}}\left(x\right)}}={\displaystyle \frac{c_p^{\mathrm{0}}\left(T\right)}{\mathrm{1005}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}}}-\mathrm{1}.\end{array}$$

Apart from the absolute heating rates for the change in absolute temperature, the change in potential temperature due to a diabatic heating rate $\frac{\mathrm{d}q}{\mathrm{d}t}$ is of interest. For example, it is the change in potential temperature that modifies the altitude of modelled trajectories in Lagrangian chemical transport models based on isentropic coordinates rather than the change in absolute temperature (e.g. the SLIMCAT model, Chipperfield, 2006; or CLaMS model, Pommrich et al., 2014).

Taking the relation T ds=dq for the specific entropy into account, Gibbs' Eq. (8) may be rewritten as

$$\begin{array}{}\text{(40)}& {\displaystyle \frac{\mathrm{d}q}{T}}={\displaystyle \frac{c_p\left(T\right)}{T}}\mathrm{d}T-R_{\mathrm{a}}{\displaystyle \frac{\mathrm{d}p}{p}}.\end{array}$$

Comparing the right-hand side of this equation to the total derivative of the new reference potential temperature θ_ref (see Appendix E for the detailed computation and Eq. E6 for the result) Eq. (40) amounts to

$$\begin{array}{}\text{(41)}& {\displaystyle \frac{\mathrm{d}q}{T}}=c_p\left({\mathit{\theta }}_{\mathrm{ref}}\right){\displaystyle \frac{\mathrm{d}{\mathit{\theta }}_{\mathrm{ref}}}{{\mathit{\theta }}_{\mathrm{ref}}}}.\end{array}$$

Consequently, the following two diabatic heating rates,

$$\begin{array}{}\text{(42)}& \begin{array}{rl}{\displaystyle \frac{\mathrm{d}{\mathit{\theta }}_{\mathrm{ref}}}{\mathrm{d}t}}& ={\displaystyle \frac{{\mathit{\theta }}_{\mathrm{ref}}}{c_p^{\mathrm{0}}\left({\mathit{\theta }}_{\mathrm{ref}}\right)T}}{\displaystyle \frac{\mathrm{d}q}{\mathrm{d}t}}={\mathrm{HR}}_{\mathrm{ref}}\left({\displaystyle \frac{\mathrm{d}q}{\mathrm{d}t}}\right),\\ {\displaystyle \frac{\mathrm{d}{\mathit{\theta }}_{\mathrm{1005}}}{\mathrm{d}t}}& ={\displaystyle \frac{{\mathit{\theta }}_{\mathrm{1005}}}{\left(\mathrm{1005}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}\right)\cdot T}}{\displaystyle \frac{\mathrm{d}q}{\mathrm{d}t}}={\mathrm{HR}}_{\mathrm{1005}}\left({\displaystyle \frac{\mathrm{d}q}{\mathrm{d}t}}\right),\end{array}\end{array}$$

for the potential temperatures θ_ref and θ₁₀₀₅ may be defined. Denoting again by x an arbitrary diabatic heating rate, the relative difference between the heating rates (Eq. 42) is

$$\begin{array}{}\text{(43)}& {\displaystyle \frac{{\mathrm{HR}}_{\mathrm{1005}}\left(x\right)-{\mathrm{HR}}_{\mathrm{ref}}\left(x\right)}{{\mathrm{HR}}_{\mathrm{ref}}\left(x\right)}}={\displaystyle \frac{{\mathit{\theta }}_{\mathrm{1005}}}{{\mathit{\theta }}_{\mathrm{ref}}}}{\displaystyle \frac{c_p^{\mathrm{0}}\left({\mathit{\theta }}_{\mathrm{ref}}\right)}{\mathrm{1005}\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}}}-\mathrm{1}.\end{array}$$

Figure 11 (a) Monthly averaged temperatures profiles for 52^∘N. (b) The relative differences between the absolute heating rates, defined in Eq. (39). (c) The relative differences between the heating rates (Eq. 43) for the potential temperatures θ₁₀₀₅ and θ_ref. (d) The resulting relative deviation of the potential temperatures ${\mathit{\theta }}_{\mathrm{1005}}^{\ast },{\mathit{\theta }}_{\mathrm{ref}}^{\ast }$ after 24 h of heating with constant diabatic heating $\frac{\mathrm{d}q}{\mathrm{d}t}$ and the resulting heating rates HR_ref, HR₁₀₀₅ at constant pressure.

In order to judge the magnitudes of the relative differences (Eqs. 39 and 43), the monthly averaged temperature profiles from ERA-Interim (Dee et al., 2011) data for 52^∘ N geographical latitude are used; see Fig. 11a. The relative differences of the absolute heating rates (Eq. 39) are shown in Fig. 11b, and the difference appears to be small. However, the relative differences of the heating rates (Eq. 43) in Fig. 11c are much larger, as relative deviations exceeding 50 % are reached in the upper stratosphere and lower mesosphere (at pressures below 1 hPa). Additionally, the temperatures that resulted after 24 h of heating with a constant heating rate $\frac{\mathrm{d}q}{\mathrm{d}t}$ as given in the (averaged) data set were computed, where a constant pressure is assumed for simplicity. As may be anticipated from the small deviations in Fig. 11b, the differences in the final absolute temperatures by using the absolute heating rates AHR₁₀₀₅ or AHR_ref are smaller than 0.044 K. However, the differences in the potential temperatures ${\mathit{\theta }}_{\mathrm{1005}}^{\ast },{\mathit{\theta }}_{\mathrm{ref}}^{\ast }$, computed with the heating rates HR_ref, HR₁₀₀₅, are much larger (Fig. 11d) and amount to about 3 % at 10 hPa and about 15 % at 1 hPa. For transport calculations done in isentropic coordinates, these differences are of the same order of magnitude as the deviations resulting from the use of the temperature-dependent instead of the constant c_p. It remains to be decided on individual application whether this additional effect in the calculation is significant.

A standard diagnostic for the speed of the stratospheric circulation is the time lag of the upward-propagating seasonal signal in tropical stratospheric water vapour (the so-called tape recorder, Mote et al., 1996). Here, differences between calculations (done in isentropic coordinates) based on different current meteorological reanalysis data sets amount to about 10 %–30 % for the signal's upward propagation below about 10 hPa (Tao et al., 2019), such that the additional deviation from using the temperature-dependent c_p is comparably small. However, in cases of smaller inter-model differences the additional c_p-related uncertainty needs to be assessed.

Note that the determination of absolute temperatures $T_{\mathrm{1005}}^{\ast },T_{\mathrm{ref}}^{\ast }$ which correspond to the resulting potential temperatures ${\mathit{\theta }}_{\mathrm{1005}}^{\ast },{\mathit{\theta }}_{\mathrm{ref}}^{\ast }$ after 24 h differs by less than 0.014 K (not shown).

8 Summary and conclusions

Under the assumption that dry air is an ideal gas, a re-assessment of computing the potential temperature was introduced that accounts for the hitherto unconsidered temperature dependence of air's specific heat capacity. The new reference potential temperature θ_ref was introduced, which is thermodynamically consistent and based on a state-of-the-art parameterisation of the ideal-gas specific heat capacity of dry air from the National Institute of Standards and Technology (NIST). This reference potential temperature was compared to a potential temperature θ_real wherein the real-gas behaviour of dry air is considered. In the range of temperatures from 180 to 300 K and the range of pressures from 1000 to 0.5 hPa, covering the atmospheric conditions of roughly the entire troposphere and stratosphere, the relative differences between θ_ref and θ_real are smaller than 0.03 % and may be considered negligible. Consequently, θ_ref even provides a reasonable approximation to the potential temperature of the real gas.

The difference between the newly derived reference potential temperature θ_ref and the conventionally determined potential temperature ${\mathit{\theta }}_{c_p}$ (with constant c_p=1005 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$, as recommended by the World Meteorological Organisation, WMO, 1966) increases with altitude, e.g. Δθ≥1 K at pressures p≤60 hPa.

Derivation of a potential temperature that is consistent with thermodynamics and that accounts for the ideal-gas properties of dry air requires the integration of Gibbs' equation and the subsequent solution of the resulting nonlinear equation. With a constant c_p, both analytical steps are straightforward, resulting in the conventional expression (13) as suggested by WMO (1966). However, if instead the temperature dependence of air's specific heat capacity c_p(T) is considered, the integrals as well as the equations are not analytically solvable, and, thus, the solution must be approximated. Both approximations were performed and described in detail. The integral was treated with the basic approximation and the solution of the nonlinear equation was approximated by the second iterate of Newton's method. As an alternative to Newton's classical method, a modified formulation of Householder's iteration method is provided, featuring accelerated convergence properties.

The suggested approximation steps to obtain a reference potential temperature have two main sources of error: the error ${\mathit{\theta }}_{\mathrm{ref}}-{\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}$ inherent in the integral's basic approximation and the error ${\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}-{\mathit{\theta }}^{\left(k\right)}$ of the kth Newton iterate. The latter error approaches zero as k→∞, whereas the error resulting from the basic approximation remains well below 0.1 % (along the US Standard Atmosphere) for values of θ_ref of up to ∼2000 K, hence up to stratopause altitudes. To keep this low error level also for θ_ref>2000 K, the approximation may require an extension by means of a higher-order polynomial.

One of the foremost implications of the re-assessed potential temperature's definition concerns the use of θ as a vertical coordinate for the sorting, grouping, and comparison of (measured) data, e.g. along or across isentropes. Thereby, the re-assessed potential temperature constitutes a more accurate consideration of the air's actual properties. This particularly concerns for example the specific heat capacity which is conventionally assumed as constant and for which various values are given depending on the textbook consulted (offering a range from 1000 to 1010 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$; see Table 1).

Significant errors and biases may arise if, for instance, the conventional derivation of θ (WMO, 1966) is used together with values for air's specific gas constant (R_a) or air's specific heat capacity (c_p), which better comply with the most recent state of knowledge. Moreover, the use of the standard pressure 1013.25 instead of 1000 hPa as defined by WMO (1966) and consistently used herein as ground level pressure (p₀) may cause an additional deviation of the resulting θ. Thus, the re-assessment of θ's definition could largely diminish such errors and biases and improve the comparability of data.

In addition to the vertical sorting of data, implications of the new reference potential temperature were discussed for several other applications in which the potential temperature is used. On the one hand, results may appear mostly unaffected by using θ_ref instead of the conventional θ₁₀₀₅, such as the values of the Brunt–Väisälä frequency or the temperature change in air parcels due to diabatic heating. On the other hand, it was illustrated that any formula which involves the potential temperature needs to be carefully reviewed to see if its derivation relies on the assumed constancy of the specific heat capacity. If this is the case, substituting θ_ref for all occurrences of θ within the particular formula may lead to a wrong computation.

In contrast, examples were shown where the computation of Ertel's potential vorticity and the rate of change of potential temperature in response to diabatic heating yields different results by the use of θ_ref instead of θ₁₀₀₅. The differences increased with altitude; hence they become more important for applications within the stratosphere and above.

It should be emphasised that all these examples were based on assuming particular profiles of temperature and pressure together with other assumptions. Moreover, only a limited number of examples could be investigated, while the applications of potential temperature are numerous. Consequently, a well-founded, individual decision is required for each application of the potential temperature as to whether it is worth applying the more rigorous calculation in the particular context.

On the one hand, such a re-assessment could take into account the current state of knowledge regarding the accuracy of thermodynamic variables and substance-related properties. On the other hand, this way, the conceptional abstractness already inherent in θ is not further complicated by a misleading selection of parameters or reputed constants. There is no doubt that the conventional method is suitable for the description of most processes occurring within the troposphere. However, at stratospheric or even mesospheric altitudes, the neglect of the temperature dependence of the ideal-gas heat capacity in the conventional definition increasingly distorts the resulting absolute values as well as the vertical course of the potential temperature. Ultimately, it seems obvious to profit from the computing capacities available today and from the known higher accuracy of physical variables and atmospheric parameters to carry out a reappraisal of the potential temperature, a useful (but not always consistently used) meteorological quantity.

Appendix A: Derivation of the specific heat capacity from thermodynamics

In the following, the derivation of the air's specific heat capacities C_V, C_p (capital letters indicate molar units) at constant volume and pressure, respectively, is summarised, mainly following the textbook exposition by Kondepudi and Prigogine (1998). We start with the ideal-gas law

$$\begin{array}{}\text{(A1)}& pV=NRT,\end{array}$$

with p the pressure, V the volume of the system, N the amount of gas within the volume, T the temperature, and R the universal gas constant. Additionally, the first law of thermodynamics is

$$\begin{array}{}\text{(A2)}& \mathrm{d}U=\mathrm{d}Q-p\mathrm{d}V,\end{array}$$

with U the internal energy of the system and dQ specifying the change in heat. Insertion of the total derivative of the internal energy U in Eq. (A2), and assuming the system as thermodynamically closed, i.e. the molar amount N remains conserved (dN=0), leads to

$$\begin{array}{}\text{(A3)}& \mathrm{d}Q-p\mathrm{d}V={\left.{\displaystyle \frac{\partial U}{\partial T}}\right|}_{V,N}\mathrm{d}T+{\left.{\displaystyle \frac{\partial U}{\partial V}}\right|}_{T,N}\mathrm{d}V\end{array}$$

and subsequently

$$\begin{array}{}\text{(A4)}& \mathrm{d}Q={\left.{\displaystyle \frac{\partial U}{\partial T}}\right|}_{V,N}\mathrm{d}T+\left(p+{\left.{\displaystyle \frac{\partial U}{\partial V}}\right|}_{T,N}\right)\mathrm{d}V.\end{array}$$

If the system's volume is held constant, Eq. (A4) represents the definition of the constant-volume heat capacity C_V in molar units, i.e.

$$\begin{array}{}\text{(A5)}& \mathrm{d}Q={\left.{\displaystyle \frac{\partial U}{\partial T}}\right|}_{V,N}\mathrm{d}T=C_V(p,T)\mathrm{d}T.\end{array}$$

Alternatively, assuming the system's pressure as constant, its volume is variable with total derivative

$$\begin{array}{}\text{(A6)}& \mathrm{d}V={\left.{\displaystyle \frac{\partial V}{\partial T}}\right|}_{p,N}\mathrm{d}T+{\left.{\displaystyle \frac{\partial V}{\partial p}}\right|}_{T,N}\underset{=\mathrm{0}}{\underbrace{\mathrm{d}p}}={\left.{\displaystyle \frac{\partial V}{\partial T}}\right|}_{p,N}\mathrm{d}T\end{array}$$

and therefore

$$\begin{array}{}\text{(A7)}& \begin{array}{rl}\mathrm{d}Q& ={\left.{\displaystyle \frac{\partial U}{\partial T}}\right|}_{V,N}\mathrm{d}T+\left(p+{\left.{\displaystyle \frac{\partial U}{\partial V}}\right|}_{T,N}\right)\mathrm{d}V\\ & ={\left.{\displaystyle \frac{\partial U}{\partial T}}\right|}_{V,N}\mathrm{d}T+\left(p+{\left.{\displaystyle \frac{\partial U}{\partial V}}\right|}_{T,N}\right)\left({\left.{\displaystyle \frac{\partial V}{\partial T}}\right|}_{p,N}\mathrm{d}T\right)\\ & =\left[{\left.{\displaystyle \frac{\partial U}{\partial T}}\right|}_{V,N}+\left(p+{\left.{\displaystyle \frac{\partial U}{\partial V}}\right|}_{T,N}\right){\left.{\displaystyle \frac{\partial V}{\partial T}}\right|}_{p,N}\right]\mathrm{d}T\\ & =C_p(p,T)\mathrm{d}T,\end{array}\end{array}$$

defining the isobaric molar heat capacity C_p. In general, this quantity depends on pressure as well as on temperature. However, if the gas is assumed as ideal, an important conclusion from the statistical description of an ideal gas is the fact that the internal energy U must be independent of the pressure (see e.g. Fay, 1965).

Using this result, together with Eq. (A7) and the ideal-gas law (A1), it follows

$$\begin{array}{}\text{(A8)}& \begin{array}{rl}{C}_p& ={\left.{\displaystyle \frac{\partial U}{\partial T}}\right|}_{V,N}+\left(p+{\left.{\displaystyle \frac{\partial U}{\partial V}}\right|}_{T,N}\right){\left.{\displaystyle \frac{\partial V}{\partial T}}\right|}_{p,N}\\ & ={\left.{\displaystyle \frac{\partial U}{\partial T}}\right|}_{V,N}+p{\left.{\displaystyle \frac{\partial V}{\partial T}}\right|}_{p,N}\\ & ={\left.{\displaystyle \frac{\partial U}{\partial T}}\right|}_{V,N}+{\left.{\displaystyle \frac{\partial }{\partial T}}\left(pV\right)\right|}_{p,N}\\ & ={\left.{\displaystyle \frac{\partial U}{\partial T}}\right|}_{V,N}+{\left.{\displaystyle \frac{\partial }{\partial T}}\left(NRT\right)\right|}_{p,N}\\ & ={\left.{\displaystyle \frac{\partial U}{\partial T}}\right|}_{V,N}+NR.\end{array}\end{array}$$

In the previous computations, there is no restriction on the temperature dependence of the internal energy U(T). Therefore, even by assuming ideal-gas behaviour, the specific heat capacity C_p in Eq. (A8) is in general a function of temperature.

Appendix B: Sensitivity of the conventional definition of θ to perturbations of c_p

This section explores, from a mathematical perspective, the sensitivity of the potential temperature formulation (13) based on a constant specific heat capacity. Considering the specific heat capacity c_p as a variable, the sensitivity of ${\mathit{\theta }}_{c_p}$ (Eq. 13) to a small perturbation δ of c_p is described by its Taylor expansion

$$\begin{array}{}\text{(B1)}& \begin{array}{rl}{\mathit{\theta }}_{c_p+\mathit{\delta }}& ={\mathit{\theta }}_{c_p}+{\displaystyle \frac{\partial {\mathit{\theta }}_{c_p}}{\partial c_p}}\mathit{\delta }+\mathcal{O}\left({\mathit{\delta }}^{\mathrm{2}}\right)\\ & ={\mathit{\theta }}_{c_p}-{\mathit{\theta }}_{c_p}{\displaystyle \frac{R_{\mathrm{a}}}{c_p^{\mathrm{2}}}}\ln \left({\displaystyle \frac{p_{\mathrm{0}}}{p}}\right)\mathit{\delta }+\mathcal{O}\left({\mathit{\delta }}^{\mathrm{2}}\right).\end{array}\end{array}$$

For any constant value of the specific heat capacity c_p and for a minor perturbation δ, the second summand within the expansion (B1) remains small for small values of $\ln \left(\frac{p_{\mathrm{0}}}{p}\right)$. If the interval between the two pressure levels is very narrow, i.e. p≈p₀, the expression $\ln \left(\frac{p_{\mathrm{0}}}{p}\right)$ approximately equals ln (1)=0. Contrarily, if the pressure approaches very low values, i.e. p→0 Pa, the logarithmic expression diverges to negative infinity, i.e. $\ln \left(\frac{p_{\mathrm{0}}}{p}\right)\to -\mathrm{\infty }$, implying that the impact of the second summand intensifies with decreasing pressure, i.e. for increasing altitudes. Moreover, this may explain why the deviation between θ₁₀₀₀ and θ₁₀₁₀, as illustrated in Fig. 2b, remains comparatively small within the troposphere and systematically increases with rising altitude, i.e. decreasing pressure levels.

Appendix C: Approximate computation of the reference potential temperature

This section summarises the detailed steps of approximating the function F(x), defined in Eq. (23), by $\widehat{F}\left(x\right)$, defined in Eq. (28) (Sect. C1), as well as the approximations of the solutions of the resulting nonlinear equations by Newton's method (Sect. C2).

C1 Reformulating the function F(x)

Proceeding from the definition of a function h(x),

$$\begin{array}{}\text{(C1)}& h\left(x\right)=\underset{T_{\mathrm{1}}}{\overset{x}{\int }}{\displaystyle \frac{c_p\left(T^{\prime }\right)}{T^{\prime }}}\mathrm{d}{T}^{\prime },\end{array}$$

with T₁=180 K, the function F(x) may be rearranged as

$$\begin{array}{}\text{(C2)}& \begin{array}{rl}F\left(x\right)& =\underset{x}{\overset{T}{\int }}{\displaystyle \frac{c_p\left(T^{\prime }\right)}{T^{\prime }}}\mathrm{d}{T}^{\prime }-R_{\mathrm{a}}\ln \left({\displaystyle \frac{p}{p_{\mathrm{0}}}}\right)\\ & =h\left(T\right)-h\left(x\right)-R_{\mathrm{a}}\ln \left({\displaystyle \frac{p}{p_{\mathrm{0}}}}\right).\end{array}\end{array}$$

The advantage of this reformulation of F(x) is the inclusion of h(x), consisting of an integral with fixed lower bound and a sole variable upper bound. This way, the function h(x) is numerically solvable, and subsequently h(x) can be substituted by an approximation f(x) that is defined as

$$\begin{array}{}\text{(C3)}& f\left(x\right)=b_{\mathrm{0}}+b_{\mathrm{1}}\ln \left(x-b_{\mathrm{2}}\right)+b_{\mathrm{3}}x+b_{\mathrm{4}}{x}^{\mathrm{2}}.\end{array}$$

Notably, if c_p is constant, this function reduces to an exact primitive of the integrand $\frac{c_p}{T^{\prime }}$ with $b_{\mathrm{3}}=b_{\mathrm{4}}=\mathrm{0}$. Moreover, in this case, the resulting root-finding problem 0=F(x) is exactly solvable and finally leads to the known conventional definition (13) of the potential temperature.

As a further step, the function h(x) is numerically approximated, while c_p(T) in Eq. (C1) is replaced by the ideal-gas limit of air's specific heat capacity $c_p^{\mathrm{0}}\left(T\right)$. The integration interval [T₁, x] with $T_{\mathrm{1}}\le x\le \mathrm{2000}\mathrm{K}$ is traversed in steps of at most 0.001 K, while each step of the integration process is carefully approximated by using Simpson's rule.

By solving a least-squares problem, the coefficients in Eq. (C3) for the approximation of h(x) by the function f(x) are estimated as

$$\begin{array}{}\text{(C4)}& \begin{array}{rl}{b}_{\mathrm{0}}& =-\mathrm{4072.2121328563667},\\ b_{\mathrm{1}}& =\mathrm{797.09247926609601},\\ b_{\mathrm{2}}& =\mathrm{29.587047521428016},\\ b_{\mathrm{3}}& =\mathrm{0.41981158226925142},\\ b_{\mathrm{4}}& =-\mathrm{5.1008025097060311}\times {\mathrm{10}}^{-\mathrm{5}}.\end{array}\end{array}$$

In Fig. C1a the function h(x) is graphed together with the approximation f(x), as well as the respective deviations h(x)−f(x) in Fig. C1b. Evidently, the absolute error inherent to the approximations is comparatively small as, over the entire temperature range above 190 K, the approximation error never exceeds ±1 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$. Exclusively at temperatures below 190 K, the approximation error rapidly rises above 1 $\mathrm{J}{\mathrm{kg}}^{-\mathrm{1}}{\mathrm{K}}^{-\mathrm{1}}$, bearing in mind that such absolute temperatures are only occasionally found in the atmosphere within a relatively narrow altitude interval at the cold point tropopause. Moreover, the difference between f(x) and h(x) appears negligible as the profiles almost ideally coincide (cf. Fig. C1a).

Figure C1 (a) Numerically evaluated function h(x) together with its approximation f(x); (b) the absolute approximation error h(x)−f(x).

Table C1 Values of the new reference potential temperature θ_ref, together with the first two iterates θ⁽¹⁾, θ⁽²⁾ using Newton's method and the first iterate ${\mathit{\theta }}_{\mathrm{Householder}}^{\left(\mathrm{1}\right)}$ using Householder's method for five pairs of temperature and pressure along the US Standard Atmosphere. The computed values are rounded to three digits.

C2 Finalised approximation of the reference potential temperature

As discussed in Sect. 5.1, the new formulation of the potential temperature based on the temperature-dependent specific heat capacity c_p(T) requires solving the root-finding problem 0=F(x), where the function F(x) is defined in Eq. (23). However, since F(x) contains an integral that complicates the root-finding process, this integral is substituted by the difference f(T)−f(x), where f is given in Sect. C1. Therefore, F(x) is replaced by the function $\widehat{F}\left(x\right)$ as defined in Eq. (28), and the zero of the equation $\mathrm{0}=\widehat{F}\left(x\right)$ is denoted as ${\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}$.

The equation $\mathrm{0}=\widehat{F}\left(x\right)$ is still not analytically solvable, so Newton's method is once more required. Using again x₀=θ₁₀₀₅ as the initial guess, cf. Eq. (26), the iteration sequence for Newton's method is given by the recursion

$$\begin{array}{}\text{(C5)}& \begin{array}{rl}{x}_{k+\mathrm{1}}& =x_k-{\displaystyle \frac{\widehat{F}\left(x_k\right)}{{\widehat{F}}^{\prime }\left(x_k\right)}}\\ & =x_k-{\displaystyle \frac{f\left(T\right)-f\left(x_k\right)-R_{\mathrm{a}}\ln \left(\frac{p}{p_{\mathrm{0}}}\right)}{-f^{\prime }\left(x_k\right)}}\\ & =x_k-{\displaystyle \frac{R_{\mathrm{a}}\ln \left(\frac{p}{p_{\mathrm{0}}}\right)-f\left(T\right)+f\left(x_k\right)}{f^{\prime }\left(x_k\right)}}.\end{array}\end{array}$$

Instead of this standard formulation of Newton's method (C5), Householder's formulation

$$\begin{array}{}\text{(C6)}& \begin{array}{rl}{x}_{k+\mathrm{1}}& =x_k-{\displaystyle \frac{\widehat{F}\left(x_k\right)}{{\widehat{F}}^{\prime }\left(x_k\right)}}-{\displaystyle \frac{{\widehat{F}}^{\prime \prime }\left(x_k\right)}{\mathrm{2}{\widehat{F}}^{\prime }\left(x_k\right)}}{\left[{\displaystyle \frac{\widehat{F}\left(x_k\right)}{{\widehat{F}}^{\prime }\left(x_k\right)}}\right]}^{\mathrm{2}}\\ & =x_k-{\displaystyle \frac{R_{\mathrm{a}}\ln \left(\frac{p}{p_{\mathrm{0}}}\right)-f\left(T\right)+f\left(x_k\right)}{f^{\prime }\left(x_k\right)}}\\ & -{\displaystyle \frac{f^{\prime \prime }\left(x_k\right)}{\mathrm{2}{f}^{\prime }\left(x_k\right)}}{\left[{\displaystyle \frac{R_{\mathrm{a}}\ln \left(\frac{p}{p_{\mathrm{0}}}\right)-f\left(T\right)+f\left(x_k\right)}{f^{\prime }\left(x_k\right)}}\right]}^{\mathrm{2}}\end{array}\end{array}$$

may be used, which allows for reducing the computation time due to its accelerated convergence speed. For completeness, the required derivatives f^′ and f^′′ in the recursion formulas (C5) and (C6) are

$$\begin{array}{}\text{(C7)}& \begin{array}{rl}{f}^{\prime }\left(x\right)& ={\displaystyle \frac{b_{\mathrm{1}}}{x-b_{\mathrm{2}}}}+b_{\mathrm{3}}+\mathrm{2}{b}_{\mathrm{4}}x,\\ f^{\prime \prime }\left(x\right)& =\mathrm{2}{b}_{\mathrm{4}}-{\displaystyle \frac{b_{\mathrm{1}}}{{\left(x-b_{\mathrm{2}}\right)}^{\mathrm{2}}}}.\end{array}\end{array}$$

The final step on the way to formulate a new expression for the potential temperature requires defining one of the iterates x_k as appropriate enough for the approximations that result from applying the different methods:

• the standard of Newton's method (C5), simply referred to as Newton's method in the sequel, or

• Householder's method (C6).

While the mathematical expressions in Eqs. (C5) and (C6) are of increasing complexity, the convergence rate of the approximating sequence increases with rising mathematical complication. The preferred method is determined by the accuracy required; i.e. better accuracy is necessarily associated with elevated computational effort for the approximation method. A discussion of the approximation errors is found in Appendix D.

Table C1 collects values of the new reference potential temperature θ_ref, together with the first two iterates θ⁽¹⁾, θ⁽²⁾ using Newton's method (C5) and the first iterate ${\mathit{\theta }}_{\mathrm{Householder}}^{\left(\mathrm{1}\right)}$ using Householder's method (C6) for five pairs of temperature and pressure along the US Standard Atmosphere, cf. Fig. 1, which allows verification of computations. The first height is chosen midway along the linearly decreasing temperature profile within the troposphere, while the other heights correspond to the kinks of the temperature profile.

Appendix D: Approximation error for the reference potential temperature

The following aims at a comprehensive investigation of the errors inherent with approximating the ultimate reference potential temperature θ_ref. As discussed in Sect. 5.2, the total error is a combination of the basic error ${\mathit{\theta }}_{\mathrm{ref}}-{\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}$ and the approximation error that results from the approximation sequence ${\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}-{\mathit{\theta }}^{\left(k\right)}$, where θ^(k) denotes the kth iterate of the approximation sequence which is computed in accordance with either Newton's or Householder's method. The formulations of Newton's (Eq. C5) and Householder's (Eq. C6) method require replacing the function F(x) by $\widehat{F}\left(x\right)$, and the approximation sequences θ^(k) converge to ${\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}$ for k→∞. Consequently, the approximation error ${\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}-{\mathit{\theta }}^{\left(k\right)}$ tends to zero for k→∞.

Figure D1 Total relative error along the US Standard Atmosphere arising from the iteration process by declaring (a) the first iterate θ⁽¹⁾ or (b) the second iterate θ⁽²⁾ as the final approximation to the reference potential temperature θ_ref. Red curves: iterates computed using Newton's method (C5); blue curves: iterates computed using Householder's method (C6). Note the different range of the abscissae.

The analysis of the approximation error is initially based on the pressure and temperature profiles of the US Standard Atmosphere. Figure D1 shows the total relative errors $\left({\mathit{\theta }}_{\mathrm{ref}}-{\mathit{\theta }}^{\left(\mathrm{1}\right)}\right)/{\mathit{\theta }}_{\mathrm{ref}}$ of the first iterate (Fig. D1a) and $\left({\mathit{\theta }}_{\mathrm{ref}}-{\mathit{\theta }}^{\left(\mathrm{2}\right)}\right)/{\mathit{\theta }}_{\mathrm{ref}}$ of the second iterate (Fig. D1b), computed with Newton's or Householder's method. The first iterate still causes the approximation to have significant errors, especially at altitudes above 35 km. However, the second iterate with either Newton's or Householder's method yields results with negligible approximation error. Hence, the total error of the approximation procedure is dominated by the unavoidable basic error and may be deduced from the provided figures whenever the total error profile nearly congruently follows the profile of the basic error (cf. Figs. D1b and 5a).

It may be noted that Householder's method achieves a significantly lower error level than Newton's method due to its accelerated rate of convergence. Compared to the first iterate approximations, computation up to the second iterate (cf. Fig. D1b) achieves, in general, a considerable improvement for both methods, and both second iterate approximations approach the basic error quite closely (cf. Fig. D1b). As is also evident from Fig. D1b, compared to Householder's method, the second iterate with Newton's method results in a smaller total relative error $\left({\mathit{\theta }}_{\mathrm{ref}}-{\mathit{\theta }}^{\left(\mathrm{2}\right)}\right)/{\mathit{\theta }}_{\mathrm{ref}}$ relative to the ultimate reference potential temperature (indicated by a smaller distance to the dashed zero-line above 45 km altitude). Nevertheless, the relative approximation error, $\left({\mathit{\theta }}_{\mathrm{ref}}^{\mathrm{approx}}-{\mathit{\theta }}^{\left(\mathrm{2}\right)}\right)/{\mathit{\theta }}_{\mathrm{ref}}$, is larger compared to the second iterate with Householder's method. So, luckily, the second iterate with Newton's method provides a better approach to the reference potential temperature than that with Householder's method.

Figure D2 Relative error of the second iterates θ⁽²⁾ with (a) Newton's method and (b) Householder's method for the ranges of pressure and temperature from 1000 to 0.5 hPa and from 180 to 300 K, respectively. (c) The absolute error arising from the first iterate θ⁽¹⁾ with Householder's method. The white solid line indicates the p–T profile from the US Standard Atmosphere. Note the different ranges of the Δθ scales.

As with the discussion of the basic error in Sect. 5.2, the analysis of the total error should include all possible combinations of pressure and temperature in order to take into account fluctuations in the real atmosphere that deviate from the profile of the US Standard Atmosphere. Therefore, the extended analysis of the approximation error is summarised in Fig. D2. The upper panels illustrate the total relative error of the second iterate for Newton's (Fig. D2a) and Householder's method (Fig. D2b). As previously shown, further iteration with either method does not improve the approximation quality. The contour patterns in these panels show a remarkable similarity to the contours for the relative error of the basic approximation in Fig. 5c. Also here (upper panels of Fig. D2), two regions are highlighted by the contours, i.e. at ∼100 hPa and in a pressure range from ∼5 to 1 hPa, featuring the same impact on Δθ∕θ_ref of identical strength as the basic error. This result may not be surprising, since the second iteration step with both methods, Newton's and Householder's, was already proven to approach the approximation comparatively well, without worsening the total error level (cf. Fig. D1b).

Consequently, concerning the required number of iterations and the method to use, the second iteration of Newton's method can be recommended to deliver appropriate results, with a relative error of less than 0.3 %, up to the stratopause level (∼50 km). Householder's method features an accelerated convergence rate, and its use up to its first iterate θ⁽¹⁾ may be already appropriate for certain applications. According to the total error of Householder's method up to its first iterate θ⁽¹⁾ (Fig. D2c), the resulting relative error remains below 7 % to a pressure level of ∼50 hPa and Δθ stays below 0.3 % to pressures of ∼2 hPa. Thus, Fig. D2 may serve as guidance to decide how many iterations with one or the other method best meet the individual accuracy requirements.

Appendix E: The derivative of the reference potential temperature

As discussed in Sect. 5.1, the new reference potential temperature is defined as the zero of the function

$$\begin{array}{}\text{(E1)}& F(x,p,T)=\underset{x}{\overset{T}{\int }}{\displaystyle \frac{c_p\left(T^{\prime }\right)}{T^{\prime }}}\mathrm{d}{T}^{\prime }-R_{\mathrm{a}}\ln \left({\displaystyle \frac{p}{p_{\mathrm{0}}}}\right)\end{array}$$

for given values of pressure p and temperature T; see Eq. (23). More precisely, for varying p, T, a function $(p,T)\mapsto {\mathit{\theta }}_{\mathrm{ref}}(p,T)$ is implicitly defined by the equation

$$\begin{array}{}\text{(E2)}& F\left({\mathit{\theta }}_{\mathrm{ref}}\right(p,T),p,T)=\mathrm{0}.\end{array}$$

According to the implicit function theorem (e.g. Protter and Morrey, 1985, chap. 7), Eq. (E2) is uniquely solvable for θ_ref(p, T); i.e. the function $(p,T)\mapsto {\mathit{\theta }}_{\mathrm{ref}}(p,T)$ actually exists as a differentiable function of (p, T), if the condition $\frac{\partial F}{\partial \mathit{\theta }}\ne \mathrm{0}$ holds. According to Eq. (E1), this partial derivative equals

$$\begin{array}{}\text{(E3)}& {\displaystyle \frac{\partial F}{\partial \mathit{\theta }}}({\mathit{\theta }}_{\mathrm{ref}},p,T)=-{\displaystyle \frac{c_p\left({\mathit{\theta }}_{\mathrm{ref}}\right)}{{\mathit{\theta }}_{\mathrm{ref}}}},\end{array}$$

being strictly negative, since the specific heat capacity is always positive. Moreover, the implicit function theorem states that the derivatives of the implicit function θ_ref(p, T) are given by

$$\begin{array}{}\text{(E4)}& \begin{array}{rl}& \left[{\displaystyle \frac{\partial {\mathit{\theta }}_{\mathrm{ref}}}{\partial p}}(p,T),{\displaystyle \frac{\partial {\mathit{\theta }}_{\mathrm{ref}}}{\partial T}}(p,T)\right]\\ & =-{\left({\displaystyle \frac{\partial F}{\partial \mathit{\theta }}}({\mathit{\theta }}_{\mathrm{ref}},p,T)\right)}^{-\mathrm{1}}\\ & \cdot \left[{\displaystyle \frac{\partial F}{\partial p}}({\mathit{\theta }}_{\mathrm{ref}},p,T),{\displaystyle \frac{\partial F}{\partial T}}({\mathit{\theta }}_{\mathrm{ref}},p,T)\right]\\ & ={\displaystyle \frac{{\mathit{\theta }}_{\mathrm{ref}}}{c_p\left({\mathit{\theta }}_{\mathrm{ref}}\right)}}\left[-{\displaystyle \frac{R_{\mathrm{a}}}{p}},{\displaystyle \frac{c_p\left(T\right)}{T}}\right]\\ & =\left[-{\displaystyle \frac{R_{\mathrm{a}}}{c_p\left({\mathit{\theta }}_{\mathrm{ref}}\right)}}{\displaystyle \frac{{\mathit{\theta }}_{\mathrm{ref}}}{p}},{\displaystyle \frac{{\mathit{\theta }}_{\mathrm{ref}}}{T}}{\displaystyle \frac{c_p\left(T\right)}{c_p\left({\mathit{\theta }}_{\mathrm{ref}}\right)}}\right].\end{array}\end{array}$$

Note that these partial derivatives coincide with the partial derivatives of ${\mathit{\theta }}_{c_p}$ in the case of a constant specific heat capacity. Using the partial derivatives (Eq. E4), the total differential of θ_ref may be written as

$$\begin{array}{}\text{(E5)}& \begin{array}{rl}\mathrm{d}{\mathit{\theta }}_{\mathrm{ref}}& ={\displaystyle \frac{\partial {\mathit{\theta }}_{\mathrm{ref}}}{\partial p}}\mathrm{d}p+{\displaystyle \frac{\partial {\mathit{\theta }}_{\mathrm{ref}}}{\partial T}}\mathrm{d}T\\ & =-{\displaystyle \frac{R_{\mathrm{a}}}{c_p\left({\mathit{\theta }}_{\mathrm{ref}}\right)}}{\displaystyle \frac{{\mathit{\theta }}_{\mathrm{ref}}}{p}}\mathrm{d}p+{\displaystyle \frac{{\mathit{\theta }}_{\mathrm{ref}}}{T}}{\displaystyle \frac{c_p\left(T\right)}{c_p\left({\mathit{\theta }}_{\mathrm{ref}}\right)}}\mathrm{d}T\end{array}\end{array}$$

or

$$\begin{array}{}\text{(E6)}& c_p\left({\mathit{\theta }}_{\mathrm{ref}}\right){\displaystyle \frac{\mathrm{d}{\mathit{\theta }}_{\mathrm{ref}}}{{\mathit{\theta }}_{\mathrm{ref}}}}=c_p\left(T\right){\displaystyle \frac{\mathrm{d}T}{T}}-R_{\mathrm{a}}{\displaystyle \frac{\mathrm{d}p}{p}}.\end{array}$$