Restriction on the timestep to be used in stochastic Lagrangian models of turbulent dispersion

By considering two analytical solutions of G. I. Taylor (1921) for dispersion in homogeneous turbulence, we derive a quantitative upper limit for the timestep dt to be used in the stochastic Lagrangian model; a more severe upper limit will probably exist in inhomogeneous turbulence. For practical purposes, there is no lower limit to the timestep.     

Instability in Lagrangian stochastic trajectory models,and a method for its cure

A number of authors have reported the problem of unrealistic velocities (“rogue trajectories”) when computing the paths of particles in a turbulent flow using modern Lagrangian stochastic (LS) models, and have resorted to ad hoc interventions. We suggest that this problem stems from two causes: (1) unstable modes that are intrinsic to the dynamical system constituted by the generalized Langevin equations, and whose actual triggering (expression) is conditional on the fields of the mean velocity and Reynolds stress tensor and is liable to occur in complex, disturbed flows (which, if computational, will also be imperfect and discontinuous); and, (2) the “stiffness” of the generalized Langevin equations, which implies that the simple stochastic generalization of the Euler scheme usually used to integrate these equations is not sufficient to keep round-off errors under control. These two causes are connected, with the first cause (dynamical instability) exacerbating the second (numerical instability); removing the first cause does not necessarily correct the second, and vice versa. To overcome this problem, we introduce a fractional-step integration scheme that splits the velocity increment into contributions that are linear (U _i ) and nonlinear (U _i U _j ) in the Lagrangian velocity fluctuation vector U, the nonlinear contribution being further split into its diagonal and off-diagonal parts. The linear contribution and the diagonal part of the nonlinear contribution to the solution are computed exactly (analytically) over a finite timestep Δt, allowing any dynamical instabilities in the system to be diagnosed and removed, and circumventing the numerical instability that can potentially result in integrating stiff equations using the commonly applied explicit Euler scheme. We contrast results using this and the primitive Euler integration scheme for computed trajectories in a drastically inhomogeneous urban canopy flow.     

Simulation of the global bio-geophysical interactions during the Last Glacial Maximum

 The bio-geophysical feedbacks during the Last Glacial Maximum (LGM, 21 000 y BP) are investigated by use of an asynchronously coupled global atmosphere-biome model. It is found that the coupled model improves on the results of an atmosphere-only model especially for the Siberian region, where the inclusion of vegetation-snow-albedo interaction leads to a better agreement with geological reconstructions. Furthermore, it is shown that two stable solutions of the coupled model are possible under LGM boundary conditions. The presence of bright sand desert at the beginning of a simulation leads to more extensive subtropical deserts, whereas an initial global vegetation cover with forest, steppe, or dark desert results in a northward spread of vegetation of up to some 1000 km, mainly in the western Sahara. These differences can be explained in the framework of Charney’s theory of a “self-induction” of deserts through albedo enhancement. Moreover, it is found that the tropical easterly jet is strengthened in the case of the “green” Sahara, which in turn leads to a modification of the Indian summer monsoon. Received: 20 June 1997/Accepted: 14 January 1998     

On the Formation of an Elevated Nocturnal Inversion Layer in the Presence of a Low-Level Jet: A Case Study

We report on observed nocturnal profiles, in which an inversion layer is located at the core of a low-level jet, bounded between two well-mixed layers. High-resolution vertical profiles were collected during a field campaign in a small plain in the Israeli desert (Negev), distant 100 km from the eastern shore of the Mediterranean Sea. During the evening hours, the synoptic flow, superposed on the late sea breeze, forms a low-level jet characterized by a maximum wind speed of 12 m s ⁻¹ at an altitude of 150 m above the ground. The strong wind shear at the jet maximum generates downward heat fluxes that act against the nocturnal ground cooling. As a result, the typical ground-based nocturnal inversion is “elevated” towards the jet centre, hence a typical early morning thermal profile is observed a few hours after sunset. Since the jet is advected into the region, its formation does not depend on the presence of a surface nocturnal inversion layer to decouple the jet from surface friction. On the contrary, here the advected low-level jet acts to hinder the formation of such an inversion. These unusual temperature and wind profiles are expected to affect near-ground dispersion processes.     

