Toward Unification of the Multiscale Modeling of the Atmosphere

As far as the representation of deep moist con-vection is concerned, only two kinds of model physics are used at present: highly parameterized as in the conventional general circulation models (GCMs) and explicitly simulated as in the cloud-resolving models (CRMs). Ideally, these two kinds of model physics should be unified so that a continuous transition of model physics from one kind to the other takes place as the resolution changes. With such unification, the GCM can converge to a global CRM (GCRM) as the grid size is refined. This paper suggests two possible routes to achieve the unification. ROUTE I continues to follow the pa-rameterization approach, but uses a unified parameterization that is applicable to any horizontal resolutions between those typically used by GCMs and CRMs. It is shown that a key to construct such a unified parameterization is to eliminate the assumption of small fractional area covered by convective clouds, which is commonly used in the conventional cumu-lus parameterizations either explicitly or implicitly. A preliminary design of the unified parameterization is presented, which demonstrates that such an assumption can be eliminated through a relatively minor modification of the existing mass-flux based parameterizations. Partial evaluations of the unified parameterization are also presented. ROUTE II follows the " multi-scale modeling framework (MMF) " approach , which takes advantage of explicit representation of deep moist convection and associated cloud-scale processes by CRMs. The Quasi-3-D (Q3-D) MMF is an attempt to broaden the applicability of MMF without necessarily using a fully three-dimensional CRM. This is accomplished using a network of cloud-resolving grids with large gaps. An outline of the Q3-D algorithm and highlights of preliminary results are reviewed.


