## Abstract

Advanced operational four-dimensional variational data assimilation (4DVAR) schemes include a linearized version of moist convective parameterization and its adjoint. At the Meteorological Service of Canada, work is underway to implement 4DVAR for both global and regional operational data assimilation. Moreover, the Kain– Fritsch moist convective parameterization scheme is currently under operational testing for global and regional weather forecasting. Consequently, tangent linear and adjoint versions of this convective scheme have been developed. Sources of nonlinearities and accuracy of the tangent linear approximation of the convective scheme itself were examined. The procedure to test this latter aspect uses Monte Carlo simulations based on background error covariances from the operational three-dimensional variational data assimilation (3DVAR) system at the Canadian Meteorological Centre. It is shown that for a critical level of amplitudes of vertical perturbations of temperature or moisture greater than typically 0.1 K or 0.1 g kg^{−1}, the tangent linear approximation becomes inaccurate (e.g., typical perturbation response having the wrong sign and amplitude errors larger than 100%). For such perturbation amplitudes, there is a rapid increase of convective points where the tangent linear convective approximation is very strongly in error. Deactivation of the Kain–Fritsch scheme becomes frequent and a significant source of invalid tangent linear approximation for input perturbations exceeding typically 0.3 K or 0.3 g kg^{−1}. Potential implications of this study for linearized moist convection in the context of 4DVAR and moist singular vector computation are discussed.

## 1. Introduction

In order to extend the usefulness of linearized physical processes within variational data assimilation and moist singular vector computations, a parameterization of deep convection needs to be considered for global and regional models when convection cannot be explicitly represented. Because of the nature of the physical properties of moist convection, threshold conditions naturally appear. In addition to the intrinsic nonlinearities of such schemes, the activation of threshold conditions represents a major source of mathematical difficulty in formulating a linearized convection scheme because the differentiability condition is violated. In advanced four-dimensional variational data assimilation (4DVAR) used for operational applications [e.g., the European Centre for Medium-Range Weather Forecasts (ECMWF) 4DVAR], a linearized moist convective parameterization scheme is used. The reader can find useful details on the basic assumptions used in the ECMWF linearized physics in Mahfouf (1999). Other operational centers like the National Centers for Environmental Prediction (NCEP) also rely on linearization of their convective parameterization scheme for surface precipitation assimilation within their global three-dimensional variational data assimilation (3DVAR) scheme, for instance (see Treadon et al. 2002). At the Meteorological Service of Canada (MSC), the Kain–Fritsch (KF) moist convective parameterization scheme (Kain and Fritsch 1990) is currently under operational testing (in addition to many other upgrades) for global and regional weather forecasting.

It is thus important for MSC data assimilation purposes to evaluate the usefulness of a linearized version of the KF scheme in order to improve 4DVAR in deep convective regions and also assimilation of precipitation data now becoming available at the Canadian Meteorological Centre (CMC) [e.g., radar and Special Sensor Microwave Imager (SSM/I) data]. There is currently abundant literature available on the subject of precipitation assimilation. In order to give a broad overview of previous studies, we refer the reader to the proceedings of the ECMWF/EuroTRMM Workshop on Assimilation of Clouds and Precipitation (ECMWF/EuroTRMM 2000), and the proceedings of the ECMWF/ GEWEX Workshop on Humidity Analysis (ECMWF/ GEWEX 2002), respectively. The impact of deep convection on moist singular vectors and related subjects is currently a very active area of research with potentially useful impacts on operational products. Details on these aspects can be found in recent scientific publications by Barkmeijer et al. (2001), Errico et al. (2003), Errico and Ehrendorfer (2003), and Coutinho et al. (2004). Work is currently under way at MSC to include more physical processes in the key analysis products potentially useful for operational purposes (Laroche et al. 2002).

In this study, a detailed examination is done of the linearization properties of the KF scheme in operational testing at the Canadian Meteorological Centre. This is done by extracting all the input fields required by the scheme from a 1-day global forecast at 35-km horizontal resolution. The linearization aspects of the scheme are then examined in isolation (vertical column mode). A wide variety of convective cases normally present in a regular operational forecast is considered.

Section 2 presents the basic features of the KF convection scheme. Section 3 goes into the details of the sensitivity (Jacobian) properties of the KF scheme by a partitioning of various contributing terms. The basic approximation used to derive a tangent linear of the KF scheme is then discussed. Once the linearization is performed, we examine in section 4 the degree of validity of the tangent linear approximation using extensive Monte Carlo simulations. Section 5 gives a summary and conclusions together with important implications for the potential usefulness of the linearized scheme in more complete 4DVAR contexts and moist singular vector computations.

## 2. Features of the Kain–Fritsch scheme and design of experiments

Rather than introducing all the features of the KF scheme (which are numerous), we choose the more direct approach of presenting the concepts as they appear in the actual computations of the KF scheme as used at MSC. Because of the complexity of the KF scheme, we refer the unfamiliar reader to standard literature treating the KF scheme, for example, Kain and Fritsch (1990, hereafter referred to as KF90) and Emanuel and Raymond (1993).