A revised method of presenting wavenumber–frequency power spectrum diagrams that reveals the asymmetric nature of tropical large-scale waves

The popular method of presenting wavenumber–frequency power spectrum diagrams for studying tropical large-scale waves in the literature is shown to give an incomplete presentation of these waves. The so-called “convectively coupled Kelvin (mixed Rossby-gravity) waves” are presented as existing only in the symmetric (anti-symmetric) component of the diagrams. This is obviously not consistent with the published composite/regression studies of “convectively coupled Kelvin waves,” which illustrate the asymmetric nature of these waves. The cause of this inconsistency is revealed in this note and a revised method of presenting the power spectrum diagrams is proposed. When this revised method is used, “convectively coupled Kelvin waves” do show anti-symmetric components, and “convectively coupled mixed Rossby-gravity waves (also known as Yanai waves)” do show a hint of symmetric components. These results bolster a published proposal that these waves should be called “chimeric Kelvin waves,” “chimeric mixed Rossby-gravity waves,” etc. This revised method of presenting power spectrum diagrams offers an additional means of comparing the GCM output with observations by calling attention to the capability of GCMs to correctly simulate the asymmetric characteristics of equatorial waves.     

Sensitivity of the thermohaline circulation and climate to ocean exchanges in a simple coupled model

 A new simple, coupled climate model is presented and used to investigate the sensitivity of the thermohaline circulation and climate to ocean vertical and horizontal exchange. As formulated, the model highlights the role of thin, ocean surface layers in the communication between the atmosphere and the subsurface ocean. Model vertical exchange is considered to be an analogue to small-scale, diapycnal mixing and convection (when present) in the ocean. Model horizontal exchange is considered to be an analogue to the effects of the wind-driven circulation. For small vertical exchange in the ocean, the model exhibits only one steady-state solution: a relatively cold, mid-high-latitude climate associated with a weak, salinity-driven circulation (“off ” mode). For large vertical and horizontal exchange in the ocean, the model also exhibits only one steady-state solution: a relatively warm, mid-high-latitude climate associated with a strong, thermally-driven circulation (“on” mode). For sufficiently weak horizontal exchange but large enough vertical exchange, both modes are possible stable, steady-state solutions. When model parameters are calibrated to fit tracer distributions of the modern ocean-atmosphere system, only the “on” mode is possible in this standard case. This suggests that the wind-driven circulation in consort with diapycnal mixing suppresses the “off ” mode in the modern ocean-atmosphere system. Since both diapycnal mixing and the wind-driven circulation would be expected to increase in a cold climate with greater meridional temperature gradients and enhanced winds, vertical and horizontal exchange in the ocean are probably associated with strong negative feedbacks which tend to stabilize climate. These results point to the need to resolve ocean wind-driven circulation and to greatly improve the treatment of ocean diapycnal mixing in more complete models of the climate system. Received: 16 November 1999 / Accepted: 19 June 2000     

The coupling of optimal economic growth and climate dynamics

In this paper, we study optimal economic growth programs coupled with climate change dynamics. The study is based on models derived from MERGE, a well established integrated assessment model (IAM). We discuss first the introduction in MERGE of a set of “tolerable window” constraints which limit both the temperature change and the rate of temperature change. These constraints, obtained from ensemble simulations performed with the Bern 2.5-D climate model, allow us to identity a domain intended to preserve the Atlantic thermohaline circulation. Next, we report on experiments where a two-way coupling is realized between the economic module of MERGE and an intermediate complexity “3-D-” climate model (C-GOLDSTEIN) which computes the changes in climate and mean temperature. The coupling is achieved through the implementation of an advanced “oracle based optimization technique” which permits the integration of information coming from the climate model during the search for the optimal economic growth path. Both cost-effectiveness and cost-benefit analysis modes are explored with this combined “meta-model” which we refer to as GOLDMERGE. Some perspectives on future implementations of these approaches in the context of “collaborative” or “community” integrated assessment modules are derived from the comparison of the different approaches.     

Closure of the Panama Seaway during the Pliocene: implications for climate and Northern Hemisphere glaciation

