Inherently mass-conservative version of the semi-Lagrangian absolute vorticity (SL-AV) atmospheric model dynamical core

The semi-Lagrangian absolute vorticity (SL-AV) atmospheric model is the global semi-Lagrangian hydrostatic model used for operational medium-range and seasonal forecasts at the Hydrometeorological Centre of Russia. The distinct feature of the SL-AV dynamical core is the semiimplicit, semi-Lagrangian vorticity-divergence formulation on the unstaggered grid. A semi-implicit, semi-Lagrangian approach allows for long time steps but violates the global and local mass conservation. In particular, the total mass in simulations with semi-Lagrangian models can drift significantly if no a posteriori mass-fixing algorithm is applied. However, the global mass-fixing algorithms degrade the local mass conservation. The new inherently mass-conservative version of the SLAV model dynamical core presented here ensures global and local mass conservation without mass-fixing algorithms. The mass conservation is achieved with the introduction of the finite-volume, semi-Lagrangian discretization for a continuity equation based on the 3-D extension of the conservative cascade semi-Lagrangian transport scheme (CCS). Numerical experiments show that the new version of the SL-AV dynamical core presented combines the accuracy and stability of the standard SL-AV dynamical core with the massconservation properties. The results of the mountain-induced Rossby-wave test and baroclinic instability test for the massconservative dynamical core are found to be in agreement with the results available in the literature.


Motivation for the research
The modern atmospheric models used for long-range forecasting or climate change modeling should treat concentrations of the greenhouse gases and certain other atmospheric constituents as prognostic variables.The mass field of such constituents is characterized by the local and global mass conservation in the absence of sources and sinks and chemical transformations.The conservation properties should be maintained by the numerical method employed, since the global mass drift can introduce biases into the model feedback to the radiative forcing and the lack of the local conservation may contaminate the physical sources and sinks of the constituents masses due to the chemical transformations.
Treatment of the atmospheric constituent concentrations as the prognostic variables is a difficulty for the semi-Lagrangian (SL) models that is well known to violate both local and global mass conservation.In particular, the total mass of the atmosphere and the mass of its constituents was found to drift significantly during the long-range integration of the SL models (see, for example, Bates et al., 1993).The global mass correction approach (e.g., Priestley, 1993) used in some SL models obviously degrades the local mass-conservation properties.
Despite the abovementioned mass conservation issues, the semi-implicit, semi-Lagrangian (SISL) treatment of the atmospheric equations is very suitable for use in general circulation models because of its computational efficiency.Attempts were made to develop advection schemes and the atmospheric equations discretizations that combine massconservation properties with the efficiency and robustness Published by Copernicus Publications on behalf of the European Geosciences Union.
V. V. Shashkin and M. A. Tolstykh: Mass-conservative SL-AV of the SL approach.Zerroukat and Allen (2012) present the 3-D inherently mass-conservative transport scheme on the sphere.CSLAM (Lauritzen et al., 2010), the locally massconservative 2-D SL scheme on the cubed sphere provides the transport computations with great multi-tracer efficiency.The approach for consistent coupling between the discrete tracer transport and continuity equations in the SISL shallowwater model is implemented by Wong et al. (2013).Lauritzen et al. (2008) developed the inherently mass-conservative, limited-area SL dynamical core for HIRLAM model using floating Lagrangian vertical levels.
This article presents the cell-integrated, mass-conservative discretization of the continuity equation in the SISL framework for the semi-Lagrangian absolute vorticity (SL-AV) global atmospheric dynamical core.Unlike Lauritzen et al. (2008), fixed vertical levels are used.We consider this research as a first step towards the hydrostatic SL dynamical core with mass-conservative and consistent tracer transport (as discussed in Wong et al., 2013), although the tracer transport problems are beyond the scope of the article.