### a. Input fields

For linearization as defined in section 3, we need to identify clearly the independent (input) variables of the KF scheme. For numerical weather prediction (NWP), various strategies can be used for computing dynamical and physical terms at a given time step (Staniforth et al. 2002). Without going into such details, it is nevertheless possible to identify the input variables of the KF scheme as

nk: Number of vertical levels

Δ

*t*: Forecast model time step (900 s here)ps: Surface pressure

*T*: Temperature*q*: Specific humidity*ω*: Vertical velocity in pressure coordinate*σ*: Vertical coordinate on which computations are done*A*: Horizontal area occupied by a grid element*u,**υ*: Horizontal wind components.

Of those input quantities, we only look at input perturbations on *T* and *q* fields. Other fields above (except *ω*) only have an impact on the triggering of convection but not on its intensity, except for a slight effect of the horizontal wind components on precipitation efficiency. Consequently, no perturbations on *u, **υ,* and ps are considered.

### b. The triggering conditions

In a typical convective parameterization scheme, there are basically two types of *conditionals* (if statements in a computer code). The first group defines conditions under which the decision is made if convection will be activated. We refer to these conditions as *triggering conditionals.* The latter appear first in the KF sequence of computations. The second group of conditionals appears throughout the rest of the code and characterize various physical mechanisms that can occur depending on the state of the system (some examples are given in the next section). We briefly describe in the following the triggering conditionals since they are directly affected by active variables like temperature and moisture. The first triggering conditional requires that the convective available potential energy (CAPE) value (defined and discussed in detail below) be greater than zero. The second triggering conditional requires that there be enough low-level vertical velocity. Finally, the third triggering conditional is that the activation of the scheme be restricted to deep convective cases giving rise to cloud depths equal to or greater than 4 km in the vertical.

### c. The prognostic equations

Following KF90, the prognostic equations for temperature and moisture tendencies are written as

where *M*_{u} and *M*_{d} are updraft and downdraft mass flux (kg s^{−1}), *M*_{ud} and *M*_{dd} are detrainment rates of the updraft and downdraft (kg s^{−1}), and *M*_{e} is the effective mass of environmental air represented by each model level (kg)