The “Panama Hypothesis” states that the gradual closure of the Panama Seaway, between 13 million years ago (13 Ma) and 2.6 Ma, led to decreased mixing of Atlantic and Pacific water Masses, the formation of North Atlantic Deep water and strengthening of the Atlantic thermohaline circulation, increased temperatures and evaporation in the North Atlantic, increased precipitation in Northern Hemisphere (NH) high latitudes, culminating in the intensification of Northern Hemisphere Glaciation (NHG) during the Pliocene, 3.2–2.7 Ma. Here we test this hypothesis using a fully coupled, fully dynamic ocean-atmosphere general circulation model (GCM) with boundary conditions specific to the Pliocene, and a high resolution dynamic ice sheet model. We carry out two GCM simulations with “closed” and “open” Panama Seaways, and use the simulated climatologies to force the ice sheet model. We find that the models support the “Panama Hypothesis” in as much as the closure of the seaway results in a more intense Atlantic thermohaline circulation, enhanced precipitation over Greenland and North America, and ultimately larger ice sheets. However, the volume difference between the ice sheets in the “closed” and “open” configurations is small, equivalent to about 5 cm of sea level. We conclude that although the closure of the Panama Seaway may have slightly enhanced or advanced the onset of NHG, it was not a major forcing mechanism. Future work must fully couple the ice sheet model and GCM, and investigate the role of orbital and CO₂ effects in controlling NHG.     

“Initialization” using an iterative Matsuno style scheme in the Eta Model adjustment stage

Summary  Following a recent approach of Fox-Rabinovitz an iterative Matsuno or a “Super-Matsuno” style scheme is applied as a filter in the Eta Model. In contrast to Fox-Rabinovitz, we however apply the scheme not for all of the model’s time-differencing but for its adjustment terms only. These distinctions compared to the original Fox-Rabinovitz’-method are made easy to implement by the split time differencing approach of the Eta, and at the same time would appear clearly appropriate for the “initialization” purpose. In addition, while Fox-Rabinovitz emphasizes the use of the method within a long time-scale data assimilation framework, we are focusing on the impact of the method in a short-range forecasting environment/time-scale. After a short one hour “initialization” procedure is completed, standard model integration is continued, now very much free of noise. The Super-Matsuno style scheme is found to balance initially unbalanced external and internal modes and to significantly reduce the high-frequency noise during the first 6 time steps. In a control case noise also reduces in amplitude as integration proceeds, but at a much slower rate. The model integration results with and without “initialization” after 6 hours are however very similar. Even so, it is to be expected that small differences, given that they have resulted from the removal of spurious initial noise, have to be beneficial. Received July 5, 1999/Revised January 11, 2000     

Numerical Considerations for Lagrangian Stochastic Dispersion Models: Eliminating Rogue Trajectories,and the Importance of Numerical Accuracy

When Lagrangian stochastic models for turbulent dispersion are applied to complex atmospheric flows, some type of ad hoc intervention is almost always necessary to eliminate unphysical behaviour in the numerical solution. Here we discuss numerical strategies for solving the non-linear Langevin-based particle velocity evolution equation that eliminate such unphysical behaviour in both Reynolds-averaged and large-eddy simulation applications. Extremely large or ‘rogue’ particle velocities are caused when the numerical integration scheme becomes unstable. Such instabilities can be eliminated by using a sufficiently small integration timestep, or in cases where the required timestep is unrealistically small, an unconditionally stable implicit integration scheme can be used. When the generalized anisotropic turbulence model is used, it is critical that the input velocity covariance tensor be realizable, otherwise unphysical behaviour can become problematic regardless of the integration scheme or size of the timestep. A method is presented to ensure realizability, and thus eliminate such behaviour. It was also found that the numerical accuracy of the integration scheme determined the degree to which the second law of thermodynamics or ‘well-mixed condition’ was satisfied. Perhaps more importantly, it also determined the degree to which modelled Eulerian particle velocity statistics matched the specified Eulerian distributions (which is the ultimate goal of the numerical solution). It is recommended that future models be verified by not only checking the well-mixed condition, but perhaps more importantly by checking that computed Eulerian statistics match the Eulerian statistics specified as inputs.     