Introduction
In his landmark paper published near the beginning of the 20th century, Vilhelm Bjerknes (Fig. 1) pointed out that necessary and sufficient conditions for the rational solution of the forecasting problem are the following (Bjerknes, 1904): 1.One has to know with sufficient accuracy the state of the atmosphere at a given time; 2. One has to know with sufficient accuracy the laws according to which one state of the atmosphere develops from another.
These conditions may sound obvious.In pointing out 2., however, he distinguished the laws for changes "from degree to degree in meridian and from hour to hour in time" from those for changes "from millimeter to millimeter and from second to second."Numerical modeling of the atmosphere Published by Copernicus Publications on behalf of the European Geosciences Union.
has been a struggle for finding such laws, and still is.Bjerknes specifically expressed his concern on the representation of frictional stress, which "depends on very small-scale velocity differences."We see that similar problems exist for all microphysical processes.The progress of our ability to represent cloud microphysical processes in climate models has been especially slow (Randall et al., 2003).This paper is based on the Vilhelm Bjerknes Medal Lecture presented by the first author at the 2010 EGU General Assembly in Vienna.The rationale for the main theme of the paper is the following assessment of the current status of numerical modeling of the atmosphere: As far as the representation of deep moist convection is concerned, only two kinds of model physics are used at present: highly parameterized and explicitly simulated.Correspondingly, besides those models that explicitly simulate turbulence such as Direct Numerical Simulation (DNS) and Large Eddy Simulation (LES) models, we have two discrete families of atmospheric models as shown in Fig. 2: one is represented by the conventional general circulation models (GCMs) and the other by the cloud-resolving models (CRMs).In this figure, the abscissa is the horizontal resolution and the ordinate is a measure for the degree of parameterization, such as the reduction in the degrees of freedom, increasing downwards.These two families of models have been developed with applications to quite different ranges of horizontal resolution in mind.
Naturally, there have been a number of studies examining a broader applicability of each family as shown by the horizontal arrows in Fig. 2: applicability of CRMs to lower resolutions and that of GCMs to higher resolutions.Weisman et al. (1997), for example, examined the applicability of a CRM to squall-line simulations for midlatitude-type environments with resolutions between 1 and 12 km, and concluded that with resolutions coarser than 4 km, the evolution is characteristically slower and the resultant mature mesoscale circulation is stronger.
With respect to the applicability of GCMs to higher resolutions, the work of Williamson (1999) is particularly intriguing.The paper shows that, when the horizontal resolution of the NCAR CCM2 is increased for both the dynamics and physical parameterizations, the upward branch of the Hadley circulations increases in strength and there is no sign of convergence.When the horizontal resolution is increased for the dynamics but not for the parameterizations, the solution converges.But the converged state is similar to that obtained with the coarse resolution for both so that the increased resolution for the dynamics is wasted.Together with other evidence, he concludes, "the results raise a serious questionare the parameterizations correctly formulated in the model ?[...] The parameterizations should explicitly take into account the scale of the grid on which it is based."A similar question on parameterization was raised by Skamarock and Klemp (1993) in the context of adaptive grid refinement.Also, analyzing the impact of horizontal resolution increases on the error growth of ECMWF forecasts, Buizza (2010) sug- gested that rather than resolution, it is model improvements that might lead to better predictions and longer predictability limits.
Strictly speaking, truncation of a continuous system can be justified only when the resulting error can be made arbitrarily small by using a higher resolution.Our problem is, therefore, more demanding than just a convergence; the GCM should converge to a physically meaningful high-resolution model such as a CRM applied to the global atmosphere.This requires that both the dynamics and physics of GCMs formally converge to those of the CRM as shown schematically by the dashed curve in Fig. 2. If the GCM and CRM share the same dynamics core, which must necessarily be nonhydrostatic, we expect that the convergence is not an issue as far as the model dynamics is concerned.Unfortunately, the same is not true for the conventional formulations of model physics, especially when cloud processes are involved.
Figure 3 schematically illustrates the qualitative difference of model physics between the two families of models.For a given observed large-scale condition, we can identify the apparent heat source, Q 1 , and the apparent moisture sink, Q 2 , from the residuals in large-scale heat and moisture budgets (e.g., Yanai et al. 1973).Here the heat source and moisture sink refer to the source of the sensible heat, c p T , and the sink of the latent heat, Lq, respectively.The left panel of Fig. 3 schematically shows typical profiles of Q 1 , Q 2 , and Q 1 − Q 2 for disturbed tropical conditions.The difference Q 1 −Q 2 gives the apparent moist static energy source.Here the moist static energy is defined by c p T + Lq + gz, where gz is the geopotential energy.As shown in the figure, the profile of Q 1 − Q 2 typically has negative values in the lower troposphere and positive values in the middle to upper troposphere, suggesting the dominant role of vertical eddy transport of moist static energy.Since GCMs must produce this type of profiles when deep moist convection is dominant, we call this type the "GCM-type".
In contrast, the local cloud microphysical processes produce practically no moist static energy source/sink except near the freezing level.This is because moist static energy is conserved under moist-adiabatic processes and thus there is no significant source/sink of moist static energy except where the ice phase is involved.Within updrafts and downdrafts/precipitation, there are sources immediately above the freezing level due to freezing and sinks immediately below that level due to melting, respectively.This is illustrated in the right panel of Fig. 3. Cloud microphysics in CRMs is expected to produce this type of profiles, which we call the "CRM-type".
As Arakawa (2004) emphasized, it is important to recognize that any space/time/ensemble average of the CRM-type profiles does not give a GCM-type profile.This means that the cumulus parameterization problem is more than a statistical theory of cloud microphysics.Also, it is not a purely physical/dynamical problem because it is needed as a consequence of mathematical truncation.Finally, it is not a purely mathematical problem since the use of a higher resolution or an improved numerical method does not automatically improve the result.A complete theory for cumulus parameterization must address all of these aspects in a consistent manner including the transition between the GCM-type and CRM-type profiles.
To see the resolution dependence of model physics, Jung and Arakawa (2004) performed budget analyses of data simulated by a CRM with different space/time resolutions and with and without the cloud-microphysics component of the CRM.A highlight of their results is summarized in Fig. 4, which shows ensemble-averaged profiles of the "required" moist static energy source.Here, the required source is defined as the source needed for the result of the low-resolution run without cloud microphysics to agree with that of highresolution run with cloud microphysics.The space/time intervals used in the analysis are (32 km, 60 min) for the green profile, (2 km, 2 min) for the red profile, and in-between for other profiles.We note that the green and red profiles are of the GCM-and CRM-type, respectively.The existence of the in-between profiles strongly suggests that the required model physics smoothly transitions from one type to the other as the resolution changes.If the model physics of GCMs is reformulated to produce such transition, GCMs and CRMs are unified to a single family of models that can cover a wide range of horizontal resolutions.
We can think of two routes to achieve the unification of the two families of models.ROUTE I continues to follow the parameterization approach, but uses a unified parameterization with which the GCM converges to a global CRM (GCRM) as the grid size is refined.ROUTE II, on the other hand, replaces the parameterization of deep moist convection with a partial simulation of cloud-scale processes by CRMs and formulates the coupling of GCM and CRMs in such a way that the coupled system formally converges to a GCRM as the GCM grid size is refined.Thus ROUTE I and ROUTE II share the same destination.On the other hand, they have different kinds of departure points corresponding to the two kinds of model physics pointed out above, one of which uses the physics of GCMs and the other that of CRMs.Naturally there can be many variations within each of ROUTE I and ROUTE II and, therefore, they represent two families of routes.
The rest of this paper is organized as follows: Sect. 2 discusses ROUTE I, which follows the unified parameterization approach, while Sect. 3 discusses ROUTE II, which follows the coupled GCM/CRM approach.Finally, Sect. 4 presents summary and further discussions.