Brief model overview
The SL-AV is a global semi-Lagrangian hydrostatic atmospheric model (Tolstykh, 2010).The model includes the dynamical core developed at the Institute of Numerical Mathematics, Russian Academy of Sciences, in cooperation with the Hydrometeorological Centre of Russia and the subgridscale physics package from ALADIN/LACE NWP model.The main feature of the SL-AV dynamical core is the finitedifference, semi-implicit semi-Lagrangian formulation on the unstaggered grid with the horizontal divergence and the vertical component of the absolute vorticity as the prognostic variables.The horizontal grid is regular latitude-longitude with the options for variable latitude resolution and the use of the reduced lat.-long.grid.In the vertical, sigma σ = p/p s (p is the pressure and p s is the surface pressure) coordinates are used.
Medium-range and seasonal forecast versions of the SL-AV are operational at the Hydrometeorological Centre of Russia.New versions of the model are being developed now.In particular, the nonhydrostatic version for the mediumrange forecast and the hydrostatic mass-conservative version for long-range forecast and climate simulations are considered.

Article structure
Section 2 presents the formulation of the inherently massconservative SISL dynamical core, beginning with the set of atmospheric governing equations (dry, adiabatic) used (Sect.2.1).The inherently mass-conservative dynamical core makes use of absolute vorticity, divergence, and thermodynamical equations approximations along with the semiimplicit system formulation and many other discretizations from the standard nonconservative SL-AV model dynamical core.The nonconservative dynamical core is reviewed in Sect.2.2.Section 2.3 describes mass-conservative discretization of the continuity equation introduced to obtain massconservative dynamical core.Section 3 presents the results of numerical experiments.
2 Inherently mass-conservative SL dynamical core formulation

Governing equations
The governing equations for the SL-AV model dynamical core in the absence of humidity are the adiabatic primitive equations written in the σ vertical coordinate as follows: -The momentum equation in the vector form (Bates et al., 1993) with the advected Coriolis term (Rochas, 1990): Since the prognostic variables are the horizontal divergence and the vertical component of the absolute vorticity, the momentum equations are used only to derive the absolute vorticity and divergence equations (see below).
-The first equation of thermodynamics and the continuity equation with the orographic terms (Ritchie and Tanguay, 1996) to make spurious orographic resonance less severe: -Lastly, the hydrostatic balance equation In the above, V = (u, v) is the horizontal velocity vector, is Earth's angular velocity vector, is Earth's angular velocity, r is the vector joining Earth's center and the given point at the surface, (..) H is the horizontal projection of a vector, is the geopotential, p s is the surface pressure, ∇ is the horizontal gradient operator, T is the temperature, R is the ideal gas constant, c p is the specific heat capacity at constant pressure, σ is the vertical velocity in the σ coordinate system, s is the surface geopotential, T is the constant reference temperature, D = div 2 (u, v) is the horizontal divergence at the σ plane, and a is Earth's radius.
The absolute vorticity equation is obtained analytically from the component form of the momentum Eq. (1): where ζ is the relative vorticity; f = 2 sin ϕ is the Coriolis parameter; and (λ, ϕ) are the longitude and the latitude, respectively.The equation for the horizontal divergence D is obtained in the discrete form in Sect.2.2.The formulation of the mass-conservative dynamical core also requires the continuity equation to be rewritten in the integral form: where δV (t) represents an arbitrary 3-D reference volume moving with the air.

