Anomaly Format of Atmospheric Governing Equations with Climate as a Reference Atmosphere

Abstract

To reduce numerical instability and increase forecast accuracy of a numerical weather prediction (NWP) model, one approach is to subtract a reference atmosphere from atmospheric governing equations. In the past, scientists have proposed one-dimensional, two-dimensional, and three-dimensional static (in time) reference atmospheres with respect to temperature and pressure. These three reference atmospheres were first reviewed, and their corresponding perturbation equations were derived. Then, a new four-dimensional (space and time) all-variable (temperature, pressure, wind, moisture, etc.) reference atmosphere was defined using observed climatic states. Unlike the previous three approaches, the perturbations derived from this new method are actual anomalies relative to climate and directly a part of individual weather systems in both structure and strength. By subtracting climatic states, anomaly equations were derived and analyzed. Finally, the benefits and challenges of the anomaly-equation-based NWP model were discussed. Theoretically, an anomaly model should reduce model systematic errors (bias) and should avoid model climate drift to significantly enhance a model’s performance. An example of tropical cyclone track forecasts using the Beta advection model (vorticity) was demonstrated. The separation of model physics into climatic and anomalous physics is a significant challenge if a pure anomaly-equation-based NWP model is desired. Fortunately, a model including both anomaly and climatic equations should work with current full physics. In an anomaly climate mixed model, the anomaly part needs to be predicted and the climate parts are precalculated constants. It is hoped that this study will inspire model developers to explore the approach, which could be a possible new direction in developing next-generation NWP models. A high-resolution reanalysis is also key to the success of this new approach.

1. Introduction

Atmospheric governing equations are used in numerical weather prediction (NWP) models. Because nonlinear differential equations are difficult to solve analytically, they are solved numerically by writing them into a discrete format following a chosen numerical scheme. Consequently, during equation integration in discrete format, small numerical errors can accumulate with time and lead to large errors in forecasts. To improve accuracy in solving these equations, many approximations have been made under various situations [1,2,3,4,5,6,7,8,9,10,11,12,13,14]. One main issue of solving the atmospheric governing equations is the “small difference between two large-value terms” [15]. Since weather change is driven by those small differences or disturbances (departing from large-value balances such as hydrostatic and geostrophic balances), relatively small differences in those large-value terms often result in big errors in final forecasts. Therefore, it is desirable that large-value background terms (satisfying dominant atmospheric balances) be avoided or eliminated during model calculation to improve forecast accuracy as well as to reduce numerical instability in the calculations.

Along this line of thinking, researchers have proposed reference atmosphere to decompose atmospheric governing equations into perturbation equations [4,16,17,18,19,20,21,22,23,24]. Some of these approaches have already been applied to both grid and spectral general circulation models [20,22,25,26,27,28]. During the last half century, the reference atmosphere concept has been expanded from one (1D) to two (2D), and then three (3D) dimensional spaces. However, these past-proposed reference atmospheres are only with respect to a few variables (mainly temperature and pressure through hydrostatic relationship) and static (with no time evolution). Their perturbations also have no direct relation to weather systems. In this paper, we propose a new four-dimensional (4D) reference atmosphere for all atmospheric variables. The new reference atmosphere changes not only in three-dimensional space but also in time. More importantly, perturbations are directly a part of weather systems both in structure and strength. The advantages and challenges of this new approach are discussed, which provide model developers with possibly a new way of thinking in developing next-generation NWP models.

In this paper, in Section 2, first the previously proposed 1D, 2D, and 3D reference atmospheres are outlined and derive their corresponding perturbation equations derived. In Section 3, the new 4D reference atmosphere is proposed, related perturbation equations are derived, and its application to NWP model is discussed. A summary is given in Section 4.

2. Perturbation Equations under 1D, 2D, and 3D Reference Atmospheres

2.1. Atmospheric Governing Equations