Identification of the problem
To identify the problem to be addressed in ROUTE I, we use the Arakawa and Schubert (1974, AS hereafter) parameterization as an example of the starting point.AS says, "Consider a horizontal area -large enough to contain an ensemble of cumulus clouds but small enough to cover a fraction of a large-scale disturbance.The existence of such an area is one of the basic assumptions of this paper."In reality, the GCM grid cells are not large enough and, at the same time, not small enough.
A major part of this problem can be solved if the parameterization is modified in such a way that it can be applied to a wide range of resolution including that of CRMs.AS and many other mass-flux based parameterizations (e.g., Tiedtke, 1989;Kain and Fritsch, 1990;Emanuel, 1991;Gregory and Rowntree, 1990;Zhang and McFarlane, 1995) exclude this possibility by assuming σ 1, either explicitly or implicitly, where σ is the fractional area covered by all convective clouds in the grid cell.With this assumption, the temperature and water vapor to be predicted are essentially those for the cloud environment.Then, as illustrated by the red arrows in Fig. 5a, relevant physical processes are "cumulusinduced" subsidence in the environment and detrainment of cloud air into the environment.(Here it is important to note that the "cumulus-induced subsidence" is only a hypothetical subsidence.It is a component of the subgrid-scale eddy, which has its own mass budget closed within the same grid cell.This does not mean that the true subsidence is confined within that cell.The true subsidence is the sum of the blue and red vertical arrows in Fig. 5a, which normally tend to compensate each other.In such a case, the true subsidence occurs in other grid cells, which may well be far away, whose position is determined by the grid-scale dynamics, not by the parameterization.This point is often misunderstood in the literature.) As the grid size becomes smaller, however, the cloud may eventually occupy the entire grid cell so that there is no "environment" within the same cell.As Fig. 5b indicates, the probability density distribution of σ becomes bimodal in this limit, consisting of σ = 1 and σ = 0 (Krueger, 2002;Krueger and Luo, 2004).It is then clear that a key to open ROUTE I is to include a transition to this limit by eliminating the assumption of σ 1.
To visualize the problem raised above, we have analyzed datasets simulated by a CRM as Krueger (2002) and Krueger and Luo (2004) did.The simulations are performed by applying the 3-D anelastic vorticity equation model of Jung and Arakawa (2008) to a horizontal domain of 512 km×512 km with a 2 km grid size.Other experimental settings are similar to the benchmark simulations used by Jung and Arakawa (2010).Two 24-h simulations are made, one with and the other without background shear.Figure 6 shows snapshots of vertical velocity w at 3km height at the end of these simulations.As is clear from the figure, the two simulations represent quite different cloud regimes.For the analysis presented in the rest of this section, datasets are taken from the last 2-hour period of each simulation with 20-min intervals.
To analyze the resolution dependency of the statistics of data, we divide the original 512 km domain into sub-domains of equal size.The selected side lengths of the sub-domains are d n = 512 km/2 n−1 , n = 2,3,4,...,9.The original domain can then be identified by n = 1. Figure 7 shows the original and examples of the sub-domains.In the analysis presented here, grid points that satisfy w>0.5 m/s are considered as "cloud points".This simple criterion for "cloud points" is used only for the diagnosis and the initial evaluations, not in the formulation of the unified parameterization presented in Sect.2.2.
Figure 8 shows σ at 3 km height averaged over all cloudcontaining (i.e., σ =0) sub-domains against the sub-domain size d n for the shear case (a) and non-shear case (b).The associated standard deviation is also shown.It is clear that σ 1 can be a good approximation only when coarse resolutions, say, d n ≥32 km, are used.The average σ tends to increase as d n decreases and becomes 1 for d n = 2 km, which is the grid size of the CRM used.The standard deviation is very large for high resolutions, but it is expected since there is no reason to believe that σ is a unique function of d n .In spite of the large standard deviation, the tendency toward a bimodal distribution consisting of σ = 1 and σ = 0 can be seen for high resolutions.We see that the average σ is generally larger for the shear case than for the non-shear case.The number of cloud-containing sub-domains is, however, smaller for the shear case (not shown).