Basic (nonconservative) SL-AV dynamical core formulation
The SL-AV model uses the time-stepping scheme based on SETTLS (Hortal, 2002) time approximation in combination with the semi-implicit approach and the pseudo-secondorder decentering (Temperton et al., 2001).The discrete time form of a generic equation, is as follows: where ψ is an arbitrary variable; L and N are the linear and nonlinear operators, respectively; ψ (n+1) e = 2ψ n −ψ n−1 , t is the time step; is the small decentering parameter; and the notation ψ * means the value of ψ calculated at the departure point of the upstream semi-Lagrangian trajectory.The ψ variable can be one of prognostic variables ζ , T , or ln p s ; the time-discrete forms of the corresponding equations are as follows: -The absolute vorticity equation: where the f D-type terms are treated as a product of separately calculated f and D. The Coriolis parameter is calculated analytically, f = 2 sin ϕ, where ϕ is either the known arrival point latitude ϕ n+1 or departure point latitude ϕ n * calculated via the trajectorysearching algorithm.
-The thermodynamic equation linearized around the reference temperature T : -The continuity equation (the notation (.) * 2 implies that horizontal 2-D interpolation is used to calculate departure point values): where , and N T stands for the nonlinear terms of the thermodynamic equation.
The vertical part of D n+1 3 , i.e., ∂ σ n+1 /∂σ , is contained in the time-discrete p s Eq. ( 10).It can be excluded from Eq. ( 10) by integrating it from the model top σ = σ top to the model bottom σ = 1 using the boundary conditions σ (1) = σ σ top = 0 and treating ln p s as pseudo-3-D variable constant in the vertical: A similar technique is applied to derive the expression for σ n+1 used in the energy conversion term of the thermodynamic equation ( 9).The continuity equation ( 10) is integrated from the model top to each σ level and Eq. ( 11) is then used to eliminate p s : The terms in the braces {.} are equal to the sub-integral term in the braces from Eq. ( 11).
The time-discrete equation for the horizontal divergence is obtained with the application of the horizontal divergence operator div 2 (a 1 , a 2 ) = 1 a cos ϕ ∂a 1 ∂λ + ∂a 2 cos ϕ ∂ϕ to the component form of the momentum Eq. ( 1) linearized around T and written in the time-discrete form, similar to (7): where the vector (A n u , A n v ) is the combination of known timestep n quantities from the right-hand side of the time-discrete momentum equation.
The elimination of (T n+1 , ln p n+1 s , n+1 ) in the divergence equation (13) using Eqs.( 9), ( 11), (12), and ( 14) leads to the equation for D n+1 : where D n+1 is the vector of dimension NLEV with components D n+1 k , k = 1 . . .NLEV representing the horizontal divergence at level σ k as a function of (λ, ϕ).(Note that our considerations are still analytical in horizontal.)The vector H n is a combination of known time-level n values.The matrix M of size NLEV × NLEV results from approximation of the integrals in Eqs. ( 11) and ( 14), and the notation ∇ 2 MD n+1 means that the horizontal ∇ 2 operator is applied to each component of vector MD n+1 .
To obtain the D n+1 , the problem Eq. ( 16) is reduced to NLEV horizontal Helmholtz equations using the eigenvalue transformation M = P P T (see Bates et al., 1993, for details).The 2-D Helmholtz equations are solved on the regular latitude-longitude grid using the algorithm from Tolstykh (2002).
Given the divergence, the n + 1 time step updates of other prognostic variables (ζ n+1 , T n+1 , ln p n+1 s ) can be calculated using Eqs.( 8), (9), and (11).The horizontal wind components u and v at time step n+1 are restored from known ζ n+1 and D n+1 using the algorithm from Tolstykh and Shashkin (2012).The algorithm solves the direct problem, using fourth-order finite differences in latitude and Fourier representation in longitude.
To summarize the description above, the structure of computations at the n + 1th time step in the dynamical core is as follows: 1.The coordinates of the upstream trajectories departure points are computed using (u n , v n ), (u n−1 , v n−1 ) via the algorithm from Rochas (1990).

The horizontal wind at
is reconstructed from ζ n+1 and D n+1 using the solver described in Tolstykh and Shashkin (2012).