On the moments approximation method for constructing a Lagrangian Stochastic model

The exact Eulerian velocity probability density function (pdf) of a turbulent field is generally unknown, and one normally has available only partial information in the form of low order moments. We compare two alternative Lagrangian Stochastic (LS) approaches formed from this partial information, (i) the moments approximation approach (Kaplan and Dinar, 1993); and (ii) the well-mixed model (Thomson, 1987) that corresponds to the maximum missing information pdf formed from the available information. We show that the moments approximation model does not in general satisfy the well-mixed constraint, and can give an inferior prediction of dispersion.     

Numerical convergence of the dynamics of a GCM

 Atmospheric general circulation models (GCMs) are characterized by many features but especially by: (1) the manner of discretizing the governing equations and of representing the variables involved at a given resolution, and (2) the manner of parameterizing unresolved physical processes in terms of those resolved variables. These two aspects of model formulation are not independent and it is difficult to untangle their intertwined effects when assessing model performance. The attempt here is to separate these aspects of GCM behaviour and to ask, “Given a perfect parameterization of the physical processes in a model, what resolution is needed to capture the dominant dynamical aspects of the atmospheric climate?” Alternatively, “At what resolution do the dynamics of a GCM converge”? The perfect parameterization approach assumes that the calculation of the physical terms returns the “correct” result at all resolutions. In the idealized case, a time-independent forcing is one of the simplest that satisfies this condition. However, experiments show that it is difficult for the dynamics of a GCM to balance a time-independent forcing with atmosphere-like flows and structures. The model requires, and the atmosphere presumably includes, physical feedback mechanisms which act so as to maintain the kinds of flows and structures that are observed. Resolution experiments are performed with a simplified forcing function for the thermodynamic equation which combines a dominant time-independent specified forcing with a weak linear relaxation feedback. These experiments show that the dynamics of the GCM have essentially converged at T32 and certainly by T63 which is the next resolution considered. This is shown by the constancy of structures, variances, covariances, transports and energy budgets with increasing resolution. Experiments with an alternative forcing proposed by Held and Suarez, which has the form of a linear relaxation, show somewhat less evidence of convergence at these resolutions. In both cases the “physics” are known by assumption. However, the form and nature of the forcing is different, as is the behaviour with resolution. The implication for the real system is that the resolution required for simulating the dynamical aspects of climate is rather modest. The simulated climate does, however, apparently depend on the ability to correctly and consistently parameterize the physical processes in a GCM, involving both forcing and feedback mechanisms, as a function of resolution. Received 19 January 1996/Accepted 22 August 1996     

Flow splitting at the bifurcation between two valleys: idealized simulations in comparison with Mesoscale Alpine Programme observations

Summary This paper investigates the characteristics of channelled airflow in the vicinity of a junction of three idealized valleys (one valley carrying the incoming flow and two tributaries carrying the outflow), using a two-dimensional single-layer shallow water model. Particular attention is given to the flow splitting occurring at the junction. Nondimensionalized, the model depends on the valley geometry, the Reynolds number, which is related to the eddy viscosity, and on the difference of the hydrostatic pressure imposed at the exit of the tributaries. At the spatial scale considered in this study, the Rossby number relating the inertial and Coriolis forces is always larger than 1, implying that the effect of earth rotation can be neglected to a first approximation. The analysis of the flow structure within the three valleys as well as the calculation of the split ratio (fraction of the air flow diverted into one of the two downstream valleys with respect to the total mass flux in the upstream valley) show that (i) the flow pattern depends strongly on the Reynolds number while the split ratio is comparatively insensitive; (ii) the valley geometry and the difference between the upstream and downstream hydrostatic pressures affect the flow pattern, the location of the split point and the split ratio; (iii) the relative contribution of flow deflection by the sidewalls and the blocking/splitting mechanism differs between the settings of a “Y-shape” valley and a “T-shape” valley. Quantitative comparison of the present results with numerical simulations of realistic cases and with observations collected in the region of the Rhine and Seez valleys (Switzerland) (“Y-shape” valley) and in the region of the Inn and Wipp valleys (Austria) (“T-shape” valley) during the Mesoscale Alpine Programme (MAP) field experiment shows good agreement provided that the normalized valley depth NΔH/U_u significantly exceeds 1, i.e., when “flow around” is expected. A structural disagreement between the idealized simulations and the observed wind field is found only when NΔH/U_u ≃ 1, that is, in the “flow over” regime. This shows that the dimensionless valley depth is indeed a good indicator for flow splitting, implying that the stratification is a key player in reality.     