Assuming that the earth is a perfect sphere with a radius of a, the hydrostatic (dry) atmospheric governing equations can be written, in latitude-longitude coordinates, as follows (which can be seen in dynamic meteorology textbook or numerical weather prediction articles such as [29,30,31]:

$$\frac{du}{dt}=-\frac{C_p\theta }{a cos \phi }\frac{\partial \pi }{\partial \lambda }+fv+F_u$$

(1)

$$\frac{dv}{dt}=-\frac{C_p\theta }{a}\frac{\partial \pi }{\partial \phi }-fu+F_v$$

(2)

$$\frac{dw}{dt}=-C_p\theta \frac{\partial \pi }{\partial z}-g+F_w$$

(3)

$$\frac{d\pi }{dt}=-\frac{1}{\left(\gamma -1\right)}\pi ·D_3+\frac{1}{\left(\gamma -1\right)}\frac{F_{\theta }}{\theta }$$

(4)

$$\frac{d\theta }{dt}=\frac{F_{\theta }}{\pi }$$

(5)

$$T=\pi \theta$$

(6)

where $\lambda$ is longitude ($\overrightarrow{i}$ direction); $\phi$ is latitude ($\overrightarrow{j}$ direction); $r=z+a$, z is vertical height above the earth surface ($\overrightarrow{k}$ direction); t is time; $\overrightarrow{V}=u\overrightarrow{i}+v\overrightarrow{j}+w\overrightarrow{k}$ three-dimensional vector wind; $D_3=\nabla ·\overrightarrow{V}=\frac{1}{r cos \phi }\frac{\partial u}{\partial \lambda }+\frac{1}{r cos \phi }\frac{\partial v cos \phi }{\partial \phi }+\frac{1}{r^2}\frac{\partial r^{2}w}{\partial r}$ three-dimensional wind convergence; T is temperature; $\theta$ potential temperature, $\pi ={\left(\frac{p}{p_0}\right)}^{\frac{R}{C_p}}$; p is air pressure at height z; $p_0$ air pressure at z = 0; $\gamma =C_p/R$; heat capacity $C_p=1004 {\mathrm{JK}}^{-1}{\mathrm{kg}}^{-1}$; gas constant $R=8.314 {\mathrm{JK}}^{-1}{\mathrm{mL}}^{-1}$; $f=2\Omega sin \phi$ the Coriolis force; $\Omega =7.292\times {10}^{-5}\mathrm{rad}/\mathrm{s}$ the angular velocity of the earth rotation; $F_u$, $F_v$, and $F_w$ are frictions; $F_{\theta }$ is heating and cooling; and g = 9.8 ${\mathrm{ms}}^{-2}$ the earth gravity at earth surface.

Equations (1)–(3) describe air motion (wind) in three-dimensional space, Equation (4) describes the air pressure (mass) change, Equation (5) describes thermal state (temperature change), and Equation (6) is the relationship between temperature and potential temperature. Note that the above equations are for dry atmosphere only, since moisture is not considered.

Because there are five independent governing equations (Equations (1)–(5)), the five variables (u, v, w, $\pi$, and $\theta$ or T) can be solved through these prognostic equations to predict basic atmospheric states at a future time t. The perturbation format of these equations will be derived below using 1D, 2D, and 3D reference atmospheres, respectively.

2.2. One-Dimensional Static Temperature as a Reference Atmosphere

One main characteristic of the earth atmosphere is the vertical layered structure of temperature. Therefore, Chao [16,17], Zeng [19,32], and Zeng et al. [20] proposed a reference atmosphere with a one-dimensional static temperature, where temperature varies only with height $\overline{T}\left(z\right)$, and reference pressure $\overline{\pi }\left(z\right)$ satisfies the hydrostatic balance:

$$C_p\overline{\theta }\left(z\right)\frac{\partial \overline{\pi }\left(z\right)}{\partial z}=-g$$

(7)

The reference potential temperature is:

$$\overline{\theta }\left(z\right)=\overline{T}\left(z\right)/\overline{\pi }\left(z\right)$$

(8)

$\overline{T}\left(z\right)$, $\overline{\pi }\left(z\right)$ and $\overline{\theta }\left(z\right)$ are globally averaged values and could be obtained from historical observation data or current observation.

Under this assumption, the perturbed variables are:

$$T^{\prime }\left(\lambda ,\phi ,z,t\right)=T\left(\lambda ,\phi ,z,t\right)-\overline{T}\left(z\right)$$

(9)

$${\pi }^{\prime }\left(\lambda ,\phi ,z,t\right)=\pi \left(\lambda ,\phi ,z,t\right)-\overline{\pi }\left(z\right)$$

(10)

$${\theta }^{\prime }\left(\lambda ,\phi ,z,t\right)=\theta \left(\lambda ,\phi ,z,t\right)-\overline{\theta }\left(z\right)$$

(11)

Using Equations (7)–(11), the atmospheric governing equations (Equations (1)–(6)) become:

$$\frac{du}{dt}=-\frac{C_p{\theta }^{\prime }}{a cos \phi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }-\frac{C_p\overline{\theta }}{a cos \phi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }+fv+F_u$$

(12)

$$\frac{dv}{dt}=-\frac{C_p{\theta }^{\prime }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }-\frac{C_p\overline{\theta }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }-fu+F_v$$

(13)

$$\frac{dw}{dt}=-C_p{\theta }^{\prime }\frac{\partial {\pi }^{\prime }}{\partial z}-C_p\overline{\theta }\frac{\partial {\pi }^{\prime }}{\partial z}-C_p{\theta }^{\prime }\frac{\partial \overline{\pi }}{\partial z}+F_w$$

(14)

$$\frac{d{\pi }^{\prime }}{dt}=-\frac{1}{\left(\gamma -1\right)}{\pi }^{\prime }·D_3-\frac{1}{\left(\gamma -1\right)}\overline{\pi }·D_3+\frac{1}{\left(\gamma -1\right)}\frac{F_{\overline{\theta }}+F_{{\theta }^{\prime }}}{\overline{\theta }+{\theta }^{\prime }}+w\frac{g}{C_p\overline{\theta }}$$

(15)

$$\frac{d{\theta }^{\prime }}{dt}=\frac{F_{\overline{\theta }}+F_{{\theta }^{\prime }}}{\overline{\pi }+{\pi }^{\prime }}-w\frac{\partial \overline{\theta }}{\partial z}$$

(16)

$$T^{\prime }={\pi }^{\prime }{\theta }^{\prime }+\overline{\pi }{\theta }^{\prime }+{\pi }^{\prime }\overline{\theta }$$

(17)

where wind (u, v, w) is full wind with no decomposition. Equations (12)–(17) are the perturbation form of the governing equations after eliminating the one-dimensional static-temperature reference atmosphere. Note that the physical process term $F_{\theta }$ is also decomposed into reference and perturbed terms, i.e., $F_{\theta }=F_{\overline{\theta }}+F_{{\theta }^{\prime }}$.

If we consider $\overline{\pi }\gg {\pi }^{\prime }$ and $\overline{\theta }\gg {\theta }^{\prime }$, Equations (12)–(17) could be further simplified into:

$$\frac{du}{dt}=-\frac{C_p\overline{\theta }}{a cos \phi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }+fv+F_u$$

(18)

$$\frac{dv}{dt}=-\frac{C_p\overline{\theta }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }-fu+F_v$$

(19)

$$\frac{dw}{dt}=-C_p\overline{\theta }\frac{\partial {\pi }^{\prime }}{\partial z}-C_p{\theta }^{\prime }\frac{\partial \overline{\pi }}{\partial z}+F_w$$

(20)

$$\frac{d{\pi }^{\prime }}{dt}=-\frac{1}{\left(\gamma -1\right)}\overline{\pi }·D_3+\frac{1}{\left(\gamma -1\right)}\frac{F_{\overline{\theta }}+F_{{\theta }^{\prime }}}{\overline{\theta }}+w\frac{g}{C_p\overline{\theta }}$$

(21)

$$\frac{d{\theta }^{\prime }}{dt}=\frac{F_{\overline{\theta }}+F_{{\theta }^{\prime }}}{\overline{\pi }}-w\frac{\partial \overline{\theta }}{\partial z}$$

(22)

$$T^{\prime }=\overline{\pi }{\theta }^{\prime }+{\pi }^{\prime }\overline{\theta }$$

(23)

All terms on the right side of Equations (18)–(22) are linear. It has been reported that this decomposition approach of 1D reference atmosphere was used by the ECMWF spectral model and that it improved its medium-range forecast accuracy [22,33].

2.3. Two-Dimensional Static Temperature as a Reference Atmosphere

Considering the prominent north-south temperature gradient in real atmosphere, Zhang et al. [34] extended this static-temperature reference atmosphere from one-dimensional to two-dimensional in space, i.e., temperature is not only a function of height but also latitude $\overline{T}\left(\phi ,z\right)$. So, reference pressure through the hydrostatic balance is:

$$C_p\overline{\theta }\left(\phi ,z\right)\frac{\partial \overline{\pi }\left(\phi ,z\right)}{\partial z}=-g$$

(24)

and the reference potential temperature is:

$$\overline{\theta }\left(\phi ,z\right)=\overline{T}\left(\phi ,z\right)/\overline{\pi }\left(\phi ,z\right)$$

(25)

$\overline{T}\left(\phi ,z\right)$, $\overline{\pi }\left(\phi ,z\right)$ and $\overline{\theta }\left(\phi ,z\right)$ are zonally averaged values and could be obtained from historical observation data or current observation.

Similarly, to obtain the perturbation equations, the temperature, pressure, and potential temperature can be decomposed as follows:

$$T^{\prime }\left(\lambda ,\phi ,z,t\right)=T\left(\lambda ,\phi ,z,t\right)-\overline{T}\left(\phi ,z\right)$$

(26)

$${\pi }^{\prime }\left(\lambda ,\phi ,z,t\right)=\pi \left(\lambda ,\phi ,z,t\right)-\overline{\pi }\left(\phi ,z\right)$$

(27)

$${\theta }^{\prime }\left(\lambda ,\phi ,z,t\right)=\theta \left(\lambda ,\phi ,z,t\right)-\overline{\theta }\left(\phi ,z\right)$$

(28)

Using Equations (24)–(28), the atmospheric governing equations (Equations (1)–(6)) become:

$$\frac{du}{dt}=-\frac{C_p{\theta }^{\prime }}{a cos \phi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }-\frac{C_p\overline{\theta }}{a cos \phi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }+fv+F_u$$

(29)

$$\frac{dv}{dt}=-\frac{C_p{\theta }^{\prime }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }-\frac{C_p\overline{\theta }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }-\frac{C_p\overline{\theta }}{a}\frac{\partial \overline{\pi }}{\partial \phi }-\frac{C_p{\theta }^{\prime }}{a}\frac{\partial \overline{\pi }}{\partial \phi }-fu+F_v$$

(30)

$$\frac{dw}{dt}=-C_p{\theta }^{\prime }\frac{\partial {\pi }^{\prime }}{\partial z}-C_p\overline{\theta }\frac{\partial {\pi }^{\prime }}{\partial z}-C_p{\theta }^{\prime }\frac{\partial \overline{\pi }}{\partial z}+F_w$$

(31)

$$\frac{d{\pi }^{\prime }}{dt}=-\frac{1}{\left(\gamma -1\right)}{\pi }^{\prime }\cdot D_3-\frac{1}{\left(\gamma -1\right)}\overline{\pi }\cdot D_3+\frac{1}{\left(\gamma -1\right)}\frac{F_{\overline{\theta }}+F_{{\theta }^{\prime }}}{\overline{\theta }+{\theta }^{\prime }}-v\frac{\partial \overline{\pi }}{a\partial \phi }-w\frac{\partial \overline{\pi }}{\partial z}$$

(32)

$$\frac{d{\theta }^{\prime }}{dt}=+\frac{F_{\overline{\theta }}+F_{{\theta }^{\prime }}}{\overline{\pi }+{\pi }^{\prime }}-v\frac{\partial \overline{\theta }}{a\partial \phi }-w\frac{\partial \overline{\theta }}{\partial z}$$

(33)

$$T^{\prime }={\pi }^{\prime }{\theta }^{\prime }+\overline{\pi }{\theta }^{\prime }+{\pi }^{\prime }\overline{\theta }$$

(34)

where wind (u, v, w) is full wind with no decomposition. Equations (29)–(34) are the perturbation form of the governing equations after eliminating the two-dimensional static-temperature reference atmosphere.

If we consider $\overline{\pi }\gg {\pi }^{\prime }$ and $\overline{\theta }\gg {\theta }^{\prime }$, Equations (29)–(34) could be further simplified into:

$$\frac{du}{dt}=-\frac{C_p\overline{\theta }}{a cos \phi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }+fv+F_u$$

(35)

$$\frac{dv}{dt}=-\frac{C_p\overline{\theta }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }-\frac{C_p\overline{\theta }}{a}\frac{\partial \overline{\pi }}{\partial \phi }-fu+F_v$$

(36)

$$\frac{dw}{dt}=-C_p\overline{\theta }\frac{\partial {\pi }^{\prime }}{\partial z}-C_p{\theta }^{\prime }\frac{\partial \overline{\pi }}{\partial z}+F_w$$

(37)

$$\frac{d{\pi }^{\prime }}{dt}=-\frac{1}{\left(\gamma -1\right)}\overline{\pi }·D_3+\frac{1}{\left(\gamma -1\right)}\frac{F_{\overline{\theta }}+F_{{\theta }^{\prime }}}{\overline{\theta }}-v\frac{\partial \overline{\pi }}{a\partial \phi }-w\frac{\partial \overline{\pi }}{\partial z}$$

(38)

$$\frac{d{\theta }^{\prime }}{dt}=\frac{F_{\overline{\theta }}+F_{{\theta }^{\prime }}}{\overline{\pi }}-v\frac{\partial \overline{\theta }}{a\partial \phi }-w\frac{\partial \overline{\theta }}{\partial z}$$

(39)

$$T^{\prime }=\overline{\pi }{\theta }^{\prime }+{\pi }^{\prime }\overline{\theta }$$

(40)

All terms on the right side of Equations (35)–(39) are also linear. Zhang et al. [34] compared the 1D and 2D approaches and found that the latter could increase numerical stability for a model.

2.4. Three-Dimensional Static Temperature as a Reference Atmosphere

The two-dimensional reference atmosphere was further expanded to three-dimensional atmosphere by Wood et al. [35] and Su et al. [36,37]. Reference temperature $\overline{T}\left(\lambda ,\phi ,z\right)$ is not only a function of height and latitude but also longitude. The hydrostatic balance in the three-dimensional reference atmosphere is:

$$C_p\overline{\theta }\left(\lambda ,\phi ,z\right)\frac{\partial \overline{\pi }\left(\lambda ,\phi ,z\right)}{\partial z}=-g$$

(41)

The potential temperature becomes:

$$\overline{\theta }\left(\lambda ,\phi ,z\right)=\overline{T}\left(\lambda ,\phi ,z\right)/\overline{\pi }\left(\lambda ,\phi ,z\right)$$

(42)

$\overline{T}\left(\lambda ,\phi ,z\right)$, $\overline{\pi }\left(\lambda ,\phi ,z\right)$ and $\overline{\theta }\left(\lambda ,\phi ,z\right)$ could be time-averaged values from historical observation data [36,37] or an observation and even a forecast at a selected time, such as model initial condition or any moment during a model integration [35]. Although $\overline{T}\left(\lambda ,\phi ,z\right)$ is still static (i.e., no time dimension), it is a real atmosphere.

To obtain the perturbation equations, the temperature, pressure, and potential temperature can be decomposed as follows:

$$T^{\prime }\left(\lambda ,\phi ,z,t\right)=T\left(\lambda ,\phi ,z,t\right)-\overline{T}\left(\lambda ,\phi ,z\right)$$

(43)

$${\pi }^{\prime }\left(\lambda ,\phi ,z,t\right)=\pi \left(\lambda ,\phi ,z,t\right)-\overline{\pi }\left(\lambda ,\phi ,z\right)$$

(44)

$${\theta }^{\prime }\left(\lambda ,\phi ,z,t\right)=\theta \left(\lambda ,\phi ,z,t\right)-\overline{\theta }\left(\lambda ,\phi ,z\right)$$

(45)

Using Equations (41)–(45), the atmospheric governing equations (Equations (1)–(6)) become:

$$\frac{du}{dt}=-\frac{C_p{\theta }^{\prime }}{a cos \phi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }-\frac{C_p\overline{\theta }}{a cos \phi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }-\frac{C_p\overline{\theta }}{a cos \phi }\frac{\partial \overline{\pi }}{\partial \lambda }-\frac{C_p{\theta }^{\prime }}{a cos \phi }\frac{\partial \overline{\pi }}{\partial \lambda }+fv+F_u$$

(46)

$$\frac{dv}{dt}=-\frac{C_p{\theta }^{\prime }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }-\frac{C_p\overline{\theta }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }-\frac{C_p\overline{\theta }}{a}\frac{\partial \overline{\pi }}{\partial \phi }-\frac{C_p{\theta }^{\prime }}{a}\frac{\partial \overline{\pi }}{\partial \phi }-fu+F_v$$

(47)

$$\frac{dw}{dt}=-C_p{\theta }^{\prime }\frac{\partial {\pi }^{\prime }}{\partial z}-C_p\overline{\theta }\frac{\partial {\pi }^{\prime }}{\partial z}-C_p{\theta }^{\prime }\frac{\partial \overline{\pi }}{\partial z}+F_w$$

(48)

$$\frac{d{\pi }^{\prime }}{dt}=-\frac{1}{\left(\gamma -1\right)}{\pi }^{\prime }·D_3-\frac{1}{\left(\gamma -1\right)}\overline{\pi }·D_3+\frac{1}{\left(\gamma -1\right)}\frac{F_{\overline{\theta }}+F_{{\theta }^{\prime }}}{\overline{\theta }+{\theta }^{\prime }}-\frac{u}{a cos \phi }\frac{\partial \overline{\pi }}{\partial \lambda }-v\frac{\partial \overline{\pi }}{a\partial \phi }-w\frac{\partial \overline{\pi }}{\partial z}$$

(49)

$$\frac{d{\theta }^{\prime }}{dt}=\frac{F_{\overline{\theta }}+F_{{\theta }^{\prime }}}{\overline{\pi }+{\pi }^{\prime }}-\frac{u}{a cos \phi }\frac{\partial \overline{\theta }}{\partial \lambda }-v\frac{\partial \overline{\theta }}{a\partial \phi }-w\frac{\partial \overline{\theta }}{\partial z}$$

(50)

$$T^{\prime }={\pi }^{\prime }{\theta }^{\prime }+\overline{\pi }{\theta }^{\prime }+{\pi }^{\prime }\overline{\theta }$$

(51)

where wind (u, v, w) is full wind with no decomposition. Equations (46)–(51) are the perturbation forms of the governing equations after eliminating the three-dimensional static-temperature reference atmosphere.

If we consider $\overline{\pi }\gg {\pi }^{\prime }$ and $\overline{\theta }\gg {\theta }^{\prime }$, Equations (46)–(51) could be further simplified into:

$$\frac{du}{dt}=-\frac{C_p\overline{\theta }}{a cos \phi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }-\frac{C_p\overline{\theta }}{a cos \phi }\frac{\partial \overline{\pi }}{\partial \lambda }-\frac{C_p{\theta }^{\prime }}{a cos \phi }\frac{\partial \overline{\pi }}{\partial \lambda }+fv+F_u$$

(52)

$$\frac{dv}{dt}=-\frac{C_p\overline{\theta }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }-\frac{C_p\overline{\theta }}{a}\frac{\partial \overline{\pi }}{\partial \phi }-\frac{C_p{\theta }^{\prime }}{a}\frac{\partial \overline{\pi }}{\partial \phi }-fu+F_v$$

(53)

$$\frac{dw}{dt}=-C_p\overline{\theta }\frac{\partial {\pi }^{\prime }}{\partial z}-C_p{\theta }^{\prime }\frac{\partial \overline{\pi }}{\partial z}+F_w$$

(54)

$$\frac{d{\pi }^{\prime }}{dt}=-\frac{1}{\left(\gamma -1\right)}\overline{\pi }·D_3+\frac{1}{\left(\gamma -1\right)}\frac{F_{\overline{\theta }}+F_{{\theta }^{\prime }}}{\overline{\theta }}-\frac{u}{a cos \phi }\frac{\partial \overline{\pi }}{\partial \lambda }-v\frac{\partial \overline{\pi }}{a\partial \phi }-w\frac{\partial \overline{\pi }}{\partial z}$$

(55)

$$\frac{d{\theta }^{\prime }}{dt}=\frac{F_{\overline{\theta }}+F_{{\theta }^{\prime }}}{\overline{\pi }}-\frac{u}{a cos \phi }\frac{\partial \overline{\theta }}{\partial \lambda }-v\frac{\partial \overline{\theta }}{a\partial \phi }-w\frac{\partial \overline{\theta }}{\partial z}$$

(56)

$$T^{\prime }=\overline{\pi }{\theta }^{\prime }+{\pi }^{\prime }\overline{\theta }$$

(57)

All terms on the right side of Equations (52)–(56) are linear. The 3D approach can effectively improve numerical precision in a model’s dynamical core [36,37]. Notably, Wood [35] suggested that using temperature and pressure at a previous moment as reference atmosphere for the next moment in a model introduced time dimension to reference atmosphere.

3. Four-Dimensional All-Variable Climate as a Reference Atmosphere

The three previous reference atmospheres are all only with respect to temperature and pressure and do not evolve in time. Recently, Smolarkiewicz et al. [24] proposed a generalized perturbation form of partial differential equations for all-scale global atmospheric flow, with respect to any assumed “ambient” state. Here, we propose to use a four-dimensional (space and time) all-variable (temperature, pressure, wind, moisture, etc.) climate as reference atmosphere to decompose atmospheric governing equations. For a variable A, its climate ${\tilde{A}}_{ }\left(\lambda ,\phi ,z,t\right)$ is defined as follows:

$$\tilde{A}\left(\lambda ,\phi ,z,t\right)={\displaystyle \sum }_{y=1}^M{A}_y\left(\lambda ,\phi ,z,t\right)/M$$

(58)

which can be derived from global reanalysis data over the past M years (y), $A_y\left(\lambda ,\phi ,z,t\right)$. M is the number of past years and normally set to 30 years or more. Depending on temporal resolution of available reanalysis data, time t can be in seconds, minutes, or hours for every day within a calendar year. The currently available temporal resolution is in hours [38].

Therefore, for a variable A, its anomalous state $A^{\prime }$ can be obtained as follows [39]:

$$A^{\prime }\left(\lambda ,\phi ,z,t\right)=A\left(\lambda ,\phi ,z,t\right)-\tilde{A}\left(\lambda ,\phi ,z,t\right)$$

(59)

For individual variables, they are:

$$u^{\prime }\left(\lambda ,\phi ,z,t\right)=u\left(\lambda ,\phi ,z,t\right)-\tilde{u}\left(\lambda ,\phi ,z,t\right)$$

(60)

$$v^{\prime }\left(\lambda ,\phi ,z,t\right)=v\left(\lambda ,\phi ,z,t\right)-\tilde{v}\left(\lambda ,\phi ,z,t\right)$$

(61)

$$w^{\prime }\left(\lambda ,\phi ,z,t\right)=w\left(\lambda ,\phi ,z,t\right)-\tilde{w}\left(\lambda ,\phi ,z,t\right)$$

(62)

$${\pi }^{\prime }\left(\lambda ,\phi ,z,t\right)=\pi \left(\lambda ,\phi ,z,t\right)-\tilde{\pi }\left(\lambda ,\phi ,z,t\right)$$

(63)

$${\theta }^{\prime }\left(\lambda ,\phi ,z,t\right)=\theta \left(\lambda ,\phi ,z,t\right)-\tilde{\theta }\left(\lambda ,\phi ,z,t\right)$$

(64)

$$T^{\prime }\left(\lambda ,\phi ,z,t\right)=T\left(\lambda ,\phi ,z,t\right)-\tilde{T}\left(\lambda ,\phi ,z,t\right)$$

(65)

$$q^{\prime }\left(\lambda ,\phi ,z,t\right)=q\left(\lambda ,\phi ,z,t\right)-\tilde{q}\left(\lambda ,\phi ,z,t\right)$$

(66)

where q is specific humidity (g/kg).

3.1. Climatic Equations

The climatic states satisfy the atmospheric governing equations (Equations (1)–(6)):

$$\frac{d\tilde{u}}{dt}=\frac{\partial \tilde{u}}{\partial t}+\frac{\tilde{u}}{a cos \phi }\frac{\partial \tilde{u}}{\partial \lambda }+\frac{\tilde{v}}{a}\frac{\partial \tilde{u}}{\partial \phi }+\tilde{w}\frac{\partial \tilde{u}}{\partial z}=-\frac{C_p\tilde{\theta }}{a cos \phi }\frac{\partial \tilde{\pi }}{\partial \lambda }+f\tilde{v}+F_{\tilde{u}}$$

(67)

$$\frac{d\tilde{v}}{dt}=\frac{\partial \tilde{v}}{\partial t}+\frac{\tilde{v}}{a cos \phi }\frac{\partial \tilde{v}}{\partial \lambda }+\frac{\tilde{v}}{a}\frac{\partial \tilde{v}}{\partial \phi }+\tilde{w}\frac{\partial \tilde{v}}{\partial z}=-\frac{C_p\tilde{\theta }}{a}\frac{\partial \tilde{\pi }}{\partial \phi }-f\tilde{u}+F_{\tilde{v}}$$

(68)

$$\frac{d\tilde{w}}{dt}=\frac{\partial \tilde{w}}{\partial t}+\frac{\tilde{u}}{a cos \phi }\frac{\partial \tilde{w}}{\partial \lambda }+\frac{\tilde{v}}{a}\frac{\partial \tilde{w}}{\partial \phi }+\tilde{w}\frac{\partial \tilde{w}}{\partial z}=-C_p\tilde{\theta }\frac{\partial \tilde{\pi }}{\partial z}-g+F_{\tilde{w}}$$

(69)

$$\frac{d\tilde{\pi }}{dt}=\frac{\partial \tilde{\pi }}{\partial t}+\frac{\tilde{u}}{a cos \phi }\frac{\partial \tilde{\pi }}{\partial \lambda }+\frac{\tilde{v}}{a}\frac{\partial \tilde{\pi }}{\partial \phi }+\tilde{w}\frac{\partial \tilde{\pi }}{\partial z}=-\frac{1}{\left(\gamma -1\right)}\tilde{\pi }·{\tilde{D}}_3+\frac{1}{\left(\gamma -1\right)}\frac{F_{\tilde{\theta }}}{\tilde{\theta }}$$

(70)

$$\frac{d\tilde{\theta }}{dt}=\frac{\partial \tilde{\theta }}{\partial t}+\frac{\tilde{u}}{a cos \phi }\frac{\partial \tilde{\theta }}{\partial \lambda }+\frac{\tilde{v}}{a}\frac{\partial \tilde{\theta }}{\partial \phi }+\tilde{w}\frac{\partial \tilde{\theta }}{\partial z}=\frac{F_{\tilde{\theta }}}{\tilde{\pi }}$$

(71)

$$\tilde{T}=\tilde{\pi }\tilde{\theta }$$

(72)

where $F_{\tilde{u}}, F_{\tilde{v}}, F_{\tilde{w}}$ are frictions due to climatic flows, and $F_{\tilde{\theta }}=F_{\tilde{\theta }a}+F_{\tilde{\theta }d}+F_{\tilde{\theta }f}$ represents annual solar cycle ($F_{\tilde{\theta }a}$), daily solar cycle ($F_{\tilde{\theta }d}$), and other climatically constant heating and cooling ($F_{\tilde{\theta }f}$). The climate is an equilibrium state which is maintained by the balance between energy source (radiative heating and other heating) and energy sink (frictions and cooling). If we expand the total derivative $(\frac{d\tilde{A}}{dt})$ into local derivative $(\frac{\partial \tilde{A}}{\partial t})$, the general form of the above climatic equations can be rewritten as follows:

$$\frac{\partial \tilde{A}}{\partial t}=-\frac{\tilde{u}}{a cos \phi }\frac{\partial \tilde{A}}{\partial \lambda }-\frac{\tilde{v}}{a}\frac{\partial \tilde{A}}{\partial \phi }-\tilde{w}\frac{\partial \tilde{A}}{\partial z}+\tilde{N}+\tilde{L}+\tilde{D}$$

(73)

where $-\frac{\tilde{u}}{a cos \phi }\frac{\partial \tilde{A}}{\partial \lambda }-\frac{\tilde{v}}{a}\frac{\partial \tilde{A}}{\partial \phi }-\tilde{w}\frac{\partial \tilde{A}}{\partial z}$ is an advection term; $\tilde{N}$, $\tilde{L}$, and $\tilde{D}$ are nonlinear, linear, and diffusion terms, respectively. Since all terms on the right side can be precalculated from climatology and are known, the local variation on the left side is similar to a slow-motion movie evolving with time and repeats yearly. In other words, the climatic equations (Equations (67)–(71)) do not need to be solved in time (i.e., no prediction is needed for climatic states) but are precalculated from climatology data. This slow-motion movie includes all “averaged” meteorological features such as seasonal monsoon, daily sea breeze, diurnal temperature variation, the intertropical convergence zone (ITCZ), Walker circulation, Hadley circulation, and storm tracks.

3.2. Anomaly Equations

Inserting the decomposition (Equations (60)–(65)) into the full-variable atmospheric governing equations (Equations (1)–(6)), and then being subtracted by the climatic equations (Equations (67)–(72)), the anomaly equations can be obtained as follows:

$$\frac{d{u}^{\prime }}{dt}=\frac{\partial u^{\prime }}{\partial t}+\frac{u}{a cos \phi }\frac{\partial u^{\prime }}{\partial \lambda }+\frac{v}{a}\frac{\partial u^{\prime }}{\partial \phi }+w\frac{\partial u^{\prime }}{\partial z}=\frac{\partial u^{\prime }}{\partial t}+\frac{\tilde{u}+u^{\prime }}{a cos \phi }\frac{\partial u^{\prime }}{\partial \lambda }+\frac{\tilde{v}+v^{\prime }}{a}\frac{\partial u^{\prime }}{\partial \phi }+\left(\tilde{w}+w^{\prime }\right)\frac{\partial u^{\prime }}{\partial z}=-\frac{C_P{\theta }^{\prime }}{a cos \varphi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }-\frac{C_P\tilde{\theta }}{a cos \varphi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }-\frac{C_P{\theta }^{\prime }}{a cos \varphi }\frac{\partial \tilde{\pi }}{\partial \lambda }-\frac{u^{\prime }}{a cos \phi }\frac{\partial \tilde{u}}{\partial \lambda }-v^{\prime }\frac{\partial \tilde{u}}{a\partial \phi }-w^{\prime }\frac{\partial \tilde{u}}{\partial z}+f{v}^{\prime }+F_{u^{\prime }}$$

(74)

$$\frac{d{v}^{\prime }}{dt}=\frac{\partial v^{\prime }}{\partial t}+\frac{u}{a cos \phi }\frac{\partial v^{\prime }}{\partial \lambda }+\frac{v}{a}\frac{\partial v^{\prime }}{\partial \phi }+w\frac{\partial v^{\prime }}{\partial z}=\frac{\partial v^{\prime }}{\partial t}+\frac{\tilde{u}+u^{\prime }}{a cos \phi }\frac{\partial v^{\prime }}{\partial \lambda }+\frac{\tilde{v}+v^{\prime }}{a}\frac{\partial v^{\prime }}{\partial \phi }+\left(\tilde{w}+w^{\prime }\right)\frac{\partial v^{\prime }}{\partial z}=-\frac{C_P{\theta }^{\prime }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }-\frac{C_P\tilde{\theta }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }-\frac{C_P{\theta }^{\prime }}{a}\frac{\partial \tilde{\pi }}{\partial \phi }-\frac{u^{\prime }}{a cos \phi }\frac{\partial \tilde{v}}{\partial \lambda }-v^{\prime }\frac{\partial \tilde{v}}{a\partial \phi }-w^{\prime }\frac{\partial \check{v}}{\partial z}-f{u}^{\prime }+F_{v^{\prime }}$$

(75)

$$\frac{d{w}^{\prime }}{dt}=\frac{\partial w^{\prime }}{\partial t}+\frac{u}{a cos \phi }\frac{\partial w^{\prime }}{\partial \lambda }+\frac{v}{a}\frac{\partial w^{\prime }}{\partial \phi }+w\frac{\partial w^{\prime }}{\partial z}=\frac{\partial w^{\prime }}{\partial t}+\frac{\tilde{u}+u^{\prime }}{a cos \phi }\frac{\partial w^{\prime }}{\partial \lambda }+\frac{\tilde{v}+v^{\prime }}{a}\frac{\partial w^{\prime }}{\partial \phi }+\left(\tilde{w}+w^{\prime }\right)\frac{\partial w^{\prime }}{\partial z}=-C_P{\theta }^{\prime }\frac{\partial {\pi }^{\prime }}{\partial z}-C_P\tilde{\theta }\frac{\partial {\pi }^{\prime }}{\partial z}-C_P{\theta }^{\prime }\frac{\partial \tilde{\pi }}{\partial z}-\frac{u^{\prime }}{a cos \phi }\frac{\partial \tilde{w}}{\partial \lambda }-v^{\prime }\frac{\partial \tilde{w}}{a\partial \phi }-w^{\prime }\frac{\partial \check{w}}{\partial z}+F_{w^{\prime }}$$

(76)

$$\frac{d{\pi }^{\prime }}{dt}=\frac{\partial {\pi }^{\prime }}{\partial t}+\frac{u}{a cos \phi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }+\frac{v}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }+w\frac{\partial {\pi }^{\prime }}{\partial z}=\frac{\partial {\pi }^{\prime }}{\partial t}+\frac{\tilde{u}+u^{\prime }}{a cos \phi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }+\frac{\tilde{v}+v^{\prime }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }+\left(\tilde{w}+w^{\prime }\right)\frac{\partial {\pi }^{\prime }}{\partial z}=-\frac{1}{\left(\gamma -1\right)}{\pi }^{\prime }\cdot D_3^{\prime }-\frac{1}{\left(\gamma -1\right)}\tilde{\pi }\cdot D_3^{\prime }-\frac{1}{\left(\gamma -1\right)}{\pi }^{\prime }\cdot {\tilde{D}}_3+\frac{1}{\left(\gamma -1\right)}\frac{F_{\tilde{\theta }}}{\tilde{\theta }+{\theta }^{\prime }}-\frac{1}{\left(\gamma -1\right)}\frac{F_{\tilde{\theta }}}{\tilde{\theta }}+\frac{1}{\left(\gamma -1\right)}\frac{F_{{\theta }^{\prime }}}{\tilde{\theta }+{\theta }^{\prime }}-\frac{u^{\prime }}{a cos \phi }\frac{\partial \tilde{\pi }}{\partial \lambda }-v^{\prime }\frac{\partial \tilde{\pi }}{a\partial \phi }-w^{\prime }\frac{\partial \check{\pi }}{\partial z}$$