Expressions for vertical eddy transport
Recall that the vertical eddy transport of moist static energy is responsible for the difference between the GCM-type and CRM-type profiles illustrated in Fig. 3.We now attempt to parameterize the vertical eddy transport of a thermodynamic variable ψ in a way applicable to any value of σ including 1.For this, we require that the parameterized eddy transport decrease as the fractional convective cloudiness σ increases and it eventually vanish at σ = 1.Such decrease of the parameterized component gives more room for explicitly simulated non-eddy transport that can directly interact with the parameterized eddy transport and other processes.The objective of this subsection is to seek a simple but reasonable possibility for parameterization that satisfies this requirement.If we do not see such a possibility, Route I may not be worth exploring because we can easily get lost in the middle of the route due to the complexities involved.
As is commonly done in the conventional cumulus parameterizations, we first assume that the cloud and environment values of ψ, denoted by ψ c and ψ, respectively, are horizontally uniform individually.The existence of convective-scale downdrafts is ignored at this stage.
Let an overbar denote the mean over the entire area of a grid cell.Then, using the assumption stated above, (1) Further let w c and w be the averages of w over the clouds and the environment, respectively.Then, and From Eq. ( 3) with Eqs. ( 1) and ( 2), we find The left hand side of Eq. ( 4) is the vertical eddy transport of ψ per unit horizontal area and density, which we simply call the "eddy transport" of ψ.As expected, the eddy transport vanishes for σ = 0 and σ = 1.Equation ( 4) is, however, not convenient to use due to the appearance of the environmental values, w and ψ, which are not well defined when σ ∼1.We thus eliminate these variables using Eqs.( 1) and (2).After some manipulations, we find If we assume σ 1, | w| w c and ψ≈ψ in Eq. ( 4), we obtain where M c ≡ρσ w c is the cloud mass flux.Equation (6) can also be derived from Eq. ( 5) assuming σ 1 and |w| w c .Most conventional mass-flux based parameterizations effectively use this expression [see, for example, Eqs. ( 35) and (36) of AS], which indicates that the cloud vertical velocity w c matters only through M c as far as the eddy transport is concerned.This is not the case if the assumption σ 1 is eliminated.The unified parameterization must determine both σ and w c .We now introduce the requirement that the parameterization converge to an explicit simulation of cloud processes as σ → 1.This means that lim σ →1 w c = w and lim It follows that both w c − w and ψ c − ψ are the order of 1 − σ (or higher) so that (w c − w) ψ c − ψ is the order of (1 − σ ) 2 (or higher) near this limit.The simplest choice to satisfy this requirement is where an asterisk denotes a limiting form expected when σ 1, typically obtained by a non-interacting plume model.
It is important to keep in mind that (w c − w) ψ c − ψ * does not depend on σ though it generally depends on the res-olution, which can be anywhere between the GCM and CRM resolutions.Substituting Eq. ( 8) into Eq.( 5), we obtain This represents the basic structure of the unified parameterization.Since (w c − w) ψ c − ψ * does not depend on σ , Eq. ( 9) shows that the dependence of wψ − w ψ on σ is through the factor σ (1 − σ ), whose maximum is 0.25 for σ = 0.5.To close the unified parameterization, we must determine σ and (w c − w) ψ c − ψ * .It should be noted that what we have parameterized in this subsection is the dependency of the eddy transport on σ , not on the resolution.As pointed out in Sect.2.1, there is no reason to believe that σ is a unique function of the resolution.For the dataset we used, a coarse resolution indicates a small value of σ as Fig. 8 shows.But a small value of σ does not necessarily indicate a coarse resolution.Judging from the large error bars in Fig. 8 for the range of resolutions 8∼32 km, we expect that a small value of σ may appear for these resolutions also.Besides, the number of sub-domains with σ = 0 tends to increase as the resolution increases.Furthermore, a large value of σ including σ = 1 may well appear with low resolutions (not in this dataset) when stratiform clouds are dominant.