Multi-Source Emission Determination Using an Inverse-Dispersion Technique

Inverse-dispersion calculations can be used to infer atmospheric emission rates through a combination of downwind gas concentrations and dispersion model predictions. With multiple concentration sensors downwind of a compound source (whose component positions are known) it is possible to calculate the component emissions. With this in mind, a field experiment was conducted to examine the feasibility of such multi-source inferences, using four synthetic area sources and eight concentration sensors arranged in different configurations. Multi-source problems tend to be mathematically ill-conditioned, as expressed by the condition number κ. In our most successful configuration (average κ = 4.2) the total emissions from all sources were deduced to within 10% on average, while component emissions were deduced to within 50%. In our least successful configuration (average κ = 91) the total emissions were calculated to within only 50%, and component calculations were highly inaccurate. Our study indicates that the most accurate multi-source inferences will occur if each sensor is influenced by only a single source. A “progressive” layout is the next best: one sensor is positioned to “see” only one source, the next sensor is placed to see the first source and another, a third sensor is placed to see the previous two plus a third, and so on. When it is not possible to isolate any sources κ is large and the accuracy of a multi-source inference is doubtful.     

Bayesian multi-model projection of climate: bias assumptions and interannual variability

Current climate change projections are based on comprehensive multi-model ensembles of global and regional climate simulations. Application of this information to impact studies requires a combined probabilistic estimate taking into account the different models and their performance under current climatic conditions. Here we present a Bayesian statistical model for the distribution of seasonal mean surface temperatures for control and scenario periods. The model combines observational data for the control period with the output of regional climate models (RCMs) driven by different global climate models (GCMs). The proposed Bayesian methodology addresses seasonal mean temperatures and considers both changes in mean temperature and interannual variability. In addition, unlike previous studies, our methodology explicitly considers model biases that are allowed to be time-dependent (i.e. change between control and scenario period). More specifically, the model considers additive and multiplicative model biases for each RCM and introduces two plausible assumptions (“constant bias” and “constant relationship”) about extrapolating the biases from the control to the scenario period. The resulting identifiability problem is resolved by using informative priors for the bias changes. A sensitivity analysis illustrates the role of the informative prior. As an example, we present results for Alpine winter and summer temperatures for control (1961–1990) and scenario periods (2071–2100) under the SRES A2 greenhouse gas scenario. For winter, both bias assumptions yield a comparable mean warming of 3.5–3.6°C. For summer, the two different assumptions have a strong influence on the probabilistic prediction of mean warming, which amounts to 5.4°C and 3.4°C for the “constant bias” and “constant relation” assumptions, respectively. Analysis shows that the underlying reason for this large uncertainty is due to the overestimation of summer interannual variability in all models considered. Our results show the necessity to consider potential bias changes when projecting climate under an emission scenario. Further work is needed to determine how bias information can be exploited for this task.     

A closure to derive a three-dimensional well-mixed trajectory-model for non-Gaussian, inhomogeneous turbulence

A three-dimensional Lagrangian stochastic (LS) model to evaluate pollutant dispersion in the atmospheric boundary layer has been developed. The model satisfies the well-mixed criterion of Thomson and allows for inhomogeneous, skew turbulence. Making use of the spherical reference frame, one of the possible solutions has been obtained. A skewed joint probability density function (PDF), which reproduces the given velocity moments (means, variances, skewness and covariances), has been built-up by a linear combination of eight Gaussian PDFs. In order to verify consistency with the well-mixed criterion, the long term results have been compared with the theoretical behaviour. A comparison between our model and Thomson's published algorithms was also carried out. By comparing wind-tunnel data and numerical predictions, a further validation of our LS model has been obtained. From an analysis of the numerical results, we can state that our model is able to evaluate dispersion in the case of complex flows where the application of previous models is unsuccessful.     