where *w* is the input large-scale vertical velocity to the KF scheme, *M*_{ue} and *M*_{de} are entrainment rates of the updrat and downdraft (kg s^{−1}), *T*_{u} and *T*_{d} are temperature of the updraft and downdraft, *q*_{vu} and *q*_{vd} are specific humidity of the updraft and downdraft, *q*_{cu} is liquid/ solid cloud water content in the updraft, *A* is total grid area, and Γ is dry adiabatic lapse rate.The four terms on the rhs of (2.1) represent, respectively, environmental compensating subsidence (noted TES), detrainment from the updraft and downdraft (TUD and TDD respectively), and evaporation/sublimation of the detrained liquid/solid water (*L*/*S*–*W*) from the updraft (TED).

The way the KF scheme operates is first by computing updraft and downdraft mass fluxes, entrainment/detrainment rates for a given (fixed) radius of the updraft Δ*a* (1.5 km here). Given input thermodynamic profiles, an iterative procedure is used to rescale all these quantities so that after a time integration *τ*_{c} of (2.1)–(2.3), the final CAPE has been reduced by at least 90% of the original value. This iterative procedure generally converges, but there exists some points where it fails to converge. Such points are automatically excluded since this 90% CAPE reduction is strongly enforced within the code. We come back on these aspects in section 3.

### d. Design of experiments

The atmospheric profiles used for describing results of experiments in the next section correspond to a 12-h forecast from a 35-km (horizontal), 46-vertical-level version of the Global Environmental Multiscale model (GEM; Côté et al. 1998) starting from CMC operational analysis at 0000 UTC 1 November 2001. The input fields described in section 2a were extracted from this forecast run, and the KF scheme is used in vertical column mode to examine various convective regimes over the globe (2000 points selected for examination).

## 3. Jacobians: Tangent linear and adjoint

We first show the general properties of the KF sensitivities to small perturbations of temperature and moisture. We present next the strategy to build the tangent linear and adjoint code of the KF scheme based on these properties.

### a. Jacobians

Based on (2.1)–(2.2), we define the input vector **x** and output vector **y**:

where

Each of the **T**, **q** fields are vertical column fields of dimension nk, that is, the number of vertical levels. The convective time tendencies above are obtained according to (2.1)–(2.2).

We are interested in the partial derivatives of the output vector **y** w.r.t. the input vector **x**. Given the KF nonlinear operator noted **f**, its first-order variation **dy** (or linear variation) is given by

We also make use of the *nonlinear variation **d***y** defined as

The Jacobian matrix ∂**f**/∂**x** is explicitly obtained by finite-difference methods, that is,

where temperature and moisture perturbations are considered small (10^{−4} K for *T* and 10^{−4 }*q*_{s} for *q*, where *q*_{s} is the saturation specific humidity). This matrix is further split into its various components based on each contributing terms appearing on the rhs of equations (2.1)–(2.3). For an example, we first show such a partitioning of the full Jacobian in Fig. 1. The heating response for the unperturbed case also appears in the upper-left panel of Fig. 1, and its four contributing terms from Eqs. (2.1) and (2.2) are also plotted. As mentioned in KF90 (section 4c), the dominant contribution comes from subsidence warming (TES). The updraft and downdraft contributions are also very clear (dash and dash– dot curves, respectively), whereas the effect from detrainment of liquid/solid water content appears as the zero dashed line in Fig. 1 since its amplitude is five orders of magnitude smaller than the contribution from subsidence in the environment. In terms of sensitivities (i.e., Jacobians), the other panels bring the additional information that, except for regions close to the top and bottom of the convective cloud where sensitivities are important, the Jacobian of subsidence in the environment dominates by five orders of magnitude the Jacobian of detrainment of liquid/solid water content and explains most of the structure of the total Jacobian. We note however that the TES term in (2.1) depends implicitly on updraft and downdraft properties through the environmental vertical velocity *w*_{e} according to (2.3).

We show with more details in Fig. 2 the heating-rate sensitivities (s^{−1}) w.r.t. temperature perturbations. Again, perturbations of temperature are made at the indicated vertical level on top of each frame, and we can more easily appreciate in Fig. 2 the vertical variation of the response of the KF scheme. The linear variation **dy** (introduced above) is generally restricted to the perturbation level between cloud base (level 44) and cloud top (level 16) except between the layer contained between levels 27 and 39 (500–850 hPa). Over this layer, the linear variation still exhibits the basic local positive/ negative response seen at other levels (such as close to the cloud top), but there appears an extension of the linear variation well above and below this perturbation level. In general these lobes maintain the sign of the linear variation. We discuss more in detail the origin of these ± lobes in section 4a. Similar general conclusions apply for the heating-rate Jacobian w.r.t. *q* perturbations (figures not shown), except that this vertical spread mentioned above appears higher in the cloud. Moistening-rate Jacobians w.r.t. temperature also exhibit vertical spread of the moistening due to local temperature perturbations (figure not shown). The only Jacobian part that is nearly diagonal (i.e., local variation of the response due to local perturbation) is the moistening-rate variation due to *q* perturbations (figure not shown).

The results just presented for a single profile are robust in the sense that for all 1000 convective points examined, these conclusions apply systematically. Based on these Jacobian sensitivities, we now describe the approach taken to develop the tangent linear and adjoint KF schemes.

### b. Tangent linear and adjoint

We have chosen to monitor the behavior of perturbations in the code and examine whether some basic approximation could be used to alleviate the large amount of computation involved in the KF scheme. In addition, this gives us the opportunity to learn more precisely about the origins of strong nonlinearities at different stages in the scheme. We found it possible to identify an approximation to the full Jacobian matrix that we discuss in the following. In order to clearly illustrate where this approximation lies in the chain of computations, we now summarize formally the basic blocks of computation involved in the KF scheme. Block 1: the triggering section is the first part computed. This part remains exactly the same in the nonlinear KF scheme as in the classical tangent linear approach [see Zou (1997) for definitions and extension avenues]. No linearization is involved here since this section can in principle be skipped by a proper flagging of convective points and saving basic-state cloud properties such as cloud base, cloud top, lifting condensation level (LCL), etc. Block 2: assuming a 1.5-km radius at the base of the sample updraft plume, compute updraft, downdraft, and environmental quantities appearing in (2.1)–(2.3), that is,

We now refer to this block as the *kernel computations.* This part of the scheme involves nonlinearities in terms of *T* and *q* input fields (KF90). We come back on the need to perform a linearization of these quantities. Block 3: then follows what is called the “stabilization loop.” This means, starting from the kernel quantities above, the prognostic equations (2.1)–(2.3) are integrated in time over *τ*_{c} (typically 1 h). Then, from the resulting environmental *T* and *q* at the end of the time integration, an updated value of CAPE is computed. If this CAPE value is less than 10% of the kernel CAPE, the computations stop since those *T* and *q* profiles at time zero lead to a sufficient number of homogeneous kernel plumes (updraft/downdraft) properties that produce enough reduction of CAPE in the given time *τ*_{c}. The time tendencies of *T* and *q* are then computed as a simple forward time-differencing scheme over *τ*_{c}, using the final values of *T* and *q* at *τ*_{c} and those at the initial time. If the computed value of CAPE is greater than 10% of the kernel CAPE, a scaling factor *α* is applied to the *kernel flux quantities **M*_{u}, *M*_{d}, *M*_{ud}, *M*_{dd}, and the time integration of (2.1)–(2.3) is performed again but now using the scaled quantities. Usually, four iterations are sufficient for convergence.

From (2.1)–(2.3) and the scaling process just described, the source of nonlinearity (at a given iteration) in these computations comes from 1) the nonlinear dependence of the scaling parameter *α* in terms of the CAPE value and defined as

where the index *j* is the iteration count, and CAPE_{0} is the value at iteration zero, that is, the kernel CAPE, and 2) the nonlinear dependence of the CAPE upon *T* and *q* fields:

Here, LFC and ETL stand, respectively, for level of free convection and equilibrium temperature level (i.e., level where the parcel becomes negatively buoyant). The quantity (*θ*_{e})_{Mix} is the equivalent potential temperature of a well-mixed parcel of air that is lifted vertically along an undiluted pseudoadiabatic profile. The exact formula used for computing this latter quantity is given by [see also Bolton 1980, Eq. (43)]:

where

Quantities appearing with a “Mix” subscript are quantities of the well-mixed parcel over a depth of at least 60 hPa. Tilde quantities (e.g., *T̃, **q̃*) mean physical constraints (e.g., supersaturation, release of latent heat, positiveness) have been applied to the given quantity. These constraints represent typical second group of conditionals referred to previously in section 2b. The constant factors have the following values:

Once (*θ*_{e})_{Mix} has been determined, its value is kept constant in the CAPE integral in the vertical. A similar formula is used to compute *θ*_{e}(*z*), but the overbar means an arithmetic mean is taken in the computational layers used in the discretization of the vertical integral.

For illustration purposes, Fig. 3 shows the vertical distribution of the integrand in the CAPE formula above, that is, CAPE_{k} (left) and the Jacobian of CAPE w.r.t. temperature perturbations at different levels *k* in the vertical, for the same deep convective point as before. Note that only the positive area of the curve on the left panel is taken into account in the CAPE integral. It is clear that the major sensitivity lies close to the LCL (level 35) and rapidly decreases outside this region. If no mixing is applied in the CAPE computation, this large low-level sensitivity becomes enormously larger (Fillion and Mahfouf 2000), showing that, since mixing the parcel before lifting is more realistic, this mixing is largely beneficial to reduce the strong sensitivity of this CAPE value. Nevertheless, as will be seen, the remaining sensitivity of such computed CAPE will remain problematic in terms of tangent linear applications.

Using the procedure described previously to compute the Jacobian, we now examine the rainfall Jacobian obtained by computing the rainfall perturbation to small (see definition of “small” above) input *T* and *q* perturbations to the KF scheme. Restricting our attention to rainfall was very useful in detecting the possibility of greatly simplifying the amount of KF code that needed linearization. This was the starting point of our approximation for building the tangent linear KF scheme. An estimate of rainfall is explicitly done by computing vertical integral of moisture change by the KF scheme. Figure 4 (left) shows the exact rainfall Jacobian w.r.t. *T* (dashed line) and *q* (solid line) perturbations for several typical tropical convective points (similar results were obtained for extratropical points including the profile used in Fig. 1). Right-hand-side panels in Fig. 4 show the same rainfall Jacobian measure but where the kernel quantities (see above) were not perturbed as compared to the basic-state values; that is, only perturbations in the stabilization loop have been retained. It is seen that the structure of the approximate Jacobian is very similar in shape to the full Jacobian. It is important to mention that the kernel CAPE, that is, CAPE_{0}, must be allowed to vary in the computation of the perturbed run with the KF scheme, otherwise these Jacobians would differ in amplitude (not in shape) from the exact Jacobian especially close to the cloud base. This result is very general in terms of profile type chosen since this approximation was examined for 1000 profiles out of the 2000 available for different geographical regions inside and outside the Tropics. Moreover, examination of the details of the heating-rate and moistening-rate Jacobians in the form of Fig. 2 but using our approximation shows very slight differences (figure not shown) for the majority of profiles examined. We come back to the impact of these differences in section 4 when finite-amplitude perturbations are used to compare our tangent linear KF scheme (TLKF) with a perturbation approach TLKF (called “full” TLKF) not using any approximations. The conclusion here is that it is highly accurate to only consider linearization of the kernel CAPE and the stabilization section. This allows us to significantly reduce the construction of the tangent linear KF scheme; that is, only about 25% rather than 75% of the code needs to be linearized, the first 25% of the code being the triggering section that does not require any linearization (classical tangent linear). This extensive examination of the approximation to the full Jacobian was found to be very accurate in terms of rainfall Jacobian and for heating- and moistening-rate profiles.

For clarity and to help understand the results in the next section, we summarize the computations involved in the KF code mentioned above (ignoring the triggering section and fixing the number of iterations of the stabilization loop to a given value niter (niter = 10 in our experiments)]:

Perform stabilization loopdo

*j*= 1, niterSolve equations (2.1)–(2.2)

Compute CAPE

_{j}using (3.2)*α*_{j+1}=*α*_{j}CAPE_{0}/(CAPE_{0}− CAPE_{j})Scale

*M*_{u},*M*_{d},*M*_{ud},*M*_{dd},*M*_{e}with*α*_{j+1}

end do.

## 4. Tangent linear approximation

It was shown in the previous section that an accurate and efficient approximation to the full rainfall Jacobian is to neglect perturbations on the kernel quantities except for CAPE. In order to generate the TLKF, differentiation of the kernel CAPE and the stabilization loop was performed. This represents our approximate TLKF code. The adjoint of the TLKF was also coded. As already mentioned, the major source of nonlinearity that needs to be taken into account for the TLKF comes from the CAPE computation for the kernel and within the stabilization loop.