Partial evaluations of the unified parameterization
Before proceeding to the closure problem, we present in this subsection partial evaluations of the unified parameterization, again using the dataset mentioned in Sect.2.1.The evaluation is designed to check the formal structure of Eq. ( 9) using water vapor mixing ratio q as an example of ψ.This equation with ψ = q becomes wq − w q = σ (1 − σ )[(w c − w)(q c − q)] * . (10) Since the solution of the unified parameterization can converge to only one of many realizations that can be simulated by the CRM, the dependence of the eddy transport on σ must inevitably be statistical especially with high resolutions.
Even when high resolution is used, our primary interest may still be in the effects of cloud processes on larger scales.We thus include evaluations of Eq. ( 10) applied to the fields of wq −w q and σ (1 − σ ) averaged over large horizontal areas.
Table 1 shows the correlation coefficients between wq − w q and σ (1 − σ ) for various combinations of sub-domain and averaging sizes for the shear case (a) and the non-shear case (b).The sub-domain size is necessarily less than or equal to the averaging size.We see that, except when subdomain and averaging sizes are both small, the correlations are almost uniformly high for both cases, suggesting that the σ -dependence of the eddy transport shown in Eq. ( 10) is basically valid.For the example marked in bold font in Table 1a, the dots in Fig. 9 show the σ -dependence of the diagnosed wq − w q averaged over all sub-domains that share σ in the same sub-range.(Here datasets are taken from the last 12-h period of the simulation to increase the sample size.)The dashed curve in Fig. 9, on the other hand, shows the right hand side of Eq. ( 10) when [(w c − w)(q c − q)] * is chosen for the best fit to the dots for small and medium values of σ .The excellent fit of the two lines again indicates that the simple choice we have made for the σ -dependence of the eddy transport is basically valid.For high values of σ , however, the diagnosed values are significantly larger than the expected values shown by the dashed curve.Also, as one may expect from the relatively low correlation coefficient 0.62 presented in Table 1a, the standard deviations associated with these averaged values are very large with a magnitude roughly equal to the mean value itself.These two aspects of the results suggest interesting physical/statistical problems to be considered in future, which are pointed out in Sect.2.5.
While the σ -dependence in Eq. ( 10) seems justifiable except for large values of σ , it remains to be seen whether this equation can correctly produce the dependence of the eddy transport on the resolution.The σ -dependency in Eq. ( 10) is explicit, but the resolution-dependence is not.We can test the resolution-dependence by independently diagnosing the left and right hand sides from the dataset, for each sub-domain size, and then comparing the results to see if the resolutiondependence of the right hand side agrees with that of the left hand side.Following this approach, another test is performed focusing on the situations with large values of σ .To allow a compact presentation of the results, we define a weighted ensemble average X by the weighted average of X over all sub-domains of the same size with the weight σ , where X is a variable defined for each sub-domain.We use the weight σ to emphasize the contributions from the sub-domains with large values of σ .The weighted ensemble average of Eq. ( 10) is given by (11) Fig. 9.The σ -dependence of the diagnosed (wq − w q) averaged over all sub-domains that share σ in the same sub-range for the case in bold face in Table 1a.The dashed curve shows σ (1−σ ) multiplied by a constant chosen for the best fit to the dots for small and medium values of σ .
In deriving this result, we have ignored the difference of [(w c − w)(q c − q)] * between the sub-domains of the same size because it is independent of σ .This quantity depends on the cloud model that we choose for the idealized case σ 1 and, therefore, it cannot be directly diagnosed from the dataset.We thus eliminate this quantity by using the weighted ensemble average of Eq. ( 8) with ψ = q given by From Eqs. ( 11) and ( 12), we obtain The open circles in Fig. 10 show the values of wq − w q at a height of 3 km as estimated from the dataset using the right hand side of Eq. ( 13) for each sub-domain size d n .These values have peaks in the mesoscale range because they are small for coarse resolutions with σ 1 and for high resolutions with σ ∼1.For comparison, the values of wq − w q directly diagnosed from the dataset are shown by the solid circles.Amazingly, the resolution dependence of the estimated values is very similar to that of the directly diagnosed values.The magnitudes of the former are, however, systematically smaller than those of the latter.This is not surprising in view of the various idealizations used in the derivation of Eqs. ( 8) and ( 9), such as neglecting convective downdrafts and possible coexistence of different types of clouds and different phases of cloud development.Moreover, the criterion w>0.5 m/s we have adopted for cloud points influences the estimated values but not the directly diagnosed values.In any case, the results presented in this subsection provide evidence that the formal structure of the unified parameterization is basically valid even from the resolution-dependence point of view.13).Closed circles: directly diagnosed with the left hand side of Eq. ( 13).