An Evaluation of the Flux–Gradient Relationship in the Stable Boundary Layer

Data collected during the SHEBA and CASES-99 field programs are employed to examine the flux–gradient relationship for wind speed and temperature in the stably stratified boundary layer. The gradient-based and flux-based similarity functions are assessed in terms of the Richardson number Ri and the stability parameter z/Λ_*, z being height and Λ_* the local Obukhov length. The resulting functions are expressed in an analytical form, which is essentially unaffected by self-correlation, when thermal stratification is strong. Turbulence within the stably stratified boundary layer is classified into four regimes: “nearly-neutral” (0 < z/Λ_* < 0.02), “weakly-stable” (0.02 < z/Λ_* < 0.6), “very-stable” (0.6 < z/Λ_* < 50), and “extremely-stable” (z/Λ_* > 50). The flux-based similarity functions for gradients are constant in “nearly-neutral” conditions. In the “very-stable” regime, the dimensionless gradients are exponential, and proportional to (z/Λ_*)^3/5. The existence of scaling laws in “extremely-stable” conditions is doubtful. The Prandtl number Pr decreases from 0.9 in nearly-neutral conditions and to about 0.7 in the very-stable regime. The necessary condition for the presence of steady-state turbulence is Ri < 0.7.     

Turbulence Scale Dependence of the Richardson Constant in Lagrangian Stochastic Models

We investigate the relative dispersion properties of the well-mixed class of Lagrangian stochastic models. Dimensional analysis shows that, given a model in the class, its properties depend solely on a non-dimensional parameter, which measures the relative weight of Lagrangian-to-Eulerian scales. This parameter is formulated in terms of Kolmogorov constants, and model properties are then studied by modifying its value in a range that contains the experimental variability. Large variations are found for the quantity, g* = 2gC₀^{− 1}, where g is the Richardson constant.     

Probabilistic Model for Concentration Fluctuations in Compact-Source Plumes in an Urban Environment

A comprehensive model for the prediction of concentration fluctuations in plumes dispersing in the complex and highly disturbed wind flows in an urban environment is formulated. The mean flow and turbulence fields in the urban area are obtained using a Reynolds-averaged Navier-Stokes (RANS) flow model, while the standard k-ϵ turbulence model (k is the turbulence kinetic energy and ϵ is the viscous dissipation rate) is used to close the model. The RANS model provides a specification of the velocity statistics of the highly disturbed wind flow in the urban area, required for the solution of the transport equations for the mean concentration and concentration variance (both of which are formulated in the Eulerian framework). A physically-based formulation for the scalar dissipation time scale t _d , required for the closure of the transport equation for , is presented. This formulation relates t _d to an inner time scale corresponding to “internal” concentration fluctuation associated with relative dispersion, rather than an outer time scale associated with the entire portion of the fluctuation spectrum. The two lowest-order moments of concentration ( and ) are used to determine the parameters of a pre-chosen functional form for the concentration probability density function (clipped-gamma distribution). Results of detailed comparisons between a water-channel experiment of flow and dispersion in an idealized obstacle array and the model predictions for mean flow, turbulence kinetic energy, mean concentration, concentration variance, and concentration probability density function are presented.     

Regime Diagrams for <Emphasis Type="Italic">K</Emphasis>-Theory Dispersion

In atmospheric dispersion, the “non-Gaussian” effects of gravitational settling, the vertical gradient in diffusivity and the surface deposition do not enter uniformly but rather break up parameter space into several discrete regimes. Here, we describe regime diagrams that are constructed for K-theory dispersion of effluent from a surface line source in unsheared inhomogeneous turbulence, using a previously derived Fourier–Hankel method. This K-theory formulation differs from the traditional one by keeping a non-zero diffusivity at the ground. This change allows for turbulent exchange between the canopy and the atmosphere and allows new natural length scales to emerge. The axes on the regime diagrams are non-dimensional distance defined as the ratio of downwind distance to the characteristic length scale for each effect. For each value of the ratio of settling speed to the K gradient, two to four regimes are found. Concentration formulae are given for each regime. The regime diagrams allow real dispersion problems to be categorized and the validity of end-state concentration formulae to be judged.