Mass-conservative SL discretization of the continuity equation
The mass of the air contained in the elementary volume dV = a 2 cos ϕdλdϕdσ in the hydrostatic atmosphere is m = p s (λ, ϕ)dV /g.The total mass of the atmosphere is M = p s dV /g = 1 − σ top p s dS/g, where the first integral is assumed over the all atmosphere, the second integral is over the sphere, and g is gravitational acceleration.
To get the semi-implicit mass-conservative discrete equation for p s , the integral form of the continuity equation ( 6) is linearized around p ref (λ, ϕ): Following the strategy of the SL methods, the arrival cell δV t n+1 supposed to coincide with some grid cell V and the departure cell δV (t n ) = V * is then determined with the SL trajectory-searching algorithm.Given the arrival and departure cells, Eq. ( 19) is discretized in time using the same approach (7) as for the nonconservative continuity equation: The arrival cell integral of a function is treated here as the cell-averaged value of the function multiplied by the arrival cell volume.
The mass-conservation properties of the continuity equation approximation (20) depend on the scheme used for the computation of the departure volume integrals and the approximation of ∇(p ref V ) terms.As for the departure volume computations, we use 3-D extension of the conservative cascade scheme (CCS) by Nair et al. (2002).The CCS 3-D implies the approximation of the departure volume geometry by the polyhedron with the sides parallel to the coordinate planes Oλϕ, Oλσ , and Oϕσ .Following the ideology of the cascade approach, the form of the polyhedron in CCS 3-D is chosen in a way to allow for the splitting of the 3-D integration into the three consecutive 1-D integrations.Piecewise parabolic subgrid reconstruction (Colella and Woodward, 1984) without limiters and filters is used for the 1-D integral approximation.
The Much like to the nonconservative p s equation ( 10), the mass-conservative one (20) contains the (∂ σ/∂σ ) n+1 term.As in the nonconservative case, Eq. ( 20) is integrated from the model top to the model bottom using the boundary conditions σ (1) = σ σ top = 0 to eliminate the vertical velocity σ .The vertical integration in the case of Eq. ( 20) is equal to the sum over the vertical column of the arrival cells V k , k = 1 . . .NLEV spreading from the model top to the model bottom.The resulting mass conservative p s equation can be written as where S is the square of the base of the vertical column, V k * is the departure cell corresponding to the arrival cell V k , and Equations ( 20) and ( 21) are used to derive the expression for σ n+1 to be used in the energy conversion term of the thermodynamic Eq. ( 9) consistent with the mass-conservative  20) is summed over the vertical column of cells V k , k = 1 . . .K, and K = 1 . . .NLEV − 1, and Eq. ( 21) is used to eliminate p n+1 s .The resulting equation for σK+1/2 is the terms in braces {. ..} are equal to the term in the braces in Eq. ( 21).The computational procedure of the n+1th time step of the presented mass-conservative dynamical core is as follows: 1.The coordinates of the departure points of the upstream trajectories are computed.
3. The Helmholtz problem ( 16) is solved and divergence D n+1 is obtained.Note that nonconservative continuity equation ( 10) is still implicitly used in the Helmholtz equation system ( 16).
5. The horizontal wind V n+1 is reconstructed from known D n+1 and ζ n+1 using the solver described in Tolstykh and Shashkin (2012).

Numerical experiments
We test the presented mass-conservative version of the SL-AV model dynamical core (further denoted as SLAV-MC) with the mountain-induced Rossby-wave and the Jablonowski and Williamson (2006a) baroclinic instability test cases.The tests are carried out using four regular grids with 400 × 250, 640 × 400, 800 × 500, and 1200 × 750 grid points in longitude and latitude, and the corresponding horizontal grid spacings are 0.9 • × 0.72 • , 0.5625 • × 0.45 • , 0.45 • × 0.36 • , and 0.3 • × 0.24 • in longitude and latitude, respectively.In the vertical, we use the set of 28 equally spaced levels with σ top = 10 −3 .
The implicit ∇ 4 horizontal diffusion (see Tolstykh, 1997, for details) is applied for ζ , D, and T in all tests.The resolution-dependent diffusion coefficients K ζ , K D , and K T are presented in Table 1.