Determination of σ and relaxed adjustment
The closure of conventional cumulus parameterizations determines the apparent source of thermodynamic prognostic variables.For ψ, it is given by S ψ − ∂ρ wψ − w ψ /ρ∂z, where S ψ is the true source of ψ per unit mass due to subgrid cloud processes and ρ is the density.From this together with S ψ , also determined by the parameterization, we can calculate the eddy transport wψ − w ψ .Let the value of wψ − w ψ calculated in this way for a full adjustment to a quasi-equilibrium state be wψ − w ψ adj .For wψ − w ψ adj to be consistent with (w c − w) ψ c − ψ * in view of Eq. ( 5), Then, for given wψ − w ψ adj and (w c − w) ψ c − ψ * , σ must be such that We see that the condition 0≤σ ≤1 is automatically satisfied by Eq. ( 15), as long as wψ − w ψ adj and (w c − w) ψ c − ψ * have the same sign, with σ → 0 as wψ − w ψ adj → 0 and σ → 1 as wψ − w ψ adj → ∞.
The unified parameterization uses σ determined in this way with a selected thermodynamical variable for ψ.Since our objective is to determine σ for use in formulating the eddy transport, it is good to select a quasi-conservative variable for ψ, for which S ψ is small, such as the moist static energy used in Figs. 3 and 4.
As far as the basic reasoning is concerned, this approach of determining σ is in parallel to Emanuel (1991) in the sense that the following two informations are combined: vertical profiles of cloud properties determined by a cloud model and the total vertical transport necessary for the adjustment to a quasi-equilibrium.
Finally, elimination of (w c − w) ψ c − ψ * between Eqs. ( 9) and ( 14) gives This is the equation we use in practical applications.
As expected, the right hand side of Eq. ( 16) becomes wψ − w ψ adj when σ 1.When σ is finite, however, Eq. ( 16) shows that the actual eddy transport is less than wψ − w ψ adj due to the factor (1 − σ ) 2 , indicating that the adjustment should be applied only partially when σ is large.This is true regardless of the resolution.The situation is then similar to the relaxed Arakawa-Schubert parameterization (Moorthi and Suarez, 1992) though in our case the degree of relaxation depends on σ .

Anticipated impact of the unified parameterization and remaining problems
In this section, we have presented a new framework for cumulus parameterization in which the result of parameterization converges to an explicit simulation of cloud processes as the resolution increases.Thus the new parameterization unifies parameterizations in GCMs and those in CRMs as far as the representation of deep moist convection is concerned.the unified parameterization, the error of the GCM in satisfying the CRM equations can be made arbitrarily small by using a higher resolution.In this way, multiscale numerical methods such as multiply-nested grids and adaptive mesh refinement (AMR) methods, for example, can be used with no problem in model physics.We emphasize that this drastic broadening of the applicability of cumulus parameterization can be achieved by a relatively minor modification of the conventional mass-flux based parameterizations.However, a good cloud model to determine (w c − w) ψ c − ψ * and a reasonable closure to determine the magnitude of wψ − w ψ adj are prerequisites to the success of determining σ and, therefore, the success of the entire unified parameterization.We can start with the methods used in conventional cumulus parameterizations to determine these quantities.But it is very much likely that these methods need to be revised as we accumulate experience with the unified parameterization.Also, the dynamics and formulations of cloud microphysics, turbulence and radiation in the model for use in explicit simulations must be such that they are applicable to a wide range of resolutions.
There are two other problems to be considered.First, as Fig. 9 indicates, the eddy transport can be larger than the value expected from the parameterization.We interpret this error of the parameterization as a consequence of the assumption that the cloud and environment values of thermodynamic variables are horizontally uniform individually.With this assumption, the parameterization ignores the eddy transport due to the coexistence of different types of clouds and different phases of cloud development and that due to the internal structures of clouds and the environment.Since the error appears in the range of large σ , it is likely that the neglect of internal structure of huge cumuliform or stratiform clouds is responsible for that error.Second, as indicated by the large standard deviations associated with the averaged eddy transports in Fig. 9, the results of the parameterization have significant uncertainties.There can be two main sources for this: one is the non-deterministic nature of the closure, as discussed by Xu and Arakawa (1992), and the other is estimating cloud properties with huge dimensions by a simple cloud model.This may suggest the necessity of including a stochastic component in the parameterization.We note, however, both of these problems may not be so serious in the unified parameterization because explicit simulations produce the dominant part of transport for large values of σ which automatically generate non-eddy transport and random processes.
When it is successfully implemented, the practical merits of the unified parameterization will be great.But we should remember that it has a limit as a "parameterization", which requires a number of idealizations to reduce the degrees of freedom.When sufficient computer resources are available, therefore, we should pursue the other approach, ROUTE II, for more realistic numerical weather prediction and climate simulations as discussed in the next section.

