Slovenia-PECS Implementation Contract 4000106730/12/NL/KML
Multivariate relationships between the aerosols, moisture and winds in
four-dimensional data assimilation for the global monitoring for environment and
security
Final report
Nedjeljka ˇ Zagar and ˇ Ziga Zaplotnik
Department of Physics
Faculty of mathematics and physics University of Ljubljana
Ljubljana 2018
Multivariate relationships between the aerosols, moisture and winds in four-dimensional variational data assimilation for the global monitoring for environment and security by Nedjeljka ˇZagar ˇZiga Zaplotnik and Nedjeljka ˇZagar
iii
Executive summary
The increasing amount of remotely sensed data on atmospheric trace constituents has been recently provided by satellites in addition to numerous observations of the atmospheric vertical temperature and moisture profiles (in form of radiances). The projects GEMS and MACC have developed a comprehensive monitoring and forecasting systems for trace atmo- spheric constituents important for climate and air quality as a part of Europe’s Global Mon- itoring for Environment and Security (GMES, http://www.esa.int/esaLP/LPgmes.html).
Both MACC and GEMS build on the global NWP system operated by the European Centre for Medium-Range Weather Forecasts (ECMWF). In these projects, the ECMWF data assimilation system based on the four-dimensional variational data assimilation (4D- Var) is combined with the expertise of the diverse research groups engaged in atmospheric composition modelling to build an integrated monitoring and forecasting system.
At the same time, a lack of wind measurements over many parts of the globe has been recognized as the main shortcoming of the current global observing system. New wind observations will soon become available thanks to the first spaceborne Doppler wind li- dar Aeolus (https://www.esa.int/Our−Activities/Observing−the−Earth/Aeolus) that will provide the global coverage with wind profiles twice per day. However, this is a demon- stration mission of limited duration and with no follow-on mission the gap in the wind observations will remain. The observation gap is especially large in the tropics where data assimilation for numerical weather prediction (NWP) faces challenges that far exceeds those in the midlatitudes. Resulting large uncertainties in tropical analyses have impeding effect on our knowledge of weather and climate.
This provides a strong motivation to explore data assimilation methodology in order to utilize existing observations of trace constituents. Very little is known about the coupling between the aerosol, moisture and winds in data assimilation. In the ECMWF system developed for GMES, the aerosols are treated as passive scalars without impact on analysis increments for other variables. In contrast, moisture cannot be regarded as a pure tracer as it actively impacts dynamics. Moisture observations influence the quality of wind analyses in 4D-Var assimilation, but details of the process are difficult to understand in a full scale NWP system. An intermediate modelling framework which would contains the most important interactions between the aerosols, moisture and dynamics in 4D-Var would be beneficial. Such a framework has been developed within this project.
Building on the previous work, the project has developed a new model of intermediate complexity to simulate the aerosol, moisture, temperature and wind interactions in 4D- Var. The Moist Aerosol Dynamics Data Assimilation system (MADDAM) mimics in many details the operational systems for NWP such as the ECMWF. MADDAM uses the spectral methods and the transformed moisture variable for the assimilation to account for the non-Gaussian properties of the humidity. In MADDAM, observations of temperature field produce analysis increments in the wind field and vice versa. In contrast, moisture and aerosol variables are analyzed univariately so that their feedbacks on the wind field are due to internal model dynamics which adjusts assimilated observations, the background-error covariances and the dynamical and physical processes represented by the model equations.
iv The now model is applied for the estimation of the ability of 4D-Var assimilation to re- trieve the wind-field information from the time series of aerosol and moisture observations.
Numerous experiments based on the the observing system simulation experiment (OSSE) approach reveal potential of aerosol observations in wind tracing.
In the atmosphere near saturation, it is essential to have reliable information about the temperature field along with humidity observations in order to extract useful wind information using 4D-Var. The reason is the susceptibility of the condensation process to slightest changes in saturation humidity which is temperature dependent. It is demon- strated that useful wind information can be retrieved from moisture observations even if the flow is highly nonlinear and precipitation takes place, provided there are significant hu- midity gradients. Adding humidity to temperature observations systematically improves wind analyses with respect to the case when only temperature data was assimilated, es- pecially in a longer 4D-Var assimilation window. Although these results are not obtained with the full-scale NWP system, the highlighted factors affecting the wind retrieval from moisture data are likely explain significant analysis increments from humidity observations in operational systems.
In the case of aerosol observations in saturated atmosphere, significant retrieval of wind information from aerosols requires that all background thermodynamical variables are simulated well. A successful retrieval depends on the level on nonlinearity. The more nonlinear the problem, the slower the convergence of the variational minimization process.
A set of OSSE experiments addressed the impact of observation sampling and their errors in relation to the assumptions applied in data assimilation modelling and dynamical and physical processes in the model. The most important conclusions are as follows. The 4D-Var of dense aerosol observations (available at every or every second grid point) assim- ilated with hourly time interval can in ideal circumstances of the perfect-model 4D-Var improve wind analysis for 10% to 15% with respect to the background. The 24-hour assim- ilation window contributes a few percent of the improvement on the top of that obtained using the 12-hour window. For comparison, the assimilation of humidity observations in the same environment is about twice more efficient and a combined assimilation of dense humidity and temperature data improves wind analysis for 5-10% compared to the case with only temperature observations. These values should be compared with the 60-70%
improvement of the wind analysis by the assimilation using wind observations. The aerosol observations barely impacts wind analyses when added to other observation types. Reduc- ing the temporal sampling of aerosol and moisture observations gradually reduces their impact on the wind analysis as the advection process becomes more difficult to reconstruct by 4D-Var in spite of the high spatial sampling.
The modelling choices regarding the background-error covariance representation for the moisture and aerosol forecast errors effect their information content in the assimilation.
The applied short background error correlation lengths were shown appropriate and in the agreement with the suggested need for the spatially dense moisture and aerosol observations to reconstruct winds. Sensitivity to the correlation length reduces as the assimilation time window becomes longer since 4D-Var is more efficient to deduce information about the advection.
MADDAM illuminated the way moist processes decrease the efficiency of 4D-Var to
v extract winds from aerosol and humidity observations. The main obstacle for the wind extraction from the time series of tracers, regardless of their spatial distribution, is the saturation process in 4D-Var. Nevertheless, it is shown that 4D-Var of dense aerosol obser- vations can provide useful information about the wind field, more so in unsaturated regions.
The assimilation of dense humidity observations is however far more useful. This is illus- trated by Figure 1 which compares zonal and meridional wind analysis errors normalized by their respective background errors separately in the unsaturated and saturated parts of the model atmosphere. The efficiency of dense temperature observations to retrieve wind field in 4D-Var in both dry and moist atmosphere is largely increased by adding moisture observations on top of temperature data.
Figure 1: Analysis root-mean-square errors normalized by the background errors for the zonal wind and meridional wind in 4D-Var assimilation using the 24-hour assimilation window. The normalized errors are shown as a function of the spatial observation density for different observation types: aerorol (c), moisture (q), temperature (t) and moisture and temperature together (tq). Results of the wind tracing in unsaturated regions are shown in red whereas the results in saturated areas are presented in blue.
A number of presented features in the results calls for further studies and a more detailed comparison with the operational systems such as ECMWF in order to propose concrete advancements in data assimilation methodology that will improve weather and climate prediction.
Contents
1 Motivation and Goals 3
2 Forecast model and data assimilation methodology 8
2.1 Formulation of moist dynamics model with aerosols . . . 9
2.1.1 Equations for dynamics . . . 9
2.1.2 Moist processes . . . 11
2.1.3 Aerosol equation . . . 12
2.1.4 MAD numerical formulation . . . 13
2.1.5 Model setup . . . 15
2.1.6 Forecast experiments using MAD . . . 15
2.2 Data assimilation modelling . . . 22
2.2.1 Four-dimensional variational data assimilation . . . 22
2.2.2 Data assimilation modelling of dynamics . . . 23
2.2.3 Data assimilation modelling of moisture . . . 25
2.2.4 Aerosol data assimilation modelling . . . 30
2.2.5 MADDAM tangent linear and adjoint model . . . 31
2.3 Setup of the assimilation experiments . . . 33
2.3.1 Observing system simulation experiments . . . 33
2.3.2 Model domain and numerical setup . . . 34
2.4 Single observation experiments . . . 34
2.4.1 A single zonal wind observation in 12-hour 4D-Var . . . 35
2.4.2 A single temperature observation . . . 35
2.4.3 A single humidity observation . . . 38
3 Factors affecting the 4D-Var wind tracing 41 3.1 Wind tracing in homogeneous dry atmosphere . . . 41
3.2 Impact of nonlinearities on wind tracing in dry atmosphere . . . 46
3.3 Wind tracing in saturated atmosphere . . . 49
3.3.1 Wind tracing from humidity observations . . . 49
3.3.2 Wind tracing from aerosol data . . . 51
3.4 Wind tracing in multiscale moist flow . . . 55
3.5 Summary of findings . . . 61
CONTENTS 2 4 Quantification of the potential for wind tracing from 4D-Var assimilation
of aerosols and moisture 62
4.1 Preparation of OSSEs . . . 62
4.2 An example of the added value of tracer observations in the tropics . . . 65
4.3 Impact of observation density and the length of the assimilation window . . 68
4.4 Impact of observation frequency . . . 71
4.5 Impact of observation accuracy . . . 74
4.6 Impact of the background correlation length . . . 75
4.7 Differences between wind tracing in dry and moist atmosphere . . . 77
5 Conclusions and Recommendations 80
References 83
Chapter 1
Motivation and Goals
The increasing amount of remotely sensed data on atmospheric trace constituents has been provided by satellites in recent years in addition to numerous observations of the atmospheric temperature and moisture (in form of radiances). The latter are thinned to a great extent. The thinning is applied to a much smaller extent to wind measurements as their amount is much smaller. In fact, wind observations, especially the vertical wind profiles remain the main missing component of the global observing system (e.g. Baker et al., 2014).
The projects GEMS (Global and regional Earth-system (Atmosphere) Monitoring using Satellite and in-situ data, http://gems.ecmwf.int/) and MACC (Monitoring Atmospheric Composition and Climate, http://www.gmes-atmosphere.eu/) have recently developed a comprehensive monitoring and forecasting systems for trace atmospheric constituents im- portant for climate and air quality as a part of Europe’s Global Monitoring for Environment and Security (GMES, http://www.esa.int/esaLP/LPgmes.html). Both MACC and GEMS build on the global NWP system operated by the European Centre for Medium-Range Weather Forecasts (ECMWF). In these projects, the ECMWF data assimilation system based on the four-dimensional variational data assimilation (4D-Var) is combined with the expertise of the diverse research groups engaged in atmospheric composition modelling to build an integrated monitoring system (Hollingsworth et al., 2008).
At the same time, a lack of wind measurements over many parts of the globe has been a standing shortcoming of the current global observing system (Baker et al., 2014). New wind observations will soon become available thanks to the first spaceborne Doppler wind lidar (Stoffelen et al., 2005) that will provide the global coverage with wind profiles twice per day. However, this is a demonstration mission that will last a couple of years only. A gap in the observing system will remain.
The gap is especially large in the tropics where data assimilation for NWP faces chal- lenges that far exceeds those in the midlatitudes. Figure 1.1 shows the radiosonde obser- vations that were received at ECMWF for the assimilation purposes on 3 January 2017 at 00 UTC. The distribution of radisonde sites is very unequally distributed and constrained to the continental areas, mainly over the Maritime Continent and South America. Wide regions (the eastern Pacific, Indian Ocean) remain virtually void of any direct wind obser- vations. The belt within 30◦ off the equator represents one half of the Earth surface, but
4 the radiosonde measurements over the tropical land areas make only a small percentage of the global radiosonde observations.
0°N
30°S
60°S 30°N 60°N
0°N
30°S
60°S 30°N 60°N 0°E
30°W 60°W 90°W 120°W
150°W 30°E 60°E 90°E 120°E 150°E
0°E 30°W 60°W 90°W 120°W
150°W 30°E 60°E 90°E 120°E 150°E
0°N
30°S
60°S 30°N 60°N
0°N
30°S
60°S 30°N 60°N 0°E
30°W 60°W 90°W 120°W
150°W 30°E 60°E 90°E 120°E 150°E
0°E 30°W 60°W 90°W 120°W
150°W 30°E 60°E 90°E 120°E 150°E
Total number of obs = 688 03/Jan/2017; 00 UTC
ECMWF Data Coverage (All obs DA) - Temp
687 LAND 1 SHIP 0 DROPSONDE 0 MOBILE
Magics 2.29.0 (64 bit)
Figure 1.1: Global distribution of radiosonde observations at a randomly chosen recent date (3 January 2017, 00 UTC). Measurements from 688 radiosonde stations around the globe observed the winds. Source:
ECMWF.
In contrast to the extratropics, a predominant balance between the mass and wind fields in the tropics does not exist. As a consequence, the wealth of indirect temperature measurements provided by satellites is not as useful to constrain the wind field as it is in the extratropics, not even in the perfect-model case (ˇZagar et al., 2005). In parts of the troposphere, atmospheric motion vectors (AMVs) data represent an important source of wind observations, but in the upper troposphere and lower stratosphere the only regular information on winds is provided by radiosoundings. As a result, there are occasional large discrepancies between stratospheric analyses by different operational centres (e.g.
Podglajen et al., 2014; Baker et al., 2014). Figure 1.2 from Baker et al. (2014) illustrates wind analysis uncertainties at 300 hPa in 2010.
Atmospheric composition affects circulation in different ways, either directly by chang- ing radiative fluxes or indirectly by modifying cloud properties. However, these properties become important only on the climate scale, and have little effect on the synoptic scale.
Nonetheless, not much is known about the dynamics of the coupling between the aerosol, moisture and winds in the data assimilation process. In the ECMWF DA system developed for GMES (Benedetti and Coauthors, 2009), the aerosols are treated in a conservative way - as passive scalars, subjected only to advection, convection and diffusion. The feedback of the aerosol analysis increments on the wind is turned off to avoid potential spurious wind increments due to observational biases. That means that the aerosol prognostic equation is
5
Figure 1.2: Root mean square differences in 300 hPa wind speed (m s−1) analyses produced by ECMWF and GFS in a one year period from January to December 2010. From: Baker et al. (2014).
not used as a strong constraint. Similarly, the impact of stratospheric ozone observations on the wind analysis is also turned off (Han and McNally, 2010). On the other hand, moisture observations in ECMWF 4D-Var system influence the wind field, both in the tropics and the mid-latitudes. This was first observed by Andersson et al. (1994) and later by Bormann and Th´epaut (2004), Geer et al. (2008) and Peubey and McNally (2009).
For example, Bormann and Th´epaut (2004) showed that the water vapor observations re- trieved from the MODIS satellite can be used for deriving high-latitude tropospheric wind information. In 4D-Var assimilation, mass observations affect wind analysis through the balance constraints in the background-error covariances and through the internal model adjustment during forward and adjoint model integration (ˇZagar et al., 2004b; Bonavita and Holm, 2016).
This study was motivated by the question whether the time series of spatially dense observations of aerosol concentrations may produce a positive impact on wind analysis, similar to moisture observations. The aerosol distribution patterns often involve sharp horizontal gradients suggesting their potential to describe the transport properties. Wind retrieval from aerosol considered as a perfect tracer has been a subject of a number of studies. However, the combined effects of aerosols, moisture and temperature observations on wind tracing has not been studied. The goal of this project is to highlight the potential of these combined mass-field observations in comparison to direct (missing) wind data in interaction in the past which are reviewed in 4D-Var.
The theoretical foundation for studying wind retrieval from perfect tracer observations was provided by Daley (1995, 1996) used analytical 1D and 2D transport models and an extended Kalman filter. With no sources or sinks of constituents in the model, he concluded
6 that the winds can be retrieved in case of sufficient tracer field variability (i.e. large enough spatial gradients in tracer field) and sufficiently frequent, dense (data voids small) and accurate (particularly for low constituent concentrations) observations. A variety of studies addressed the problem of wind extraction from observations of perfect trace gases in 4D-Var with reference to stratospheric ozone. Varying levels of realism included idealized models of different complexities and exclusively simulated observations (Riishøjgaard, 1996;
Allen et al., 2014, 2016) and NWP systems with simulated observations (Peuch et al., 2000; Allen et al., 2013). Studies based on perfect observations in general reported a positive wind retrieval outcome. However, positive results were less clear by using real observations in NWP environment (Semane et al., 2009) with the main limiting factor being the availability and accuracy of the observations.. Even if in principle the aerosols might be good tracers, in the large system such as ECMWF the aerosol couplings with dynamical variables are so complex that the aerosol impact on the wind field might appear unrealistic. Thus in the operational setting, the impact of the assimilation of ozone data on the wind analysis is turned off to avoid spurious wind increments due to biased observations in certain atmospheric conditions (Han and McNally, 2010)
Studying the aerosol, moisture and wind interactions within an idealized modelling framework has some important advantages. An idealized (with respect to NWP) frame- work allows to develop, implement and test new algorithms and to perform numerical experiments in the controlled environment in which various issues, difficult to grasp in real NWP case, can easier be understood. However, simplified models should still be complex enough to capture main dynamical and physical aspects of the phenomena of interest in order to explain the observed features of circulation and to be of any value of NWP. An important such model is based on the rotating non-linear shallow water equations (e.g.
Vallis, 2006) which include bot balanced (vorticity dominated) dynamics and gravity-wave dynamics as well as their interactions. Shallow water models (SWMs) were applied a number of data assimilation studies to develop new concepts and to study the value of mass-field and wind-field observations (e.g. ˇZagar, 2004, and references therein), as well as to evaluate the impact of future observing systems (e.g. ˇZagar et al., 2008).
The envisaged modelling framework shares the same basic dynamical properties as the ECMWF model and therefore can highlight most relevant aspects of the aerosol-wind feed- back. A special feature of the proposed modelling framework is a novel data assimilation methodology that takes into account the large-scale inertio-gravity (IG) waves in the trop- ics in the representation of the forecast error covariances (ˇZagar et al., 2004a). This is a step forward compared to the traditional data assimilation approach that does not consider the IG waves explicitly in the background-error term. Namely, the standard approach has been focused on the midlatitudes where the largest forecast errors develop in the storm track regions and can be understood with the help of quasi-geostrophic dynamics (Daley, 1991). In the tropics, the equatorially trapped inertio-gravity waves play a role in dynamics at all scales, and they are characterized by phase speeds that are comparable with those of the equatorial Rossby waves. ˇZagar et al. (2005, 2007, 2013) showed that the equatorial large-scale inertio-gravity waves contribute a large portion of the short-range forecast-error variance statistics.
This project reports presents the potential of incremental 4D-Var assimilation to extract
7 wind information from the observations of moisture and aerosol concentrations alone or together with temperature observations. This is achieved by developing a new, intermediate complexity model for studying the aerosol-moisture-wind interactions in 4D-Var, which albeit simplified, shares a majority of assimilation relevant properties with NWP systems.
Performed 4D-Var experiments address the spatial and temporal sampling of observa- tions, and the role of observation errors. The role of model dynamics, especially nonlin- earities, in relation to the applied assimilation window in 4D-Var is studied in detail. It is shown how nonlinearities, associated with shear instabilities, moist processes and aerosol- moisture interactions such as deposition (aerosol sink) undermine the assumption of aerosol being a perfectly advected tracer and complicate the wind tracing.
The report consists of five chapters. The following chapter presents the modelling framework. Its first part describes the development of the forecasting system involving the interaction between the moisture, tracers and dynamics, and demonstration of the complex dynamics represented by the new model. The second part of this chapter describes the variational data assimilation modelling. Factors affecting the retrieval of aerosols in 4DVar are discussion by using a number of simple experiments in chapter 3. Chapter 4 presents comprehensive results from the ensembles of experiments with tracers, moisture and dynamics in in 3D-Var and 4D-Var. In particular, multiple experiments with the assimilation of tracers and moisture observations are compared with the efficiency of data assimilation using observations of dynamical variables. Conclusions and recommendations are presented in chapter 5.
Chapter 2
Forecast model and data assimilation methodology
In this chapter we present the applied modelling framework: the moist prognostic model and the formulation of the 4D-Var with moist processes. The two components of the system are abbreviated MAD, the Moist Atmosphere Dynamics model or Moist Aerosol Dynam- ics model depending on whether the model is used to simulate the aerosols along with the moisture, and MADDAM, Moist Aerosol Dynamics Data Assimilation Model. It is a medium complexity forecasting system which applies the state-of the-art data assimilation methodology. First we present details of the forecasting system involving wind, temper- ature, specific humidity and total aerosol mixing ratio as the prognostic variables. The second part describes the novel data assimilation modelling approach.
The foundations for the new model developed within the project were set in ˇZagar et al. (2004b,a, 2008) and ˇZagar (2012). Extensive new developments performed within the project include the following:
1. tangent-linear and the adjoint of the discretized equations, 2. incremental 4D-Var formulation (inner/outer loop),
3. model extension with the aerosol conservation equation, 4. a simple scheme for large-scale precipitation,
5. a simple physical schemes for aerosol wet and dry deposition,
6. a new moisture control variable and its background-error covariance model, and 7. aerosol assimilation including a new aerosol background-error covariance model.
Numerical aspects of the model and data assimilation have been updated and improved.
The aforementioned earlier studies applied tangent-linear (TL) and adjoint (AD) of the continuous equations, and the new TL and AD models were derived for the discretized equations, i.e. using line-by-line approach. The purpose is to adequately describe moist
2.1 Formulation of moist dynamics model with aerosols 9 physical processes that include on-off switches (e.g. Zou, 1997). The new modelling frame- work has been accepted for publication in QJRMS (Zaplotnik et al., 2018).
As the existing system was developed for the Tropics, the present project has benefited from using the existing tropical framework for data assimilation. If any, this has made the addressed topic more challenging to work with. It is also easy to argue that one should focus on the Tropics where the wind-field information is missing the most.
There is hardly any moist simplified (with respect to NWP) model of complexity sim- ilar to the presented model. In particular, the model is based on the spectral methods that are used in the main European numerical weather prediction (NWP) global model, the European Centre for Medium-range Weather Forecasts (ECMWF) model and in a spectrum of mesoscale NWP models ALADIN/HIRLAM/AROME/HARMONIE that are applied operationally in the majority of European countries. The applied variational data assimilation approach closely follows the ideas and numerical solutions applied in the NWP models. After presenting the forecast model and variational data assimilation modelling, we conclude this chapter by an example of the newly developed data assimilation system.
2.1 Formulation of moist dynamics model with aerosols
The Moist Atmosphere Dynamics model (MAD) describes the dynamical response of the tropical atmosphere to the heat sources such as condensational heating and diabatic heating from thermal radiation (radiative cooling) or as a result of evaporative heating from the underlying ocean surface. The diabatic heating induces a dynamical response in form of large-scale travelling equatorially-trapped waves, which rapidly decay away from equator but can travel in the zonal direction along the equator for several thousand kilometres.
That means the diabatic heating affects circulation not only locally but it induces remote response through large-scale equatorial waves. MAD dynamical equations are based on the work of Gill (1982c,b); Heckley and Gill (1984); Davey and Gill (1987) and is a direct follow-on of ˇZagar et al. (2008) and ˇZagar (2012).
The model dynamics describes both balanced (vorticity dominated) dynamics and mainly divergent gravity-wave dynamics as well as their interactions. The prognostic model simulates also the moist atmospheric processes, i.e. non-linear moisture advection, conden- sation, evaporation and the impact of released latent heat on circulation. MAD is extended by a prognostic equation for a single mode (total) aerosol mass mixing ratio with simple parametrizations of dry and wet deposition. Especially the latter heavily influences the aerosol spatial distribution (Rasch et al., 2000) as it is the dominant sink of aerosols in the atmosphere. The model dynamics, physics and numerics are significantly simplified with respect to NWP models, however the model still shares the same basic dynamical properties and can therefore highlight most relevant aspects of the aerosol-moisture-wind feedback.
2.1.1 Equations for dynamics
The dynamical equations describe the potential temperature perturbations at the mid- troposphere level where the heating is maximal and associated horizontal motions in the
2.1 Formulation of moist dynamics model with aerosols 10 lower troposphere. In response to the half sinusoidal vertical structure of heating, the potential temperature perturbations and horizontal motions have half-sinusoidal and cos- inusoidal vertical structure, as illustrated in Figure 2.1. Therefore, motions in the upper troposphere have opposite signs, the vertical structure known as the first baroclinic mode.
To maintain the heat balance, the latent heating and the adiabatic cooling of the rising air must be appoximately balanced, and the vertical winds also need to attain half-sinusoidal vertical structure.
Figure 2.1: Simple model of tropical atmosphere following A. Gill.
Prognostic equations for the mid-level potential temperature θ and lower-level wind componentsv= (u, v) are
∂θ
∂t + (v· ∇)θ−θ0N2H
g (∇ ·v) =QLH +Q−θθ , (2.1)
∂v
∂t + (v· ∇)v+fk×v= gH θ0
∇θ−vv. (2.2)
Here, QLH denotes the latent heating whereas Q represents other diabatic forcings. Fric- tional processes for the mass and momentum are parameterized by terms. These are defined by choosing the Newtonian cooling relaxation time scale 1/θ and the Rayleigh- wind friction relaxation time scale 1/vequal 2 months. πHis the depth of the troposphere above the boundary layer, θ0 is the background potential temperature and N stands for the Brunt-V¨ais¨al¨a frequency, defined as N2 = (g/θ0)∂θ0/∂z. Other constants here have their usual meaning: ρ0 is the background air density, g is gravity and f is the Coriolis parameter. In the Tropics, the Coriolis parameter f is replaced by the equatorial β-plane approximation: f =βy, where β =∂f /∂y = 2.28·10−11 m−1s−1, andRe is the radius of the Earth.
2.1 Formulation of moist dynamics model with aerosols 11 2.1.2 Moist processes
MAD assumes exponentially decreasing vertical profile of specific humidity, which agrees well with the observed structure (Holloway and Neelin, 2009). Similar to what was done in Gill (1982b); Davey and Gill (1987), the moisture distribution can be expressed as q = ˜q(x, y, t) exp [−(z−zq)/H)], (2.3) where ˜q is the horizontal variation of the column mean specific humidity (in units kg of water vapour per kg air) and zq the vertical level of the mean humidity, thus q= ˜q(x, y, t) atz=zq. It is also assumed that the moisture flux is dominated by the bottom level winds, as the bulk of moisture is concentrated near the ground. The moisture flux error introduced by ignoring the vertical profile of winds and specific humidity and instead approximating the flux as a product of bottom level winds and column mean specific humidity is less than 15 %. The prognostic equation for the conservation of moisture (with tildes dropped) simulates the temporal evolution of specific humidity horizontal distribution at vertical levelz=zq, i.e. the distribution of column-mean specific humidity:
∂q
∂t +∇ ·(vq) =E−C. (2.4)
The parameters E and C represent the rate of evaporation and condensation. Condensa- tion occurs in the model when the atmosphere is saturated and there is a convergence of moisture. All excess moisture is precipitated. In the absence of evaporation, condensation rate can be expressed as
C=
(−∇ ·(vq) ifq≥qs
0 otherwise. (2.5)
Here, qs is saturation specific humidity, which is also assumed to have exponential ver- tical distribution. Condensation leads to the release of the latent heat QLH at the mid- tropospheric level
QLH = 1 2
Lc
cpC, (2.6)
where Lcis the latent heat of condensation and cp the specific heat of dry air at constant pressure. Factor 2 is included to relate the total column condensed moisture and total column released latent heat.
The system of prognostic equations includes some unwanted positive feedback process.
If condensation occurs due to large scale convergence in atmosphere near saturation, heat- ing is added, initiating even more convergence and thus even more heating, as depicted in Davey and Gill (1987). Rayleigh friction, Newtonian cooling and influx of surrounding drier air all partly (but not satisfactory) suppress this process. Here, saturated specific humidity is described as a function of temperature in order to damp the positive feedback.
Thus,
qs(T) = es(T) pq
Rd
Rv, (2.7)
wherees(T) is estimated from Clausis-Clapeyron equation andRdand Rv are specific gas constants for dry air and water vapour, respectively. Pressurepq =p(zq) is chosen constant
2.1 Formulation of moist dynamics model with aerosols 12 as the impact of pressure perturbations on saturation humidity is negligible. Temperature Tq=T(zq) is estimated from θq =θ(zq) =θq0+θ0q using Poisson equation
Tq =θq pq
p0
Rdcp
=θq0
1 + θ0 θ00
sin zq
H
pq p0
Rdcp
, (2.8)
where θ0 is the potential temperature perturbation, a prognostic variable. Mean potential temperature atzq is linearly extrapolated as θ0(zq) =θ00−(∂θ0/∂z) π2H−zq.
The presence of moisture does not change the gravity wave phase speed directly. It is condensational heating that reduces restoring bouyant forces in the moist zones, decreasing the amplitude of the waves and slowing them down, as shown analytically by (Gill, 1982b) using a linearized 1D SWM. Another important consequence of the bouyancy frequency value N is the meridional trapping scale for the circulation, i.e. the meridional distance at which the amplitude of the equatorial waves decays. The distance is described by the equatorial Rossby deformation radius, Re =pN H/(2β), which is around 1000 km. The associated time scale isT =Re/c≈6 hours. Therefore, the choice of reduced static stability increases the equatorial trapping and changes characteristics of propagating waves.
The fundamental moisture-dynamical coupling occurs in relation to the released heat- ing. Until the saturation is reached, presence of mid-troposphere heating will lead to lower-level convergence towards the heated region and vertical motions. The convergence leads to saturation and precipitation, which additionally intensifies convergence in the lower level. Associated temperature and wind perturbations move along the equator as well as outside the Tropics. Moisture driven condensation processes are thus the main source of nonlinearity in the model as illustrated in the next section.
2.1.3 Aerosol equation
The aerosol dynamics is described by the Eulerian advection equation with simplest possible parametrization of dry deposition and below-cloud wet deposition. The prognostic equation is written as
∂c
∂t +∇ ·(vc) =−λ c P−Kdc+S+−S−, (2.9) where the constant λ is the below-cloud wet deposition (hereafter, wet deposition), P is the precipitation rate (thus scavenging rate is Λ =λP),Kdis the dry deposition coefficient (a proxy for aerosol deposition velocity), S+ represents sources (e.g. aerosol mixing ratio forcing from the bottom boundary), while S− represents sinks. All mentioned coefficients are empirically determined. Thus, if the fluid is motionless and the precipitation rate constant, the aerosol mixing ratio decays exponentially with time.
The aerosols are, like the moisture, assumed concentrated in the lower troposphere.
Aerosol distribution is governed by internal and external processes. Internal processes include coagulation, condensation, adsorption/desorption, heterogeneous chemisty and nu- cleation mechanisms. Since this study is performed in an idealized framework and focuses on dynamical aspects of aerosol assimilation, only the dominant external aerosol processes are described and included in the aerosol prognostic equation. These processes involve
2.1 Formulation of moist dynamics model with aerosols 13 advection, dry deposition (gravitational sedimentation) and below-cloud wet deposition (particle scavenging by precipitation) (Seinfeld and Pandis, 2006; Colbeck and Lazaridis, 2014). In-cloud wet deposition is not parametrized as MAD model does not simulate clouds (liquid water). For simplicity, different aerosol species are not distinguished by their origin (e.g. dust, sea-salt etc.) or their size. The bin representation of the aerosol mass-size distribution function, often utilized in operational framework (Morcrette and Coauthors, 2009), is here not applied. Instead, the aerosol size spectra is described by a single variable, the total aerosol mass mixing ratio. However, aerosol physical and chemical processes af- fect distinct parts of the size spectra differently. For example, the value of the below-cloud scavenging rate for the same precipitation rate can vary even four orders of magnitude in different parts of the aerosol size spectra (Croft et al., 2009). Furthermore, it can vary up to two orders of magnitude within any aerosol size mode (nucleation, Aitken, accu- mulation and coarse). Nonetheless, simplified integral mass-size representation of aerosol in MAD model can be justified by the fact that our aim is only to better understand aerosol-moisture-wind coupling and their mutual dynamical feedbacks in the 4D-Var data assimilation and to address potential dynamical pitfalls of wind extraction from aerosols.
That means that we are interested in their couplings on the time scales of 12 to 48 hours.
2.1.4 MAD numerical formulation
MAD is formulated for a limited area on the sphere but it includes also options for thef−plane andβ−plane and the periodic boundary conditions. Prognostic equation are discretized by a spectral transform formulation using the Fourier series in both horizontal directions. The spectral representation of the fields in a limited-area domain with time- dependent lateral boundary conditions requires the extension zone, in order to ensure the periodicity (Haugen and Machenhauer, 1993). The elliptic truncation is applied to achieve homogeneous and isotropic spectral representation over the whole model domain. The sequence of steps involved in the numerical solutions of the forecast model is the following:
1. g−→ˆg=F(g)
2. ∂g∂x =F−1(ikˆg), ∂g∂y =F−1(ilˆg)
3. ∂g∂t =fg,∂x∂g,∂g∂y, . . .−→ c∂g∂t =F∂g∂t 4. ˆgn+1= ˆgn−1+ 2∆t
c∂g
∂t
n
5. Relaxation towards lateral boundary conditions.
6. 4th order implicit horizontal numerical diffusion.
7. Asselin filtering: ˆg˜n= ˆgn+p(ˆgn−1−2ˆgn+ ˆgn+1)
In the first step, the prognostic grid-point fields θ, u, v and q (in further explanation, g = [gij]i=1,...,Nx;j=1,...,Ny denotes each field) are represented as truncated Fourier wave series ˆg = [ˆgkl]k=0,...,Nk;l=0,...,Nl(k) (F denotes Fourier transform). Here Nx and Ny are
2.1 Formulation of moist dynamics model with aerosols 14 the numbers of zonal and meridional grid points, k and l are zonal and meridional wave numbers, while Nkand Nl(k) are the maximal numbers of Fourier waves in each direction.
The latter is dependent on k to satisfy the elliptic truncation criterion. Quadratic grid is used to avoid aliasing in the nonlinear advective terms, thus the shortest wavelength is rep- resented by 3 grid points, i.e. Nx'3(Nk+ 1). Then, the spatial derivatives are computed in spectral space (by multiplying each spectral component ˆgklbyik oril) and transformed to physical grid-point space (step 2). Temporal tendencies of prognostic variables (∂g/∂t) are computed in step 3 and transformed back to spectral space. The leapfrog time stepping is performed in step 4.
The relaxation of time-dependent fields towards time-dependent boundary fields is per- formed using the space-dependent relaxation functionα(x, y) = 1−tanh (r(x, y)/2), where r(x, y) is the distance from the grid point (x, y) to the boundary of the inner integration area, normalized by the width of the relaxation zone (step 5). The fourth-order implicit horizontal numerical diffusion is applied to the spectral coefficients to prevent an accumu- lation of energy at the smallest resolved scales (3∆x) during temporal integration (step 6).
Accumulation of energy occurs because the nonlinear interaction of two finitely truncated fields can only generate an equally truncated field, eventhough naturally those interactions can generate twice shorter waves. The natural dissipation of energy occurs at sub-milimeter scales, which is several orders of magnitude smaller than the resolution of NWP models.
Lastly, the Asselin time filtering is performed at step 7 to damp the computational mode arising due to leap-frog time scheme. The described numerical procedure is similar to the approach applied in the NWP models ALADIN and HIRLAM (e.g. Gustafsson, 1998) and its further details are available in previous studies with dry dynamics (ˇZagar et al., 2004a, 2008).
A weakness of the spectral method with the finite truncation is that it can produce neg- ative values for quantities which are physically positive definite, such as specific humidity q. The Gibbs oscillations can pose a problem for moisture at both ends, i.e. numerics can produce supersaturated humidity values and triggers “numerical rain”. Sharp gradients and discontinuities associated with precipitation and latent heat release can occur during the simulation. However, those were found to be of less importance in our experiments.
The treatment of moist processes in physical space involves the following steps:
1. At the end of time stepn:
• compute grid point specific humidity tendencies∂q/∂t.
• If tendency is positive (∂q/∂t > 0) and atmosphere already locally saturated (q ≥qs), then: condensation rateP =∂q/∂t, set ∂q/∂t= 0, specific humidity retains saturated value q=qs, else: P = 0.
2. At the beginning of time stepn+ 1:
• diagnose saturated specific humidity qs(T(θ)),
• relax the extra moisture (if q > qs) obtained from time integration and add it to condensation from the previous time step: C =C+ (q−qs)/∆t. Setq=qs,
• compute latent heating rate QLH =LcC/(2cp) and add it to potential temper- ature tendency∂θ/∂t.
2.1 Formulation of moist dynamics model with aerosols 15 2.1.5 Model setup
MAD model settings were chosen to resemble the real state of the tropical troposphere as closely as possible in a simplified model. The choice of the characteristic depth H = 4.5 km implies that the tropopause height HT is about πH ≈ 14.1 km above the bottom level HB at around 1000 m (≈ 900 hPa) (Table 2.1). The elevation of the tropopause is thus around 15 km, which slightly underestimates the tropical tropopause height in the atmoshpere. Tropospheric mid-levelHM, where condensational heating is strongest, is thus atHB+πH/2≈8 km, or roughly 400 hPa. Typical (background) potential temperature there isθ00= 333 K and the lapse ratedθ0/dz= 4 K km−1.
Moisture is assumed to be exponentially distributed between bottom level and mid- tropospheric level as depicted by Equation 2.3. Level of mean-column specific humidity is thus zq = 0.46 H = 2 km above the bottom-level or at an elevation around HB+zq = 3 km, roughly corresponding to pressurepq= 700 hPa. Due to exponential vertical profile of moisture (and thus condensation rate), the precipitation rate atz=zq is merely α= 0.42 of that at the bottom-level. Mean saturation specific humidity at zq is estimated using Equations 2.7,2.8 and is around 10.5 g kg−1, resembling well the real saturation humidity values at tropical 700 hPa pressure level.
The wet deposition process is the dominant aerosol removal process in the atmosphere, thus it strongly influences the global aerosol distributions (Rasch et al., 2000). In MAD model, the wet scavenging constant is λ = 1 mm−1. If precipitation rate is 1 mm hr−1, the wet scavenging rate equals Λ =λP ≈10−4 s−1. The scavenging rate is thus chosen in line with the calculated rates (Croft et al., 2009) used in ECHAM5-HAM model, a general circulation model (Roeckner et al., 2003) coupled to the Hamburg Aerosol Model (HAM) (Stier et al., 2005). Croft et al. (2009) reported scavenging rate values of order 10−4 s−1 for aerosols with geometric mean radius above 0.5 µm (Figure 2.2). These numbers were obtained by assuming the exponential Marshall-Palmer raindrop distribution, log-normal distribution of each aerosol size mode and by knowing the aerosol-droplet size-dependent collision efficiencies. According to Seinfeld and Pandis (2006), aerosols with radius larger than 0.5µm represent almost the entire aerosol volume and thus entire aerosol mass. Thus, the above choice of scavenging constant is also representative for aerosol mixing ratio.
2.1.6 Forecast experiments using MAD
In this subsection we present several examples of the tropical moist dynamics to demon- strate that the MAD model is an appropriate tool for studying the interactions between aerosol, moisture and dynamics. The examples include an adjustment to zonal wind per- turbation at equator and an adjustment to potential temperature perturbation (alias for discontinuousδ-shaped heating perturbation) in the ITCZ. The former is used to interpret the role of moisture in the adjustment process, while the latter is used to discuss the ef- fect of equatorial β-plane on channeling of IG waves along equator. Lastly, the impact of moisture on the aerosol via wet deposition is demonstrated.
2.1 Formulation of moist dynamics model with aerosols 16
Table 2.1: MAD model physical parameters.
Bottom-level elevation HB ≈1 km
Bottom-level dry air density ρ00 1 kg m−3
Characteristic depth H 4.5 km
Tropopause height HT =HB+πH ≈15 km
Tropospheric mid-level HM =HB+πH/2 ≈7 km
Lapse rate dθ0/dz 4 K km−1
Mid-level potential temperature θ00 333 K
Bouyancy frequency N 0.011 s−1
Horizontal phase speed c=N H 49 m s−1
Deformation radius Re=pc/(2β) 1000 km
Time scale T =Re/c 6 hours
Rayleigh friction parameters u,v 5·10−6 s−1
Newtonian cooling parameter θ 5·10−6 s−1
Latent heating rate Lc 2.5 MJ kg−1
Specific heat of dry air at constant pressure cp 1004 J kg−1 K−1
Specific gas constant of dry air Rd 287 J kg−1 K−1
Specific gas constant of water vapour Rv 461.5 J kg−1 K−1 Column-mean-humidity altitude (above HB) zq = 0.46 H 2 km
Column-mean-humidity potential temperature θq0 313 K
Column-mean-humidity pressure level pq 700 hPa
γ =δTq/δθ θθq0
00sin zHqppq
0
Rdcp
≈0.37 Precipitation coefficient at zq 1−exp(zq/H−1) 0.42
Saturation specific humidity qs ≈10.5 g kg−1
Below-cloud wet deposition constant λ 1 mm−1
Dry-deposition coefficient Kd /
2.1 Formulation of moist dynamics model with aerosols 17
Figure 2.2: Mean mass (solid lines) and number (dashed lines) below-cloud scavenging coefficients as a function of the geometric mean aerosol radius. Rainfall precipitation rates are shown in the legend. From Croft et al. (2009).
Adjustment to zonal wind perturbation at the equator
The impact of saturation on propagating wave disturbances is illustrated by comparing the adjustment process to the zonal wind perturbation in a nearly saturated atmosphere with an experiment in which the atmosphere is saturated only eastward from the initial perturbation.
The zonal wind perturbation is Gaussian-shaped with a maximal amplitude 2.5 m s−1at the equator and halfwidth 7◦ in a nearly saturated (RH = 99 %) atmosphere everywhere in the domain (denoted saturated case) or east of longtitude λ = 191◦ (denoted partly saturated case). The background is motionless with potential temperature constant across the domain and equal to θ00= 333 K.
Figure 2.3 presents the model response. In both saturated case (Figure 2.3a-c) and partly saturated case (Figure 2.3d-f), adjustment leads to IG waves, which are propagating zonally both eastward and westward. The dynamics associated with the eastward travelling IG (EIG) wave are different in each case. In the saturated case, the convergence associated with EIG wave almost immediately leads to condensation, precipitation and consequently latent heat release, warming the travelling perturbation and so reducing its amplitude. By comparing the position of the trailing edge of EIG wave at equator after 12 hours (Fig- ure 2.3c) to that in unsaturated (not shown) or partly saturated atmosphere (Figure 2.3f), it can be observed that convection and associated latent heat release slow down the equatorial waves, reducing their phase speed, similar to the observed convectively-coupled equatorial waves in the real atmosphere (Wheeler and Kiladis, 1999).
In the partly saturated case (Figure 2.3d-f), the EIG wave with convergent dynamics is
2.1 Formulation of moist dynamics model with aerosols 18
a) b) c)
d) e) f)
Figure 2.3: Adjustment to the zonal wind perturbation at the equator with amplitude 2.5 m s−1 in the atmosphere which is (a-c) everywhere near saturation (RH = 99 %) and (d-e) near saturation only in the region east ofλ = 195◦ at (a, d) 3 hours, (b, e) 6 hours, (c, f) 12 hours. Blue and red isolines denote negative and positive potential temperature perturbations with temperature spacing 0.05 K. Black isolines denote windspeed, spacing is 0.5 m s−1. Three-hour cumulative precipitation (mm) is shown by rainbow colorbar.
already meridionally elongated before it reaches the saturation area after 6 hours. As the wave travels eastward, the related convergence front maintains the condensation forming a narrow precipitation front at the leading edge similar to Pauluis et al. (2008). The released latent heat eventually damps the EIG wave amplitude (not shown). As a result, the associated convergence is too weak to trigger condensation.
Despite the vast differences in IG wave dynamics, the remaining balanced state (con- sisting predominantly of n= 1 Rossby wave) after 12 hours in the region of initial pertur- bation is almost the same in case of saturated atmosphere, partly saturated atmosphere and unsaturated atmosphere (not shown), with the amplitude of the balanced zonal wind perturbation differing by less than 10 %.
Moist adjustment to temperature perturbation in ITCZ
Figure 2.4 presents the model response to a Gaussian temperature perturbation centred at 8◦N with amplitude of 1 K and halfwidth 5◦ and 1.5◦ in the zonal and meridional directions, respectively. The response in the unsaturated case (Figure 2.4a-c) substantially differs from the nearly saturated case defined by RH = 99 % in the initial state (Figure 2.4d- h). In the unsaturated case, the perturbation excites IG waves which effectively disperse potential energy on the advective time scale. After 6 hours, there is a negligible cyclonic circulation (wind speed less than 0.1 m s−1) in the region of the initial perturbation.
In the saturated case, the positive temperature perturbation induces convergence in the lower layer towards the perturbed area. The converging winds advect moisture and lead to saturation and precipitation in the perturbed region (Figure 2.4d). The released
2.1 Formulation of moist dynamics model with aerosols 19 latent heat is warming the mid-troposphere (i.e. increasing temperature) and thus trig- gers further convergence in the lower level, initiating the positive feedback loop. For the selected parameters of the model, the process continues for about 6 hours. After 6 hours, cyclonic circulation (speeds around 2 m s−1) begins to establish whereas the precipitation continues (Figure 2.4e). The positive feedback loop is subsequently broken when the in- creasing saturation value due to heating induced temperature increase becomes too large (after 9 hours, not shown). This results in less moisture being condensed and less latent heat emitted. Thus, the convergence and so the moisture influx start to cease. The re- maining low-level heating perturbation then starts to travell westward slowly (compare the perturbation centers in at initial time in Figure 2.4d and after 36 hours in Figure 2.4h). I projects mainly on the equatorial (n= 1) Rossby wave. After 36 hours (Figure 2.4h), the teleconnected part of the original perturbation is seen south of the equator at the same longtitude.
The gravity wave front is much wider and has a larger amplitude in the saturated case as the latent heat increases temperature perturbation. The amplitude of the front decays approximately inversely with the distance from the perturbation (the statement would be exact in case of f-plane approximation), as the waves disperse on the horizontal plane.
The wave dispersion process is relatively slow as no significant diffusion and low friction have been applied in the experiment. Note also that the gravity wave front amplitude is larger on the meridional than on the zonal side of the front in Figures 2.4a,d because the initial perturbation is steeper in the meridional direction.
Figure 2.4c,f illustrates a lack of zonal symmetry in 12-hour forecast. The asymmetry is absent in thef-plane experiment with f = 10−5 s−1 as the Kelvin wave is not an eigen- solution in this case (figures not shown). As the radiated gravity waves reach the equator (Figure 2.4f), the equatorial chanelling effect begins to act, directing the wave movement along the equator (shown in Figure 2.4g,h at 24 and 36 hours, respectively). Now, the front decay is not inversely proportional to the distance from perturbation origin anymore.
Theoretically, in the absence of humidity and significant friction or diffusion, the equato- rially trapped waves would persist infinetely. However, due to applied friction Presented numerical result is in line with the 3D analytical study by Peagle (1978), who demonstrated that tropical latent heating events with precipitation rates of a few centimeters per day can produce long lasting gravity waves with cross-isobaric flow of order of 1 m s−1 outward to beyond 2000 km.
Aerosol-moisture-wind coupling
The impact of moisture on the aerosol via wet deposition is demonstrated by an example of a vastly inhomogeneous humidity and aerosol fields subject to advection by a meandering low-level tropical easterly jet with speeds up to 20 m s−1 (Figure 2.5).
In the presented domain, the aerosol mixing ratio varies from 0 to 22µg kg−1 and the horizontal mixing ratio gradients are strong. On the other hand, specific humidity field is more homogeneous, with values ranging between 6 and 10.5 g kg−1. The meandering easterly winds lead to the development of barotropic eddies. Precipitation occurs in relation to the flow convergence in saturated areas such as in Figure 2.5b,c at around λ = 210◦,
2.1 Formulation of moist dynamics model with aerosols 20
a) b) c)
d) e) f)
g) h)
Figure 2.4: Adjustment to the potential temperature perturbation in (a-c) unsaturated atmosphere, (d- h) saturated atmosphere at (a,d) 3 hours, (b,e) 6 hours, (c,f) 12 hours, (g) 24 hours and (h) 36 hours.
Blue-red colourbar denote negative and positive potential temperature perturbations with spacing 0.05 K.
Three-hour cumulative precipitation (mm) is shown by the rainbow colourbar.
φ = -5◦) or as a result of lower level inflow initiated by some heating source such as in Figure 2.5b,c at aroundλ= 187◦,φ= -5◦. In the areas of precipitation, aerosol mixing ratio exponentially decreases, as the aerosol is scavenged by precipitation. Notice in Figure 2.5d,e the white areas of low aerosol mixing ratio at λ= 187◦ at equator, at λ= 210◦,φ = -10◦ and at λ= 220◦,φ= 13◦.
The experiment represents how the aerosol prognostic equation coupled to MAD fore- cast model describes the temporal and spatial horizontal distribution of aerosol with an adequate degree of realism.
2.1 Formulation of moist dynamics model with aerosols 21 a)
b)
c)
d)
e)
f)
g)
h)
i)
j)
Figure 2.5: 24-hour simulation of aerosol mixing ratio (a-e), specific humidity field (f-j), wind field (all) and 6-hour cumulative precipitation (mm), denoted by rainbow colorbar; (a,f) at initial time, (b,g) at 6, (c,h) 12, (d,i) 18 and (e,j) at 24 hours.
2.2 Data assimilation modelling 22
2.2 Data assimilation modelling
In this section we describe the new 4D-Var assimilation with moisture and aerosol assimilation. The variational assimilation problem is formulated using the incremental approach of Courtier et al. (1994) to account for the nonlinearities due to moist processes.
The assimilation system including the moist non-linear forecast model with aerosols, and its TL and AD models constitutes our Moist Aerosol Dynamics Data Assimilation Model
−MADDAM.
2.2.1 Four-dimensional variational data assimilation
4D-Var seeks a solution (analysisxa), which minimizes the distance of the 4D trajectory to observationsyk, measured by cost function termJo, and background (previous forecast) statexb, measured by Jb. The cost function is therefore
J(x) =Jb+Jo
= 1
2(xb−x)TB−1(xb−x) +1
2
N
X
k=0
(yk−Gk(x))TR−1k (yk−Gk(x))
(2.10)
and is minimized for x = xa. Here, B is the background error covariance matrix, Rk = R are temporally invariant covariance matrices of the observation errors. Observations errors are assumed statistically independent, thereforeR becomes diagonal matrix of error covariances. The nonlinear generalised observation operator Gk =HkMt0→tk produces the model equivalent of observations by integrating state x from time t0 to time tk (denoted by an operator Mt0→tk) to obtain model state xk at time tk. For simplicity we simulated observations at model grid points, thus the observation operator isHk= 1.
The minimization problem (2.10) is solved using the incremental approach (Courtier et al., 1994) which approximated the nonlinear cost function by a sequence of quadratic cost functions around the latest guessxg. This sequence is most often denoted as the outer loop.
Each quadratic cost function is then minimized (so-called inner loop) to update the guess for the next outer loop iteration. The minimization problem is simplified by transforming the state increment δx to a new control variable χ, which transforms B matrix into an identity matrix.
x−xg=δx=Lχ , (2.11)
where L is defined such that LTB−1L = I. Similarly, the guess increment δxg is trans- formed into new variableχg, i.e. δxg =xg−xb=Lχg. The minimization problem therefore becomes
J(χ) = 1
2(χ+χg)T(χ+χg) +1
2
N
X
k=0
(dk−GkLχ)TRk−1(dk−GkLχ),
(2.12)