Mountain-induced Rossby wave
This 3-D analog of the shallow-water test case no. 5 from Williamson et al. (1992) is carried out to check the performance of the mass-conservative dynamical core in the presence of orography.The test setup presented in Jablonowski et al. (2008) is used.The initial conditions present the hydrostatically balanced smooth zonal flow, which is the stationary analytic solution to the primitive equations in the absence of the orography.Given the nonzero orography, the zonal flow breaks up and a Rossby-wave train begins its evolution.
The SLAV-MC setup for the test uses the reference surface pressure p ref = p 0 exp (− S /RT 0 ) with p 0 = 930 hPa and T 0 = 288 K (equal to the initial isothermal state of the atmosphere).This choice of p ref produces the orographic correction terms similar to Ritchie and Tanguay (1996) in the mass-conservative continuity equation ( 19) that improved the SLAV-MC accuracy near the mountain and reduced spurious orographic resonance.The reference temperature T is set to 320 K.The time step for the 400 × 250 grid simulations is 3600 s, which gives the initial zonal CFL number C ≈ 0.72.In higher-resolution simulations, the time step is chosen to keep the initial CFL number the same.The developed circulation longitudinal and meridional CFL numbers with the time steps chosen are about 3.0 and 1.8, respectively, in all simulations.
The SLAV-MC dynamical core conserves the global mass up to machine precision, whereas the standard SLAV dynamical core with the mass fixer turned off produces the monotonic global mass decrease that amounts to 0.02 % of the total atmosphere mass during the month integration of the test case initial conditions on the 640 × 400 grid.Such a mass trend cannot be considered as negligible for integration periods longer than a year.However, it does not affect the solution in a 1-month experiment.Indeed, standard and massconservative SL-AV test solutions are practically identical (the comparison is not presented); thus only the SLAV-MC solution is discussed here.
Figure 1 presents the day 25 geopotential height, temperature, and relative vorticity fields at 700 hPa from SLAV-MC simulations at the lowermost resolution (400 × 250 grid) and the finest resolution (1200 × 750 grid) used in the study.Obviously, the high-resolution solution better resolves the finer features in the vicinity of the mountain (90 • E, 30 • N), especially in the vorticity field.The large-scale structure of the solution remains the same in both high-and low-resolution runs (see the geopotential height field snapshots).The 400 × 250 grid solution is found to be in reasonable agreement with the reference solution presented in Jablonowski et al. (2008) (solution from finite-volume dynamical core of the CAM atmospheric model on a 360 × 181 regular latitude-longitude grid).