Route II: quasi-3-D multiscale modeling framework
We have developed Quasi-3-D Multi-scale Modeling Framework (Q3-D MMF) following ROUTE II for unification of GCMs and CRMs.It is a complement of ROUTE I, simulating the details of cloud processes at least partially.The Q3-D MMF is described in detail by Jung and Arakawa (2010) so that this paper gives only a brief outline of the framework and some highlights of its preliminary results.
MMF recognizes that we currently have two kinds of model physics, the GCM type and the CRM type.Correspondingly, MMF uses two grid systems, one for the GCM and the other for the CRM.Model physics is almost entirely determined by the statistics of CRM solutions, replacing the conventional parameterizations in GCMs.In contrast to many other multi-scale numerical methods, MMF gains computing efficiency by sacrificing full representation of cloud-scale 3-D processes.This is motivated by the fact that 2-D CRMs are reasonably successful in simulating thermodynamical effects of deep moist convection.The prototype MMF is called "Cloud Resolving Convective Parameterization" (Grabowski and Smolarkiewicz 1999;Grabowski 2001) or "Super Parameterization" (Khairoutdinov and Randall 2001;Randall et al. 2003).It replaces the cloud parameterization by a 2-D CRM embedded in each GCM grid cell as shown in Fig. 11a for a portion of the horizontal domain.The MMF is still called "parameterization" because it inherits the structure of conventional GCM; i.e., the CRM is forced by the GCM and the GCM recognizes only the domain-averaged values of the CRM results.
The Quasi-3-D (Q3-D) MMF we have developed is an attempt to broaden the applicability of the prototype MMF without necessarily using a fully three-dimensional CRM.The horizontal domain of the Q3-D MMF consists of a network of perpendicular sets of channels, each of which contains grid-point arrays as shown in Fig. 11b.The grey areas in the figure represent the gaps of the network.For computing efficiency, the gaps are chosen to be large by using a narrow width for the channels, barely enough to cover a typical cloud size in the lateral direction.Thus, a channel may contain only a few grid-point arrays, whose minimum number required for resolving local 3-D processes is two as in Fig. 11b.
Because the channels are so narrow, it is crucial to select a proper lateral boundary condition to realistically simulate the statistics of cloud and cloud-associated processes.Among the various possibilities, a periodic lateral boundary condition is chosen for the deviation from a background field obtained through interpolation from GCM grid points.We design the coupling of the two grid systems in such a way that the deviation vanishes as the GCM grid size approaches that of the CRM.Thus the whole system of the Q3-D MMF formally converge to a fully 3-D global CRM as schematically shown in Fig. 12. Consequently, the horizontal resolution of the GCM can be freely chosen depending on the objective of application without changing the formulation of model physics.For more details of the Q3-D algorithm, see Jung and Arakawa (2010).
To evaluate the Q3-D CRM in an efficient way, idealized experiments are performed using a small horizontal domain.First, benchmark simulations are made using a fully 3-D CRM.Then a Q3-D simulation is made for the situation corresponding to each of the benchmark simulations.The grid used in these tests is similar to that shown in the central GCM grid cell of Fig. 11b, consisting of only one pair of perpendicular channels with only two grid points across each channel.Since the horizontal domain is so small, the GCM component is made inactive in these tests.Thus the GCM grid point values are taken from the benchmark simulations after horizontal smoothing.These values are then interpo-  Figure 13 shows time series of precipitation, evaporation and sensible heat flux at the surface.The Q3-D results fluctuate more than those of the BM because the sample size of the former is much smaller than that of the latter.The