(77)

$$\frac{d{\theta }^{\prime }}{dt}=\frac{\partial {\theta }^{\prime }}{\partial t}+\frac{u}{a cos \phi }\frac{\partial {\theta }^{\prime }}{\partial \lambda }+\frac{v}{a}\frac{\partial {\theta }^{\prime }}{\partial \phi }+w\frac{\partial {\theta }^{\prime }}{\partial z}=\frac{\partial {\theta }^{\prime }}{\partial t}+\frac{\tilde{u}+u^{\prime }}{a cos \phi }\frac{\partial {\theta }^{\prime }}{\partial \lambda }+\frac{\tilde{v}+v^{\prime }}{a}\frac{\partial {\theta }^{\prime }}{\partial \phi }+\left(\tilde{w}+w^{\prime }\right)\frac{\partial {\theta }^{\prime }}{\partial z}=\frac{F_{\tilde{\theta }}}{\tilde{\pi }+{\pi }^{\prime }}-\frac{F_{\tilde{\theta }}}{\tilde{\pi }}+\frac{F_{{\theta }^{\prime }}}{\tilde{\pi }+{\pi }^{\prime }}-\frac{u^{\prime }}{acos\phi }\frac{\partial \tilde{\theta }}{\partial \lambda }-v^{\prime }\frac{\partial \tilde{\theta }}{a\partial \phi }-w^{\prime }\frac{\partial \tilde{\theta }}{\partial z}$$