Baroclinic instability test
The test case as described in Jablonowski and Williamson (2006a) consists of two parts.The first part tests the ability of the dynamical core to maintain the steady-state initial conditions with two midlatitude jets.The second part of the test consists of the same steady-state initial conditions with overlaid zonal wind speed perturbation starting the evolution of the baroclinic wave.
The time step for 400 × 250 grid simulations is 2700 s, yielding the initial maximum zonal CFL number C ≈ 1.3.In higher-resolution simulations, the time step is chosen to keep the CFL number the same.The CFL number used is at least twice higher than that used by Jablonowski and Williamson (2006a) in simulations with the SL dynamical core of the CAM3 model.The reference state used in the experiment is the constant reference surface pressure p ref = 900 hPa and the reference temperature T = 320 K.
In the first part (stationary case) of the test, the model deviation from the initial state (which is the analytic solution) is dominated by the numerical vertical integration error in the hydrostatic balance Eq. ( 14).This gives the root-meansquare l 2 error (as defined in Jablonowski and Williamson, 2006a) in p s field of about 0.2 hPa after 30 days of integration.The error is reduced twice when increasing the number of vertical levels up to 50.The initial state is also affected by the Helmholtz solver boundary conditions near the poles, which are only second-order accurate in the latitudinal direction.This produces the p s decrease of about of 3 hPa around the poles in the simulations with the 400 × 250 regular grid.The p s field remains symmetric under all test conditions.
In the second part of the test case, zonal wind speed perturbation added to the geostrophically balanced initial conditions triggers the evolution of the baroclinic wave.The maximum developed circulation longitudinal and meridional CFL numbers with the time steps used are about 3.5 and 1.8, respectively.Both the mass-conservative and standard versions of the SL-AV model remain stable for at least 30 days of simulation with the chosen time steps and diffusion coefficients.
As in the mountain-induced Rossby-wave test, the standard non-mass-conservative SL-AV solution is characterized by the monotonic loss of the total atmosphere mass of about 0.02 % per 30 days (the 640 × 400 grid simulation), whereas the mass-conservative version of the model preserves the total mass integral up to machine precision.Although the mentioned mass loss can be a problem in a longer period of integration, it does not influence the solution in a 30 day simulation.In fact, the numerical test solutions by two versions of the model are practically identical; therefore only the SLAV-MC solution is discussed here.
Figure 2 shows the comparison of the surface pressure and 850 hPa temperature fields from the SLAV-MC simulation at day 9 using grids with various resolution.We find the basic structures of the solution well resolved even in the coarsest 400 × 250 grid run.The field patterns gradually become more developed once the resolution increases (with the most remarkable development in the temperature field).No apparent phase shift can be noticed between lowerand higher-resolution runs.One can find a good agreement between the snapshots presented in Fig. 2 and the snapshots of the reference solutions given in the Jablonowski and Williamson (2006b) and ICON dynamical core solution Wan et al. (2013).
As compared to the reference solutions presented in Jablonowski and Williamson (2006b) -namely solutions from Eulerian, semi-Lagrangian, and finite-volume dynamical cores of the CAM atmospheric model (Collins et al., 2004) and the dynamical core of GME atmospheric model (Majewski et al., 2002) -the 850 hPa relative vorticity field at day 9 from the SLAV-MC solution presented at Fig. 3 looks generally smoother, and no Gibbs phenomenon as in the CAM-EUL and CAM-SLD can be observed.The shape and magnitude of the field features agree well with the reference solutions and also ICON-HDC (Wan et al., 2013) solution.One can note that there are no negative vorticity values inside the vortices in the SLAV-MC solution.
The quantitative assessment of the similarity and difference between the SLAV-MC solution and reference solution from the T340 spectral SL dynamical core of the CAM3 atmospheric model is available via the l 1 , l 2 , and l ∞ surface pressure difference norms defined in Jablonowski and Williamson (2006a).The upper row in Fig. 4 presents difference norms plots for the SLAV-MC solutions at various resolution.The lower row of the figure shows the difference norms between the lower-resolution SLAV-MC solutions and the highest resolution (1200 × 750 long.-lat.grid) SLAV-MC solution.The gray shading in Fig. 4 denotes the uncertainty of the numerical solution obtained in Jablonowski and Williamson (2006a) by comparing different reference solutions.
The difference norms shown on the upper row of Fig. 4 are all below uncertainty limit; this confirms the convergence of the SLAV-MC simulations to reference solution.Moreover, one can see that the lower-resolution SLAV-MC solutions converge to the highest resolution SLAV-MC solution (at 1200 × 750 grid) since the corresponding difference norms (lower panel of Fig. 4) also fall below the uncertainty.Finally, the difference norms between standard SLAV and SLAV-MC solutions of equal resolution are well below the uncertainty limit (not shown) that proves the similar behavior of the two versions of the dynamical core.