Summary and conclusion
As far as the representation of deep moist convection is concerned, conventional GCMs and CRMs have quite different formulations of model physics, each of which is applicable to only a limited range of horizontal resolution.These two kinds of model physics should be unified so that a continuous transition from one kind to the other takes place as the resolution changes.Then a resolution between those typically used by GCMs and CRMs can be freely chosen depending on the objective of the application.
This paper suggests two possible routes to achieve the unification.ROUTE I uses a new framework for cumulus parameterization in which the result of parameterization converges to an explicit simulation of cloud processes as the resolution increases.In this way the framework unifies parameterizations in GCMs and CRMs as far as the representation of deep moist convection is concerned.With the unified parameterization, the error of the GCM solution measured by the difference from a CRM solution can be made arbitrarily small by using a higher resolution.It is shown that a key to construct a unified parameterization is to eliminate the assumption of small fractional area covered by convective clouds, which is commonly assumed in the conventional cumulus parameterizations either explicitly or implicitly.A preliminary design of the unified parameterization is presented, which demonstrates that such an assumption can be eliminated through a relatively minor modification of the existing mass-flux based parameterizations.Partial evaluations of the unified parameterization are also presented.
It is pointed out that the unified parameterization has a limit as a "parameterization" even when it is successful.Route II, which follows the MMF approach, complements Route I through the use of CRM physics throughout.The Quasi-3-D (Q3-D) MMF is an attempt to broaden the applicability of the prototype MMF without necessarily using a fully three-dimensional CRM.A great advantage of the Q3-D MMF is that it converges to a 3-D CRM as the GCM's resolution is refined while maintaining the same CRM physics.An outline of the Q3-D algorithm and highlights of preliminary results are presented.Comparing the simulation results with the corresponding benchmark simulation performed with a 3-D CRM, it is concluded that the Q3-D CRM can reproduce most of the important statistics of the 3-D solutions, including precipitation rate and heat fluxes at the surface and vertical profiles of vertical transports of major prognostic variables.
The Q3-D MMF and GCMs with the unified parameterization still represent different families of models.As shown in Fig. 15, however, the both can converge to the same model, a GCRM, as the GCM resolution approaches the CRM resolution.Comparisons of simulated results from these models with those from a GCRM will greatly enhance our understanding of the multiscale role of cumulus convection in the global atmosphere.

Fig. 2 .
Fig. 2. Two families of atmospheric models with different model physics.

Fig. 3 .
Fig. 3. Schematic illustration of typical vertical profiles of moist static energy source under disturbed tropical conditions.Redrawn from Arakawa (2004).

Fig. 4 .
Fig. 4. Ensemble averages of the "required" source of moist static energy obtained through budget analyses of CRM-simulated data for various space/time intervals.Redrawn from Jung and Arakawa (2004).

Fig. 5 .
Fig. 5. Schematic illustration of circulations associated with clouds for (a) coarse and (b) fine resolutions.

Fig. 6 .
Fig. 6.Snapshots of vertical velocity w at 3 km height at the end of two simulations with and without background shear.

Fig. 7 .
Fig. 7.The original domain and examples of the sub-domains used for analysis.

Fig. 8 .
Fig. 8. Dependence of the mean fractional cloud cover at 3 km height on the sub-domain size and its standard deviation.

Fig. 10 .
Fig. 10.Weighted ensemble average of the eddy transport of water vapor at 3 km height.Open circles: estimated with the right hand side of Eq. (13).Closed circles: directly diagnosed with the left hand side of Eq. (13).

Fig. 11 .
Fig. 11.Examples of grid-point arrays used in the Prototype and Q3-D MMFs.

Fig. 12 .
Fig. 12. Schematic illustration of the convergence of the Q3-D MMF grid to a 3-D CRM grid.

Fig. 13 .
Fig. 13.Time series of domain-averaged precipitation, evaporation and sensible heat flux at the surface simulated by the Q3-D CRM (red lines) and by the BM (black lines).

Fig. 14 .
Fig. 14.Vertical profiles of the time-and domain-averaged vertical transports of potential temperature, water vapor mixing ratio, and the horizontal components of vorticity by the Q3-D CRM (red lines) and by the BM (black lines).

Table 1 .
Correlation coefficients between wq − w q and σ (1 − σ ) for various combinations of sub-domain size and averaging size.