(78)

$$T^{\prime }={\pi }^{\prime }{\theta }^{\prime }+\tilde{\pi }{\theta }^{\prime }+{\pi }^{\prime }\tilde{\theta }$$

(79)

Equations (74)–(79) are the perturbation form of the governing equations after eliminating the four-dimensional all-variable reference atmosphere. Similar to $F_{\theta }$, friction term $F_u$ is also decomposed into climatic and perturbed terms, i.e., $F_u=F_{\tilde{u}}+F_{u^{\prime }}$. Since the perturbations are relative to climatic states, these perturbation equations are meteorologically “anomaly equations”. Note that wind in the advection term on the left side is full wind, i.e., $\frac{d{A}^{\prime }}{dt}=\frac{\partial A^{\prime }}{\partial t}+\frac{u}{a cos \phi }\frac{\partial A^{\prime }}{\partial \lambda }+\frac{v}{a}\frac{\partial A^{\prime }}{\partial \phi }+w\frac{\partial A^{\prime }}{\partial z}=\frac{\partial A^{\prime }}{\partial t}+\frac{\tilde{u}+u^{\prime }}{a cos \phi }\frac{\partial A^{\prime }}{\partial \lambda }+\frac{\tilde{v}+v^{\prime }}{a}\frac{\partial A^{\prime }}{\partial \phi }+\left(\tilde{w}+w^{\prime }\right)\frac{\partial A^{\prime }}{\partial z}$.