Conclusions
A semi-implicit time integration scheme in conjunction with the semi-Lagrangian treatment of advection allows running of atmospheric simulations with time steps larger than time steps limited by CFL stability condition and thus building of computationally efficient models.Indeed, it was shown that the semi-Lagrangian advection can be implemented efficiently on massively parallel computer systems using up to O(10 4 ) processors (White and Dongarra, 2011).Recently, it was found that the elliptic solver necessary to implement the semi-implicit scheme can also use such systems efficiently (see Müller and Sheichl, 2013).However, the application of SISL methods in modern atmospheric models used for climate simulations is limited by the absence of inherent mass conservation requiring a global mass fixer.
We have presented here a version of SISL dynamical core for the SL-AV global model that is inherently mass conservative without use of mass correctors.The mass conservation is achieved by the introduction of the cell-integrated semi-Lagrangian discretization for the continuity equation.This discretization is based on the 3-D extension of the conservative cascade SL transport scheme (CCS) by Nair et al. (2002).Except for the new discretization of the continuity equation, approximation of the primitive equations and the semi-implicit equation system formulation in the massconservative version are the same as in the standard version, and therefore only minimal changes to the dynamical core are introduced.
The numerical experiments showed that the inherently mass-conservative version of the SL-AV dynamical core (SLAV-MC) is as accurate and stable with long time steps as the standard nonconservative version of this dynamical core.The results of SLAV-MC for the baroclinic instability test (Jablonowski and Williamson, 2006a) and mountain-induced Rossby-wave test (from Jablonowski et al., 2008) are found to be in agreement with the results available in the literature.In the baroclinic instability test case, the difference norms between SLAV-MC solutions in various resolution and the reference T340 SL solution are below the solution uncertainty calculated in Jablonowski and Williamson (2006a).The behavior of two versions of the dynamical core in the numerical www.geosci-model-dev.net/7/407/2014/Geosci.Model Dev., 7, 407-417, 2014 experiments is very similar, except that the standard version (with mass corrector turned off) produces the monotonic loss of global mass, which can be crucial in the longer-period simulations.SLAV-MC conserves the global mass up to machine precision.
The presented approach efficiently combines the advantages of the SISL method with the inherent mass conservation.Thus we believe that our research can be the basis for building an SISL dynamical core of an atmospheric general circulation model suitable for long-range forecasting and climate simulations.In particular, we plan to implement the hybrid σ -p vertical coordinate and the reduced lat.-long.grid (Fadeev, 2013), as we did for the shallow-water model (Tolstykh and Shashkin, 2012).Furthermore, consistent transport formulation similar to Wong et al. (2013) is considered.

p
s = p s − p ref .On the right-hand side of this equation, we have used the Eulerian treatment of the d dt and the fact that ∂p ref ∂t = 0 and ∂p ref ∂σ = 0.
type terms are calculated with the 2-D divergence calculation algorithm from the mass-conservative shallow-water model Tolstykh and Shashkin (2012).The algorithm used guarantees ∇ (p ref V ) dS = 0 (the integral is over the sphere) and thus in combination with CCS 3-D ensures the mass conservation of the continuity equation approximation (20) (see Tolstykh and Shashkin, 2012, for a detailed discussion).

Fig. 2 .
Fig. 2. The day 9 p s (left column) and T at 850 hPa surface (right column) in the Jablonowski and Williamson (2006a) test case simulated by the SLAV-MC dynamical core at the various resolution grids.

Fig. 3 .
Fig. 3.The day 9 850 hPa relative vorticity field in the Jablonowski and Williamson (2006a) test case simulated by the SLAV-MC dynamical core at the various resolution grids.

Fig. 4 .
Fig.4.Time evolution of the l 1 , l 2 , and l ∞ p s difference norms.Upper row: difference between the spectral T340 SL reference solution and SLAV-MC solutions at 400 × 250 (red line), 640 × 400 (green line), 800 × 500 (blue line), and 1200 × 750 (orange line) regular lat.-long.grids.Lower row: difference between the 1200 × 750 SLAV-MC solution and the lower-resolution SLAV-MC solutions (the line colors are the same as on the upper row).Gray shading presents the uncertainty of reference solutions.