Although the TLKF scheme (as for the adjoint of TLKF) in its construction is based on infinitesimal perturbations, its use in the context of finite-size perturbations is important to examine for variational data assimilation. What can be done in our context here is to examine the degree of linearity and reliability of the TLKF as compared to the full nonlinear variation of the KF scheme under finite-size perturbations. This amounts to verifying to what extent the *first-order variation ***dy** introduced in section 3a explains the *nonlinear variation **d***y**. We make use in the following of the ratio *R*_{k} = **dy**|_{k}/*d***y**|_{k}, where *k* refers to a given vertical level. Since moist convective parameterization is an integral part of forecast model for current-day operational NWP (global and regional scales), it is useful to characterize the degree of nonlinearity of such parameterizations. It is important to stress that it is the full model operator (i.e., dynamics, diffusion, radiation, etc.) that we ultimately need to examine in terms of the tangent linear assumption, for instance, in the 4DVAR context or moist singular vector computations. This aspect is out of the scope of this paper, but it is nevertheless useful to characterize the degree to which a suboperator (moist convection here) of the forecast model behaves for noninfinitesimal perturbations so as to understand and characterize the potential problems of the full tangent linear model operator.

In order to do a systematic evaluation on many deep convective profiles, we adopt the following Monte Carlo procedure (hereafter noted MC). This procedure is similar to the one used by Fillion and Errico (1997, hereafter referred to as FE97) in the context of the Kuo–Anthes and the relaxed Arakawa–Schubert moist convective parameterization schemes. Using a given number of typical midlatitude profiles (95 precisely), a random generation of *T* and *q* perturbation profiles was done using two a priori probability distributions. The same type of procedure was applied to typical tropical deep convective points, and similar results were found as those reported here. The first perturbation distribution examined is a uniform (zero mean) distribution between ± the standard deviation error on *T* and *q* (typically used in operations at CMC for 3DVAR). This perturbation is applied level by level in the vertical independently and for all profiles. Another type of perturbation (hereafter referred to as B perturbations) was used and stimulated by the fact that in normal applications of tangent linear codes, perturbations have amplitudes and vertical structures that are not independent from one level to the other but are constrained by the presence of background error covariances in 3DVAR/4DVAR contexts. To simulate input perturbations with such vertical structures, eigenvectors/eigenvalues of a Gaussian vertical background error autocovariance matrix 𝗕 for *T* and *q,* respectively, were used. Vertical correlations for *T* and *q* were specified using a Gaussian model in pressure coordinate with a fixed correlation length of 200 hPa for *T* and 100 hPa for *q.* The dominant eigenvectors/eigenvalues appear in the upper panels of Fig. 5; the thick solid line in the upper-left panel represents the vertical correlation of the Gaussian used. Standard deviation profiles used to produce the respective covariance matrices for *T* and *q* appear in the lower panels of Fig. 5. These correspond to typical values used in operations at CMC. Each perturbation vector *δ***T** or *δ***q** (noted here as *δ***x**) satisfies

where the *g*_{i} are random numbers normally distributed with zero mean and unit variance; that is, *N*(0; 1). As in FE97, a scaling factor (*S*) is applied to these perturbations, and the nonlinear regime is examined. In results described below, for each perturbation experiments, we used 10 000 MC simulations. A comparison of both perturbation types mentioned above was done for all experiments described below, and it was found that the TLKF performed slightly better under 𝗕-type rather than random-type perturbations. In consequence, we present in the following only results related to 𝗕-type perturbations.

In our estimation of the accuracy of the TLKF, we have used a criterion on the number of vertical levels in the convective cloud in which a measure of success is defined (see below). This criterion was introduced because when the TLKF begins to be inaccurate for increasingly larger input perturbation amplitudes, the failure generally starts to appear at only very few levels in the vertical. Assuming this local error of the TLKF (more specifically the error of the tangent linear approximation not the error of the TLKF scheme itself) could be monitored in practical situations (not trivial to accomplish however), we regard the remaining vertical column in the cloud as a region where the TLKF can still bring useful information. As an example, we show in Fig. 6 the mean ratio *R* (shown for temperature changes here) obtained when 𝗕-type MC perturbations are used with scaling *S* = 10^{−2} on a given midlatitude profile. Except at level 25 where the value of *R* is around −2.5, the rest of the cloud depth exhibits a mean value of *R* close to one; that is, the TLKF is generally good except at level 24. In the following results to be presented, we have used a vertical extent criteria of at least 70% of the total cloud depth (not necessarily contiguous layers).

Figure 7 shows the result of testing what we call here the “full” TLKF. This tangent linear operator is obtained by using the Jacobian matrix to compute **dy** (see section 3a) where no approximation has been used. The use of a “full” TLKF will allow us to examine the degree of validity of the tangent linear approximation and also compare the degree of accuracy of our approximate TLKF. The upper-left panel shows the intensity of convective surface rainfall computed from the KF scheme for each of these profiles. The scaling amplitude of the MC perturbations appears as the *S* value on top of each frame. The first case examined is the validity of TLKF heating-rate perturbation output (moistening-rate perturbation output is discussed below) under a scaling value of 10^{−3}. This typically means (according to the above-mentioned error specifications from error statistics) that temperature and specific-humidity perturbations are of the order of 0.001 K and 0.001 g kg^{−1}, respectively. The tolerance level of error in the verification is set to 10% throughout at least 70% of the whole cloud depth (as explained above). It is seen that the full TLKF scheme explains the nonlinear variation at this level of accuracy for almost all convective points examined. When increasing by one order of magnitude the amplitude of the MC perturbations, it is seen in the upper-right panel of Fig. 7 that even under a 50% error level tolerance, the level of accuracy of the full TLKF is significantly reduced (except for some convective points). Clearly, for such amplitudes, that is, typically 0.01 K and 0.01 g kg^{−1}, the full TLKF is less accurate. When the range of perturbation size reaches 0.1 (Fig. 7, lower-left panel), there are still many points where the TLKF shows some useful information up to around 50% of MC simulations, but equivalently it fails to provide useful information for the other half of the convective points for about half of the total simulations. This situation continues to deteriorate as we increase the size of the perturbations to 0.3 (Fig. 7, lower-left panel). For such amplitude of perturbations, the full TLKF is generally invalid 70% of the MC simulations for almost all convective points. At this level of perturbation, it is important to realize that we have reached a critical level where the convective scheme itself may be deactivated because of the perturbation size (i.e., triggering conditional active). Those situations were referred to as “pathological cases” in FE97 because of their discontinuous effect during the minimization (i.e., disappearance of the convective operator). We show in the lower-right panel of Fig. 7 the percentage of occurrence of such pathological cases as we perform the MC perturbations at *S* = 0.3. For such a case, many convective points become convectively inactive around 50% of the total MC simulations. We conclude here that for perturbation amplitude exceeding about 0.1 K and 0.1 g kg^{−1}, the TLKF has a substantial level of error that makes it almost completely inaccurate because of combined nonlinearities and threshold conditions, including possibly the complete deactivation of convection.

The previous results were obtained with the “full” TLKF. We now compare the degree of accuracy of our approximate TLKF w.r.t. the “full” TLKF. It is seen in the upper-left and lower-left panels of Fig. 8 that except for some convective points, the approximate TLKF have about the same level of accuracy as the “full” TLKF. For perturbation regimes where *S* = 0.1, 0.3, and 0.5, there is a rapid decrease of TLKF validity both for the “full” and approximate TLKF as is obvious in Fig. 8. For cases where *S* > 0.1, it can be said that our approximate TLKF is slightly less accurate than the “full” TLKF, but the latter is already a very poor approximation of nonlinear variations observed in the KF scheme. We conclude here that it may be slightly useful to take into account the variability of the kernel quantities (described in section 3) when deriving the complete TLKF, but this will only help for perturbations of size typically less than 0.1 K and 0.1 g kg^{−1}.

It is interesting to quantify the level of TLKF error in terms of the cumulative temperature and moisture changes. As mentioned in the preceding sections, the KF scheme at MSC is used with fixed time tendencies for a given time *τ* ∼ 1 h; that is, the heating and moistening rates computed from the KF scheme (once convection is first detected from the triggering conditionals) are kept constant during time *τ.* The error between the TLKF response—that is, **dy** and the nonlinear differential response *d***y** in terms of temperature and moisture changes after time *τ*—was computed, and we present in Fig. 9 the value of the maximum standard deviation error found in the vertical depth of the cloud for each convective point after MC perturbations were performed. The range of the scaling parameter used is 0.1, 0.3, and 0.5 (left to right), and results are shown in terms of errors in temperature changes after 1 h. For all scaling factor *S* considered, results show that the level of error is typically of the order of 0.2 K and generally less than 0.4 K. What is seen in the *S* = 0.3 panel in Fig. 9 is that in general, even if the TLKF scheme (full here) as discussed above is almost completely all the time in error, the amplitude of this error still remains generally less than 0.4. We conclude here that even though the TLKF is known to be inaccurate for a range of input *T* and *q* perturbation sizes, temperature perturbation errors are, however, only of the order of 0.2 K for most convective points considered.

### a. Constant mass-flux approximation

For the development of a tangent linear code for the Tiedtke (1989) moist convective parameterization scheme, Mahfouf (1999) considered an approximation where the mass flux of the *single plume* cloud model was kept constant. We examine here the degree of accuracy of this approximation in our *multiple plumes* context with the KF scheme. Our approximation for TLKF is not to keep the mass flux constant but only the *kernel* mass flux. This means we allow variations of total mass fluxes by allowing variations of the scaling parameter *α* in the stabilization loop. Remember that the role of the stabilization loop is to find the number of homogeneous plumes (of the kernel type with same properties) to reduce sufficiently the convective instability after a given time. We allow this number to vary in our approximation. The value of the scaling parameter *α* (see section 3) determines this number (a real number in fact). By allowing *α* to vary, we allow a fluctuation of the total mass fluxes. It is worth noting that because of the dependency of *M*_{e} upon the updraft and downdraft flux properties [see (2.4)] the rhs terms in (2.1), for instance, are not linear in terms of *α*; that is,

The “constant mass flux” approximation can be tested here by forcing *α* to the basic-state (i.e., unperturbed) value. We note that such an approximation reduces the Jacobian structures (such as those shown in Figs. 1 and 2) to an almost diagonal structure (i.e., nonzero elements of heating-rate perturbation only appear at the level of the local perturbation used to produce the Jacobian elements in the perturbation method). This can easily be seen from Eq. (4.1) when *M* terms are assumed constant and equal to the basic state. The only contribution to nondiagonal terms is the environmental subsidence-heating term (TES) (and similarly for the moisture tendency equation) where a vertical gradient of the input vertical profile of temperature is involved. Because of this latter effect, the vertical slope above and below the perturbation level are of opposite sign, which explains the two lobes of opposite signs in the Jacobian structures in Fig. 2 and also (but less apparent) in Fig. 1 above and below the diagonal. For illustration purposes, we show in Fig. 10 a superposition of the heating-rate Jacobians for the same convective point as Fig. 2 where the constant-mass-flux assumption is used (solid lines) as compared to our approximation of constant kernel mass flux (dotted line), with the full Jacobian being essentially superimposed on our approximate Jacobian dotted lines. It is apparent from Fig. 10 that neglecting off-diagonal elements in these Jacobians can be a source of large errors in the linear regime since at a given vertical level *k* we have a sum of contributions

From Fig. 10, the departures in partial derivatives between the two approaches at levels *i* = 28–39 contribute to produce a significant error in the constant-mass-flux TLKF response as compared to the full TLKF for all levels *k* in the cloud. If these sensitivities are maintained from infinitesimal to finite-amplitude perturbations, large errors can also arise. Interesting related problems of tangent linear failure due to change in the number of active plumes were obtained by Bélanger (1999) with the multiple-plumes relaxed Arakawa–Schubert (RAS) convective parameterization scheme. The problems are less severe here than with the RAS scheme since fluctuations of active plumes in the KF scheme are considered only with fixed vertical extension and kernel properties (the homogeneous cloud assumption), as opposed to RAS, which uses many different vertical plume extents and properties.

In order to test the constant-mass-flux approximation, the value of *α* at the last iteration of the stabilization loop for the basic-state KF computation is stored and used again in the perturbation step (without iterations of the stabilization loop). A similar set of MC experiments with a scaling factor ranging from *S* = 10^{−3} to 10^{−1} showed that the constant-mass-flux approximation fails to pass the 50% error level criteria for all range of amplitudes of perturbations and for all 95 convective points used previously. We conclude here that the constant-mass-flux approximation is not an acceptable solution in the context of the multiple plumes KF convective scheme when temperature perturbation is an important prognostic variable in the tangent linear application considered (e.g., certain form of incremental 4DVAR where tangent linear observation operators are used explicitly and moist singular vector computations). As demonstrated previously, a better approximation is to account for variations in the number of plumes while the mass flux of the kernel plume is held constant, albeit even if the full TLKF has large percentage errors at and above the 0.1-K temperature perturbation intensity.

### b. 1DVAR precipitation assimilation

We developed and tested a one-dimensional variational data assimilation (1DVAR) scheme of simulated precipitation observations using the TLKF and its adjoint in a way similar to FE97. The results (not shown) are overall very similar in terms of gradient-norm reduction and precipitation adjustment (e.g., two to three iterations needed for convergence). Similar results were obtained with the full or approximate TLKF. The only result worth mentioning is that the analysis increments obtained by 1DVAR are mostly concentrated on *moisture changes* rather than both temperature and moisture changes as in the 1DVAR results of FE97 with the relaxed Arakawa–Schubert moist convective parameterization scheme (see also Fillion and Mahfouf 2000). Since very similar background error standard deviations for temperature and specific humidity were used for both 1DVAR experiments mentioned above, we believe this latter difference is due to the relative importance of temperature and moisture contribution to heating-rate Jacobians combined with a different vertical number and placement of vertical levels.

## 5. Summary and conclusions

In this study, we have examined the linearization of the Kain–Fritsch moist convective scheme. The latter was used within a 12-h forecast with the Canadian global NWP model GEM run at 35-km resolution in the horizontal and 46 vertical levels. By extracting all the input fields necessary for the computations in the Kain– Fritsch scheme, it was possible to concentrate specifically on its nonlinearities. In the context of the classical tangent linear approximation, it was possible, by a careful look at the propagation of perturbations step by step in the KF code, to observe that the kernel quantities of the updraft and downdraft (i.e., those of the homogeneous plume model) are relatively stable to small perturbations. This is true except for the kernel CAPE. It was possible to show that the approximate rain-rate Jacobian on various convective profiles where the kernel quantities were held fixed provides a very accurate approximation to the full rain-rate Jacobian. Using this approximation, the tangent linear and adjoint of the Kain–Fritsch scheme were derived. The departures between our approximate tangent linear Kain–Fritsch scheme and a more complete approach not using this approximation have been examined.

Monte Carlo simulations in the form used by Fillion and Errico (1997) were used to examine the degree of validity of the tangent linear approximation. The MC perturbations were derived from CMC vertical error covariances for temperature and moisture. These perturbations are also more appropriate considering the presence of the background error vertical covariances in a full fledged 3DVAR/4DVAR context.

Generally, for midlatitude and tropical convective profiles, the tangent linear of the Kain–Fritsch scheme quickly becomes inaccurate when the level of temperature or moisture perturbations exceeds 0.1 K and 0.1 g kg^{−1}. It was shown that the level of error of the tangent linear scheme for temperature and moisture perturbations (with or without our approximation of constant kernel properties), is of the order of 0.1 K and 0.1 g kg^{−1} if perturbation tendencies are kept constant at each time step during typically 1-h convective time scales (as done at MSC in our application of the KF scheme). The constant mass-flux approximation used at ECMWF by Mahfouf (1999) for the Tiedtke (1989) single plume convective scheme was also tested in our multiple plume context with the Kain–Fritsch scheme. Because of nonlocal variation of the response to local temperature and moisture perturbations, this approximation for the tangent linear of the Kain–Fritsch scheme was found inaccurate at all range of perturbation sizes.

One source of nonlinearity limiting the accuracy of the tangent linear KF scheme is the CAPE computations, which are strongly sensitive to low-level perturbations. Other closure assumptions could presumably be examined to alleviate this strong source of nonlinearity. Flexibility to use various closure assumptions in convective parameterization schemes for advanced data assimilation schemes has already been introduced by Grell and Devenyi (2002). This is an avenue to be explored according to the limitations of the tangent linear convection scheme examined here and presumably occuring also in similar types of convective parameterization schemes in operational use. It can also be noted that (a) the degree of nonlinearity of a given scheme (e.g., the KF scheme) can be examined in itself without the requirement of developing a tangent linear and adjoint version of the scheme. Ways to simplify the original nonlinear scheme by one which is more linear can then be explored. (b) Minimization approaches not relying on differentiability of the cost function need to be examined and compared with the current operational practice, especially when considering the assimilation of nonconventional data such as cloud and precipitation data from satellite instruments.

The implications of this study concern the current limited reliability of the tangent linear scheme when input perturbations of temperature and moisture exceed typically 0.1 K and 0.1 g kg^{−1}. Very similar results were obtained in the past with other convective parameterization schemes such as the Kuo–Anthes scheme (Fillion and Errico 1997), the relaxed Arakawa–Schubert scheme (Fillion and Errico 1997; Bélanger 1999), and a version of Tiedtke scheme as used at ECMWF (Lopez and Moreau 2003). In those implementations of incremental 3DVAR/4DVAR making explicit use of tangent linear observation operators, for instance, temperature and moisture analysis increments during minimization can initially exceed 1 K and 1 g kg^{−1}, thus exceeding the above limitations of current-day tangent linear convective parameterization schemes. As the amplitude of analysis increments decrease during minimization iterations of 3DVAR or 4DVAR for instance, it is justified to delay the introduction of such a tangent linear physical operator (assuming it is used) after a given number of iterations during the minimization. Also, in the context of moist singular vector computations where these are by design the most intensifying structures, as soon as the tangent linear convective scheme is used in the regular chain of physical operators of the forecast model with such perturbation amplitudes at the initial time, the perturbed convective response coming out of the tangent linear scheme will be strongly in error. It would be justified and relevant for future studies on 4DVAR or moist singular vector computations to clarify the above-mentioned weaknesses of tangent linear convective parameterization schemes with a full tangent linear forecast model over geographical regions where moist convection plays a significant role. As it stands now, it is not a priori clear that a tangent linear Kain–Fritsch convective scheme can be a reliable addition to a tangent linear forecast model because of its very limited accuracy. This latter aspect needs to be carefully monitored, otherwise local instabilities could easily be triggered in the tangent linear forecast model. In the perturbation approach to build a Jacobian matrix for the tangent linear of convection, these large tendencies can be monitored and removed using a criterion on the size of eigenvalues of the Jacobian matrix of the convection scheme (Errico and Raeder 1999).

## Acknowledgments

The authors would like to acknowledge stimulating discussions with researchers in data assimilation and NWP at the Global Modeling and Assimilation Office (NASA GMAO), European Centre for Medium-Range Weather Forecasts (ECMWF), Meteorological Service of Canada (MSC), National Center for Atmospheric Research (NCAR), and National Centers for Environmental Prediction (NCEP) on early results of this study.

## REFERENCES

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

**,**

## Footnotes

*Corresponding author address:* Dr. Luc Fillion, Direction de la Recherche en Météorologie, Meteorological Service of Canada, 2121 Route Trans-canadienne, Dorval QC H9P IJ3, Canada. Email: luc.fillion@ec.gc.ca