If we consider $\overline{\pi }\gg {\pi }^{\prime }$ and $\overline{\theta }\gg {\theta }^{\prime }$ (note, wind perturbation is not necessarily smaller than background wind, and sometimes could even exceeds the background wind such as in a storm area), Equations (74)–(79) could be further simplified into:

$$\frac{d{u}^{\prime }}{dt}=\frac{\partial u^{\prime }}{\partial t}+\frac{\tilde{u}+u^{\prime }}{a cos \phi }\frac{\partial u^{\prime }}{\partial \lambda }+\frac{\tilde{v}+v^{\prime }}{a}\frac{\partial u^{\prime }}{\partial \phi }+\left(\tilde{w}+w^{\prime }\right)\frac{\partial u^{\prime }}{\partial z}=-\frac{C_P\tilde{\theta }}{a cos \varphi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }-\frac{C_P{\theta }^{\prime }}{a cos \varphi }\frac{\partial \tilde{\pi }}{\partial \lambda }-\frac{u^{\prime }}{a cos \phi }\frac{\partial \tilde{u}}{\partial \lambda }-v^{\prime }\frac{\partial \tilde{u}}{a\partial \phi }-w^{\prime }\frac{\partial \tilde{u}}{\partial z}+f{v}^{\prime }+F_{u^{\prime }}$$

(80)

$$\frac{d{v}^{\prime }}{dt}=\frac{\partial v^{\prime }}{\partial t}+\frac{\tilde{u}+u^{\prime }}{a cos \phi }\frac{\partial v^{\prime }}{\partial \lambda }+\frac{\tilde{v}+v^{\prime }}{a}\frac{\partial v^{\prime }}{\partial \phi }+\left(\tilde{w}+w^{\prime }\right)\frac{\partial v^{\prime }}{\partial z}=-\frac{C_P\tilde{\theta }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }-\frac{C_P{\theta }^{\prime }}{a}\frac{\partial \tilde{\pi }}{\partial \phi }-\frac{u^{\prime }}{a cos \phi }\frac{\partial \tilde{v}}{\partial \lambda }-v^{\prime }\frac{\partial \tilde{v}}{a\partial \phi }-w^{\prime }\frac{\partial \check{v}}{\partial z}-f{u}^{\prime }+F_{v^{\prime }}$$

(81)

$$\frac{d{w}^{\prime }}{dt}=\frac{\partial w^{\prime }}{\partial t}+\frac{\tilde{u}+u^{\prime }}{a cos \phi }\frac{\partial w^{\prime }}{\partial \lambda }+\frac{\tilde{v}+v^{\prime }}{a}\frac{\partial w^{\prime }}{\partial \phi }+\left(\tilde{w}+w^{\prime }\right)\frac{\partial w^{\prime }}{\partial z}=-C_P\tilde{\theta }\frac{\partial {\pi }^{\prime }}{\partial z}-C_P{\theta }^{\prime }\frac{\partial \tilde{\pi }}{\partial z}-\frac{u^{\prime }}{a cos \phi }\frac{\partial \tilde{w}}{\partial \lambda }-v^{\prime }\frac{\partial \tilde{w}}{a\partial \phi }-w^{\prime }\frac{\partial \check{w}}{\partial z}+F_{w^{\prime }}$$

(82)

$$\frac{d{\pi }^{\prime }}{dt}=\frac{\partial {\pi }^{\prime }}{\partial t}+\frac{\tilde{u}+u^{\prime }}{a cos \phi }\frac{\partial {\pi }^{\prime }}{\partial \lambda }+\frac{\tilde{v}+v^{\prime }}{a}\frac{\partial {\pi }^{\prime }}{\partial \phi }+\left(\tilde{w}+w^{\prime }\right)\frac{\partial {\pi }^{\prime }}{\partial z}=-\frac{1}{\left(\gamma -1\right)}\tilde{\pi }\cdot D_3^{\prime }-\frac{1}{\left(\gamma -1\right)}{\pi }^{\prime }\cdot {\tilde{D}}_3+\frac{1}{\left(\gamma -1\right)}\frac{F_{{\theta }^{\prime }}}{\tilde{\theta }}-\frac{u^{\prime }}{a cos \phi }\frac{\partial \tilde{\pi }}{\partial \lambda }-v^{\prime }\frac{\partial \tilde{\pi }}{a\partial \phi }-w^{\prime }\frac{\partial \check{\pi }}{\partial z}$$

(83)

$$\frac{d{\theta }^{\prime }}{dt}=\frac{\partial {\theta }^{\prime }}{\partial t}+\frac{\tilde{u}+u^{\prime }}{a cos \phi }\frac{\partial {\theta }^{\prime }}{\partial \lambda }+\frac{\tilde{v}+v^{\prime }}{a}\frac{\partial {\theta }^{\prime }}{\partial \phi }+\left(\tilde{w}+w^{\prime }\right)\frac{\partial {\theta }^{\prime }}{\partial z}=\frac{F_{{\theta }^{\prime }}}{\tilde{\pi }}-\frac{u^{\prime }}{a cos \phi }\frac{\partial \tilde{\theta }}{\partial \lambda }-v^{\prime }\frac{\partial \tilde{\theta }}{a\partial \phi }-w^{\prime }\frac{\partial \tilde{\theta }}{\partial z}$$

(84)

$$T^{\prime }=\tilde{\pi }{\theta }^{\prime }+{\pi }^{\prime }\tilde{\theta }$$

(85)

All terms on the far-right side of the prognostic Equations (80)–(84) are linear.

The x-y coordinate is commonly used for regional models. Under the x-y coordinate, the complete anomaly equations (Equations (74)–(79)) could be written as follows:

$$\frac{\partial u^{\prime }}{\partial t}=-[u^{\prime }\frac{\partial u^{\prime }}{\partial x}+v^{\prime }\frac{\partial u^{\prime }}{\partial y}+w^{\prime }\frac{\partial u^{\prime }}{\partial z}]-[u^{\prime }\frac{\partial \tilde{u}}{\partial x}+v^{\prime }\frac{\partial \tilde{u}}{\partial y}+w^{\prime }\frac{\partial \tilde{u}}{\partial z}]-[\tilde{u}\frac{\partial u^{\prime }}{\partial x}+\tilde{v}\frac{\partial u^{\prime }}{\partial y}+\tilde{w}\frac{\partial u^{\prime }}{\partial z}]-[C_P{\theta }^{\prime }\frac{\partial {\pi }^{\prime }}{\partial x}+C_P\tilde{\theta }\frac{\partial {\pi }^{\prime }}{\partial x}+C_P{\theta }^{\prime }\frac{\partial \tilde{\pi }}{\partial x}]+f{v}^{\prime }+F_{u^{\prime }}$$

(86)

$$\frac{\partial v^{\prime }}{\partial t}=-\left[u^{\prime }\frac{\partial v^{\prime }}{\partial x}+v^{\prime }\frac{\partial v^{\prime }}{\partial y}+w^{\prime }\frac{\partial v^{\prime }}{\partial z}\right]-\left[u^{\prime }\frac{\partial \tilde{v}}{\partial x}+v^{\prime }\frac{\partial \tilde{v}}{\partial y}+w^{\prime }\frac{\partial \tilde{v}}{\partial z}\right]-\left[\tilde{u}\frac{\partial v^{\prime }}{\partial x}+\tilde{v}\frac{\partial v^{\prime }}{\partial y}+\tilde{w}\frac{\partial v^{\prime }}{\partial z}\right]-[C_P{\theta }^{\prime }\frac{\partial {\pi }^{\prime }}{\partial y}+C_P\tilde{\theta }\frac{\partial {\pi }^{\prime }}{\partial y}+C_P{\theta }^{\prime }\frac{\partial \tilde{\pi }}{\partial y}]-f{u}^{\prime }+F_{v^{\prime }}$$

(87)

$$\frac{\partial w^{\prime }}{\partial t}=-\left[u^{\prime }\frac{\partial w^{\prime }}{\partial x}+v^{\prime }\frac{\partial w^{\prime }}{\partial y}+w^{\prime }\frac{\partial w^{\prime }}{\partial z}\right]-\left[u^{\prime }\frac{\partial \tilde{w}}{\partial x}+v^{\prime }\frac{\partial \tilde{w}}{\partial y}+w^{\prime }\frac{\partial \tilde{w}}{\partial z}\right]-\left[\tilde{u}\frac{\partial w^{\prime }}{\partial x}+\tilde{v}\frac{\partial w^{\prime }}{\partial y}+\tilde{w}\frac{\partial w^{\prime }}{\partial z}\right]-\left[C_P{\theta }^{\prime }\frac{\partial {\pi }^{\prime }}{\partial z}+C_P\tilde{\theta }\frac{\partial {\pi }^{\prime }}{\partial z}+C_P{\theta }^{\prime }\frac{\partial \tilde{\pi }}{\partial z}\right]+F_{w^{\prime }}$$

(88)

$$\frac{\partial {\pi }^{\prime }}{\partial t}=-\left[u^{\prime }\frac{\partial {\pi }^{\prime }}{\partial x}+v^{\prime }\frac{\partial {\pi }^{\prime }}{\partial y}+w^{\prime }\frac{\partial {\pi }^{\prime }}{\partial z}\right]-\left[u^{\prime }\frac{\partial \tilde{\pi }}{\partial x}+v^{\prime }\frac{\partial \tilde{\pi }}{\partial y}+w^{\prime }\frac{\partial \tilde{\pi }}{\partial z}\right]-\left[\tilde{u}\frac{\partial {\pi }^{\prime }}{\partial x}+\tilde{v}\frac{\partial {\pi }^{\prime }}{\partial y}+\tilde{w}\frac{\partial {\pi }^{\prime }}{\partial z}\right]-\left[\frac{1}{\left(\gamma -1\right)}{\pi }^{\prime }\cdot D_3^{\prime }+\frac{1}{\left(\gamma -1\right)}\tilde{\pi }\cdot D_3^{\prime }+\frac{1}{\left(\gamma -1\right)}{\pi }^{\prime }\cdot {\tilde{D}}_3\right]+[\frac{1}{\left(\gamma -1\right)}\frac{F_{\tilde{\theta }}}{\tilde{\theta }+{\theta }^{\prime }}-\frac{1}{\left(\gamma -1\right)}\frac{F_{\tilde{\theta }}}{\tilde{\theta }}+\frac{1}{\left(\gamma -1\right)}\frac{F_{{\theta }^{\prime }}}{\tilde{\theta }+{\theta }^{\prime }}]$$

(89)

$$\frac{\partial {\theta }^{\prime }}{\partial t}=-\left[u^{\prime }\frac{\partial {\theta }^{\prime }}{\partial x}+v^{\prime }\frac{\partial {\theta }^{\prime }}{\partial y}+w^{\prime }\frac{\partial {\theta }^{\prime }}{\partial z}\right]-\left[u^{\prime }\frac{\partial \tilde{\theta }}{\partial x}+v^{\prime }\frac{\partial \tilde{\theta }}{\partial y}+w^{\prime }\frac{\partial \tilde{\theta }}{\partial z}\right]-\left[\tilde{u}\frac{\partial {\theta }^{\prime }}{\partial x}+\tilde{v}\frac{\partial {\theta }^{\prime }}{\partial y}+\tilde{w}\frac{\partial {\theta }^{\prime }}{\partial z}\right]+[\frac{F_{\tilde{\theta }}}{\tilde{\pi }+{\pi }^{\prime }}-\frac{F_{\tilde{\theta }}}{\tilde{\pi }}+\frac{F_{{\theta }^{\prime }}}{\tilde{\pi }+{\pi }^{\prime }}]$$

(90)

$$T^{\prime }={\pi }^{\prime }{\theta }^{\prime }+\tilde{\pi }{\theta }^{\prime }+{\pi }^{\prime }\tilde{\theta }$$

(91)

The motion equations (Equations (86)–(88)) show that the local change of anomaly flow is controlled by six factors: interaction with anomalous environment, interaction with climatic flow, advection by climatic flow, pressure gradient force, earth rotation effect (Coriolis force, except for vertical motion), and dissipation to friction. The local pressure anomaly change (Equation (89)) is controlled by five factors: interaction with anomalous environment, interaction with climatic pressure, advection by climatic flow, air inflow and outflow, and adiabatic heating and cooling. The local temperature anomaly change (Equation (90)) is controlled by four factors: interaction with anomalous environment, interaction with background temperature, advection by climatic flow, and adiabatic heating and cooling. Since climatic states are precalculated and given, all terms related to climate on the right side are linear.

Under the assumption of $\overline{\pi }\gg {\pi }^{\prime }$ and $\overline{\theta }\gg {\theta }^{\prime }$, the simplified version is:

$$\frac{\partial u^{\prime }}{\partial t}=-[u^{\prime }\frac{\partial u^{\prime }}{\partial x}+v^{\prime }\frac{\partial u^{\prime }}{\partial y}+w^{\prime }\frac{\partial u^{\prime }}{\partial z}]-[u^{\prime }\frac{\partial \tilde{u}}{\partial x}+v^{\prime }\frac{\partial \tilde{u}}{\partial y}+w^{\prime }\frac{\partial \tilde{u}}{\partial z}]-[\tilde{u}\frac{\partial u^{\prime }}{\partial x}+\tilde{v}\frac{\partial u^{\prime }}{\partial y}+\tilde{w}\frac{\partial u^{\prime }}{\partial z}]-[C_P\tilde{\theta }\frac{\partial {\pi }^{\prime }}{\partial x}+C_P{\theta }^{\prime }\frac{\partial \tilde{\pi }}{\partial x}]+f{v}^{\prime }+F_{u^{\prime }}$$

(92)

$$\frac{\partial v^{\prime }}{\partial t}=-\left[u^{\prime }\frac{\partial v^{\prime }}{\partial x}+v^{\prime }\frac{\partial v^{\prime }}{\partial y}+w^{\prime }\frac{\partial v^{\prime }}{\partial z}\right]-\left[u^{\prime }\frac{\partial \tilde{v}}{\partial x}+v^{\prime }\frac{\partial \tilde{v}}{\partial y}+w^{\prime }\frac{\partial \tilde{v}}{\partial z}\right]-\left[\tilde{u}\frac{\partial v^{\prime }}{\partial x}+\tilde{v}\frac{\partial v^{\prime }}{\partial y}+\tilde{w}\frac{\partial v^{\prime }}{\partial z}\right]-[C_P\tilde{\theta }\frac{\partial {\pi }^{\prime }}{\partial y}+C_P{\theta }^{\prime }\frac{\partial \tilde{\pi }}{\partial y}]-f{u}^{\prime }+F_{v^{\prime }}$$

(93)

$$\frac{\partial w^{\prime }}{\partial t}=-\left[u^{\prime }\frac{\partial w^{\prime }}{\partial x}+v^{\prime }\frac{\partial w^{\prime }}{\partial y}+w^{\prime }\frac{\partial w^{\prime }}{\partial z}\right]-\left[u^{\prime }\frac{\partial \tilde{w}}{\partial x}+v^{\prime }\frac{\partial \tilde{w}}{\partial y}+w^{\prime }\frac{\partial \tilde{w}}{\partial z}\right]-\left[\tilde{u}\frac{\partial w^{\prime }}{\partial x}+\tilde{v}\frac{\partial w^{\prime }}{\partial y}+\tilde{w}\frac{\partial w^{\prime }}{\partial z}\right]-\left[C_P\tilde{\theta }\frac{\partial {\pi }^{\prime }}{\partial z}+C_P{\theta }^{\prime }\frac{\partial \tilde{\pi }}{\partial z}\right]+F_{w^{\prime }}$$

(94)

$$\frac{\partial {\pi }^{\prime }}{\partial t}=-\left[u^{\prime }\frac{\partial {\pi }^{\prime }}{\partial x}+v^{\prime }\frac{\partial {\pi }^{\prime }}{\partial y}+w^{\prime }\frac{\partial {\pi }^{\prime }}{\partial z}\right]-\left[u^{\prime }\frac{\partial \tilde{\pi }}{\partial x}+v^{\prime }\frac{\partial \tilde{\pi }}{\partial y}+w^{\prime }\frac{\partial \tilde{\pi }}{\partial z}\right]-\left[\tilde{u}\frac{\partial {\pi }^{\prime }}{\partial x}+\tilde{v}\frac{\partial {\pi }^{\prime }}{\partial y}+\tilde{w}\frac{\partial {\pi }^{\prime }}{\partial z}\right]-\left[\frac{1}{\left(\gamma -1\right)}\tilde{\pi }·D_3^{\prime }+\frac{1}{\left(\gamma -1\right)}{\pi }^{\prime }·{\tilde{D}}_3\right]+\frac{1}{\left(\gamma -1\right)}\frac{F_{{\theta }^{\prime }}}{\tilde{\theta }}$$

(95)

$$\frac{\partial {\theta }^{\prime }}{\partial t}=-\left[u^{\prime }\frac{\partial {\theta }^{\prime }}{\partial x}+v^{\prime }\frac{\partial {\theta }^{\prime }}{\partial y}+w^{\prime }\frac{\partial {\theta }^{\prime }}{\partial z}\right]-\left[u^{\prime }\frac{\partial \tilde{\theta }}{\partial x}+v^{\prime }\frac{\partial \tilde{\theta }}{\partial y}+w^{\prime }\frac{\partial \tilde{\theta }}{\partial z}\right]-\left[\tilde{u}\frac{\partial {\theta }^{\prime }}{\partial x}+\tilde{v}\frac{\partial {\theta }^{\prime }}{\partial y}+\tilde{w}\frac{\partial {\theta }^{\prime }}{\partial z}\right]+\frac{F_{{\theta }^{\prime }}}{\tilde{\pi }}$$

(96)

$$T^{\prime }=\tilde{\pi }{\theta }^{\prime }+{\pi }^{\prime }\tilde{\theta }$$

(97)

The anomaly equations have been discussed so far are all in dry version with no moisture. In real world, moisture plays an important role in daily evolving weather. Below gives a derivation of anomaly moisture equation. In x-y coordinate, the full-field moisture equation is:

$$\frac{\partial q}{\partial t}=-\left[u^{ }\frac{\partial q}{\partial x}+v^{ }\frac{\partial q}{\partial y}+w^{ }\frac{\partial q}{\partial z}\right]+Q_{source}-Q_{sink}$$

(98)

where $Q_{source}$ is moisture source such as evaporation and evapotranspiration, and $Q_{sink}$ moisture sink such as condensation and precipitation. By decomposing specific humidity, $Q_{source}$, and $Q_{sink}$ into climatic and anomalous components:

$$q^{\prime }\left(x,y,z,t\right)=q\left(x,y,z,t\right)-\tilde{q}\left(x,y,z,t\right)$$

(99)

$$Q_{source}\left(x,y,x,t\right)={\tilde{Q}}_{source}\left(x,y,z,t\right)+{Q^{\prime }}_{source}\left(x,y,z,t\right)$$

(100)

$$Q_{sink}\left(x,y,z,t\right)={\tilde{Q}}_{sink}\left(x,y,z,t\right)+{Q^{\prime }}_{sink}\left(x,y,z,t\right)$$

(101)

We have climatic moisture equation:

$$\frac{\partial \tilde{q}}{\partial t}=-\left[\tilde{u}\frac{\partial \tilde{q}}{\partial x}+\tilde{v}\frac{\partial \tilde{q}}{\partial y}+\tilde{w}\frac{\partial \tilde{q}}{\partial z}\right]+{\tilde{Q}}_{source}-{\tilde{Q}}_{sink}$$

(102)

Decomposing the full equation (Equation (98)) by Equations (99)–(101) etc., and then subtracting the climatic equation (Equation (102)), the anomaly moisture equation can be obtained:

$$\frac{\partial q^{\prime }}{\partial t}=-\left[u^{\prime }\frac{\partial q^{\prime }}{\partial x}+v^{\prime }\frac{\partial q^{\prime }}{\partial y}+w^{\prime }\frac{\partial q^{\prime }}{\partial z}\right]-\left[u^{\prime }\frac{\partial \tilde{q}}{\partial x}+v^{\prime }\frac{\partial \tilde{q}}{\partial y}+w^{\prime }\frac{\partial \tilde{q}}{\partial z}\right]-\left[\tilde{u}\frac{\partial q^{\prime }}{\partial x}+\tilde{v}\frac{\partial q^{\prime }}{\partial y}+\tilde{w}\frac{\partial q^{\prime }}{\partial z}\right]+{Q^{\prime }}_{source}-{Q^{\prime }}_{sink}$$

(103)

In addition, for saturated moist air, the potential temperature $\theta$ needs to be replaced by equivalent potential temperature ${\theta }_e$ in the anomaly equations (such as Equations (86)–(91)):

$${\theta }_e=\theta \exp (\frac{L{q}_s}{C_{p}T})$$

(104)

where L is latent heat, and $q_s$ saturated specific humidity. These moisture-related terms are connected to many complex physical processes such as evaporation, evapotranspiration, land surface features (soil moisture, vegetation, lakes, and ocean etc.), convection (condensation and latent heat), cloud micro-physics (water phase changes), and different types of precipitation (liquid and solid).

As for physical or meteorological meaning, the 4D climate-decomposition-based anomaly equations derived in this section clearly have advantages over the other three versions (1D, 2D, and 3D decomposition based) of perturbation equations derived in Section 2. The perturbations from the previous three versions are only departures from their own reference atmosphere and do not have direct connection to predicted individual weather systems, while the anomalies from the climate-decomposition-based anomaly equations are directly a part of predicted individual weather systems both in structure and strength [38].

3.3. Advantage and Challenge

Since background atmosphere has been removed in anomaly prognostic equations, an anomaly-equation-based model can avoid numerical errors accumulated from background atmosphere prediction and, therefore, might improve forecast accuracy. Theoretically speaking, anomaly model should greatly reduce model systematic errors (bias) and avoid model climate drift problem. Qian et al. [40] and Huang et al. [41] have demonstrated this advantage in hurricane track prediction using a vorticity equation.

At a non-divergence level, the full-field vorticity equation (with no friction) is:

$$\frac{\partial \varsigma }{\partial t}=-u\frac{\partial \varsigma }{\partial x}-v\frac{\partial \varsigma }{\partial y}-\beta v$$

(105)

where $\zeta =\frac{\partial v}{\partial x}-\frac{\partial u}{\partial y}$, $\beta =\left(df/dy\right)=2\Omega cos \phi /a$. Equation (105) is the full-field-based Beta advection model (BAM).

The climatic equation is:

$$\frac{\partial \tilde{\zeta }}{\partial t}=-\tilde{u}\frac{\partial \tilde{\zeta }}{\partial x}-\tilde{v}\frac{\partial \tilde{\zeta }}{\partial y}-\beta \tilde{v}$$

(106)

where $\tilde{\zeta }=\frac{\partial \tilde{v}}{\partial x}-\frac{\partial \tilde{u}}{\partial y}$.

By subtracting Equation (106) from Equation (105), the anomaly equation can be obtained:

$$\frac{\partial {\varsigma }^{\prime }}{\partial t}=-[\left(\tilde{u}+u^{\prime }\right)\frac{\partial {\varsigma }^{\prime }}{\partial x}+(\tilde{v}+v^{\prime })\frac{\partial {\varsigma }^{\prime }}{\partial y}]-[u^{\prime }\frac{\partial \tilde{\zeta }}{\partial x}+v^{\prime }\frac{\partial \tilde{\zeta }}{\partial y}]-\beta v^{\prime }$$

(107)

where ${\varsigma }^{\prime }=\frac{\partial v^{\prime }}{\partial x}-\frac{\partial u^{\prime }}{\partial y}$. Equation (107) is the anomaly-based Beta advection model. In an environment with strong tropical storms, Huang et al. [41] showed that $u^{\prime }\frac{\partial \tilde{\zeta }}{\partial x}+v^{\prime }\frac{\partial \tilde{\zeta }}{\partial y}\ll \left(\tilde{u}+u^{\prime }\right)\frac{\partial {\varsigma }^{\prime }}{\partial x}+(\tilde{v}+v^{\prime })\frac{\partial {\varsigma }^{\prime }}{\partial y}$, Equation (107) can then be simplified as:

$$\frac{\partial {\varsigma }^{\prime }}{\partial t}=-[\left(\tilde{u}+u^{\prime }\right)\frac{\partial {\varsigma }^{\prime }}{\partial x}+(\tilde{v}+v^{\prime })\frac{\partial {\varsigma }^{\prime }}{\partial y}]-[u^{\prime }\frac{\partial \tilde{\zeta }}{\partial x}+v^{\prime }\frac{\partial \tilde{\zeta }}{\partial y}]-\beta v^{\prime }$$

(108)

This simplified version of the anomaly-based Beta advection model (Equation (108)) is called the generalized Beta advection model (GBAM) and used in Qian et al. [40] and Huang et al. [41].

Qian et al. [40] and Huang et al. [41] demonstrated that using the anomaly-based GBAM could more accurately predict tropical cyclone’s track than the full-field-based BAM. For example, Huang et al. [41] reported that the GBAM steadily outperformed the BAMs for both normal and unusual tracks after examining 133 tropical cyclones. One of the reasons why the GBAM is superior to the BAM is that the steering flow level can be precisely estimated from anomalous vorticity or divergence field. The track prediction is more accurate, and it is also easier to perform sensitivity experiments to investigate controlling factors of storm track because the GBAM has two clearly distinct terms: one is the interaction with surrounding anomalous environment ($-u^{\prime }\frac{\partial {\varsigma }^{\prime }}{\partial x}-v^{\prime }\frac{\partial {\varsigma }^{\prime }}{\partial y}$), and another is the impact of background flow ($-\tilde{u}\frac{\partial {\varsigma }^{\prime }}{\partial x}-\tilde{v}\frac{\partial {\varsigma }^{\prime }}{\partial y}$), while the BAM has more complex interaction terms.

Although for simple pure-dynamic anomaly models such as the GBAM, it seems to work well, it will be a challenge to apply complete anomaly equations (such as Equations (74)–(79) or Equations (86)–(91)) to the state-of-the-art NWP model. We will try our best to elaborate some of them since it is an untouched area till now. The first task is to develop a dynamical core with anomaly equations. Since climatic states are precalculated, most of terms on the right side of anomaly prognostic equations are linear with only a few nonlinear terms. For those linear and nonlinear terms, we might already have matured numerical schemes [42,43,44] in NWP to solve them. In a decomposed anomaly climate mixed NWP model, only anomaly-related terms need to be integrated or predicted over time, while climatic terms do not need to be calculated but are precalculated constants. Although state variables are decomposed, it is still a full-field model and, therefore, works with current physics package. This approach is an analog to the perturbation model developed by Juang at NCEP [45]. The perturbation model produces a perturbed forecast, and then adds the perturbed forecast to a background forecast (a global model forecast) to have a more accurate downscaled forecast. In this case, the background is observed climate instead of a global model forecast. To realize this method, very high-resolution (both in space and time) climatology or reanalysis data would be needed, which is currently lacking.

What if a pure anomaly-equation-based NWP model such as the BAMS [40,41] used in tropical cyclone’s track prediction is desired, where climatic equations are excluded? In such a pure anomaly model, physics processes (parameterized or AI based) must be decomposed into climatic and anomalous ones too to match model dynamics. This physics separation would be a significant challenge which has not yet been explored. Physics is an important part of the anomaly equations. For example, the terms $F_{u^{\prime }}$, $F_{v^{\prime }}$, and $F_{w^{\prime }}$ in Equations (86)–(88) are frictions, which are related to “anomalous” surface layer, planetary boundary layer, turbulence, etc. physical processes. The term $F_{{\theta }^{\prime }}$ in Equations (89)–(90) is heating and cooling effects, which are related to “anomalous” radiation, and adiabatic heating and cooling processes. ${Q^{\prime }}_{source}$ and ${Q^{\prime }}_{sink}$ in Equation (103) are related to anomalous evaporation and precipitation processes. Some physical processes might be easier to decompose than others. For example, average solar radiation can be easily identified, and topography including sea-land distribution is always climatic (means that no topography is needed in an anomaly-equation-based model). While some others, including multi-scale processes, will be challenging to separate them into climatic and anomaly portions, not to mention that such a separation of physics might not be realistic in the real world. Nowadays, a longer-range NWP model (earth simulator) is often coupled with ocean, cryosphere, biosphere, and aerosol, so that chemical and biological processes are also involved in addition to physical processes. This makes such a separation even more difficult, if not impossible. One day in the future, if all model physics can be explicitly resolved by a model itself, this challenging physics decomposition issue might be overcome. For forecasts with strong flow dynamics but less physics feedback, an anomaly form of pure dynamical model (without physics) might work well, as demonstrated by the 1–2 days tropical cyclone’s track forecasts [40,41]. Such forecasts might include very short-range weather forecasts of certain types or large-scale seasonal forecasts with strong anomaly ENSO-like forcing.

4. Summary

To reduce numerical instability and increase NWP model forecast accuracy, one approach is to avoid or eliminate large-value background terms. Atmospheric scientists had proposed the 1D, 2D, and 3D static reference atmospheres with respect to temperature and pressure in the past. By subtracting a reference atmosphere from the atmospheric governing equations, perturbation equations can be constructed. The perturbation-equation-based model does run more stable numerically and have higher accuracy. These three reference atmospheres were reviewed in this paper. Perturbation equations under these three reference atmospheres were derived accordingly for reader’s reference.

However, the existing 1D, 2D, and 3D reference atmospheres are only with respect to temperature and pressure with no time evolution. The resulting perturbations have no direct connection to weather systems. Therefore, we proposed a four-dimensional (space and time) all-variable (temperature, pressure, wind, moisture, etc.) climatic state as a reference atmosphere to decompose the atmospheric governing equations. In this way, the resulting perturbations are anomalies relative to climate and directly a part of individual weather systems in both structure and strength. The anomaly equations were derived and analyzed.

Finally, the benefits and challenges of the anomaly-equation-based NWP model were discussed. Theoretically, an anomaly model should greatly reduce model systematic errors (bias) and avoid model climate drift problem due to the removal of climatic flow prediction. The anomaly-based Beta advection model (vorticity equation) did consistently outperform the full-field-based model in tropical cyclone track forecasts. However, if a purely anomaly-equation-based model is applied to a real-world NWP model, the model’s complex physics will also need to be separated into climatic physics (physics tendency) and anomalous physics (physics tendency). This separation of physics will be a significant challenge, not to mention if it is physically a valid option in real world. Fortunately, if both anomaly and climatic equations are included in an NWP model, it is still a full-field model. Therefore, model physics should still work as it does now. In such a decomposed anomaly climate mixed NWP model, only anomaly-related terms need to be integrated or predicted over time, while climatic terms do not need to be calculated but are precalculated constants. This anomaly climate mixed approach requires very high-resolution (both in space and time) climatology or reanalysis data, which is lacking now. Developing a high spatial and temporal resolution global reanalysis dataset is a key priority. This research might have opened a new door for NWP model development. For now, we encourage model developers to establish an anomaly climate mixed NWP model based on the anomaly and climatic equations derived in this study, which should advance next-generation NWP model development.