Dimension of aircraft exhaust plumes at cruise conditions : effect of wake vortices

The dispersion of aircraft emissions during the vortex phase is studied using a 3-D LES model with Lagrangian particle tracking. The simulations start with a fully rolled-up vortex pair of a type B747/A340 airplane and the tracer centred around the vortex cores. The tracer dilution and plume extent is studied for a variety of ambient and aircraft parameters until aircraft-induced effects have ceased. For typical upper tropospheric conditions, the impact of stratification is more dominant compared to turbulence intensity or vertical wind shear. Moreover, the sensitivity to the initial tracer distribution was found to be weak. Along the transverse direction, the tracer concentrations can be well approximated by a Gaussian distribution, along the vertical a superposition of three Gaussian distributions is adequate. For the studied parameter range, the vertical plume expansion ranges from 400 m to 550 m and cross-sectional area from 4.0× 104 m2 to 6× 104 m2 after six minutes. For validation, selected simulations were compared to an alternative LES model and to in-situ NO-measurements.


Introduction
From the wake-vortex dynamics point of view, it is convenient to distinguish three phases in the life cycle of an aircraft wake (Hoshizaki et al., 1975): During the jet regime the vorticity sheet around the wings rolls up into two counterrotating vortices.The vortex pair then propagates vertically downward by mutual velocity induction in the vortex regime.In the dissipation regime the organised vortices have bro-ken up in turbulence and the energy dissipates to the background level.All three phases have a distinct influence on the distribution of the exhaust from the turbines.In the jet regime, the exhaust jets mix with ambient air while they are entrained into the forming vortices.When the vortex pair descends most of the exhaust is stored in that primary wake and, hence, also sinks below flight level.However, as the primary wake moves downward in an atmosphere which is typically stably stratified at cruise flight levels, it produces a wake itself, the so-called secondary wake.Adiabatic heating leads to a baroclinic force which detrains some of the exhaust and leaves a curtain of exhaust above the descending vortex pair (Gerz et al., 1998;Holzäpfel et al., 2003;Misaka et al., 2012).When the organised vortices break up the exhaust concentrated in the primary wake is suddenly released and diffuses quickly.
The evolution of aircraft plumes during the vortex phase is important to study for several reasons.The descending vortex pair increases the vertical extent of the plume, which in turn affects the horizontal spreading by the vertically sheared (cross-) wind.Regarding the long-term evolution over many hours this initial expansion can be viewed as a "sudden" event since the time scale of (natural) vertical mixing in the free troposphere is too large to produce such deep plumes within reasonable time.Since the horizontal spreading by vertical wind shear scales linearly with the plume height, the vertical plume extension after the vortex break-up is crucial to know.
This study looks in detail at the plume mixing under a variety of atmospheric conditions and aircraft parameters.The S. Unterstrasser et al.: Dimension of aircraft plumes: effect of wake vortices presented results have implications, both for nonlinear plume chemistry and for the microphysical evolution in a contrail.This work bridges the gap between analytical jet expansion laws describing the plume mixing within the first seconds behind the aircraft and the later dispersion in the free atmosphere after aircraft-induced effects have vanished.For box model studies of the jet phase, often treating the contrail formation, analytical jet expansion laws are used (Kärcher and Yu, 2009;Wong and Miake-Lye, 2010).Especially for chemistry applications, box models treat the plume evolution also beyond the end of the jet phase (Miake-Lye et al., 1993;Danilin et al., 1994;Kärcher, 1995;Kärcher et al., 1996).Plume dilution over longer time scales (larger than O(1 h)) is often based on empirical laws derived from in situ NOmeasurements (Schumann et al., 1998) or on results from Large-Eddy Simulations (LES) (Dürbeck and Gerz, 1995;Dörnbrack and Dürbeck, 1998).
In the later stages the plume is controlled by vertical wind shear and associated turbulent mixing.The ratio of vertical and horizontal spreading also depends on stratification.The intermediate regime with the interaction of the plume and the wake vortices is more intricate.
In larger-scale models individual plumes cannot be resolved on the grid scale and nonlinear chemical effects within aircraft plumes (Vohralik et al., 2008;Vancassel et al., 2010) cannot be directly simulated.Commonly, effective emissions are introduced in order to take into account nonlinear effects (Petry et al., 1998;Paoli et al., 2011).
The focus of the present study is on the impact of wake vortices on the tracer distribution rather than on the wake vortex evolution itself.Details of the vortex dynamics are only shortly addressed, since the underlying phenomena were previously discussed and quantified in the literature.
The simulations start at a wake age when it is reasonable to assume that the jets have already cooled to the ambient tem-perature, axial velocities are negligible and the vortices are fully rolled up.The simulations are carried out until the coherent vortex structures disappeared, the downward motion stopped, and the plume in the primary wake is released.
At that time, the temporal changes in the concentration distribution become smaller and smaller as the impact of the aircraft-induced dynamics on the tracer distribution has mostly ceased.Then, the final plume dimensions can be used to initialise larger-scale models in which aircraft-induced dynamics are not incorporated.
The paper is structured as follows.In the following section, we begin with a brief description of the two employed models -EULAG and NTMIX -and the simulation setup.Section 3 discusses the plume evolution in detail for a reference case.This is followed by thorough sensitivity analyses varying stratification, ambient turbulence, linear vertical wind shear, vortex initial circulation, and initial spatial distribution.Section 4 compares the EULAG simulation results with NTMIX results and aircraft in-situ observations.More general implications of our results are discussed in Sect. 5. Conclusions are drawn in the Sect.6.

Model description and setup
The simulations have mainly been carried out with the computational model EULAG (Prusa et al., 2008), proven for large-eddy simulation of turbulent thermally-stratified flows.In addition, the alternative computational model NTMIX was applied to selected experimental setups, in order to corroborate the results obtained with EULAG.

Dv
In Eqs.
(1), v, θ, ρ, and π represent, respectively, the velocity vector, potential temperature, density, and a densitynormalised pressure.Moreover, D/Dt = ∂/∂t + v • ∇, and g = (0, 0, −g) T is the gravitational acceleration vector.The subscript b appearing with ρ and θ refers to prescribed horizontally-homogeneous basic state profiles in hydrostatic balance (Clark and Farley, 1984).In addition to the basic state, a more general ambient state denoted by the subscript e is introduced, see Prusa et al. (2008) for a discussion.In our setup, the basic and ambient state are identical.The primed variables θ and π correspond to deviations from the ambient state.The last terms on the rhs of Eqs.(1a), ( 1b) and (1d) represent sub-grid (turbulent) dissipation of momentum and diffusion of heat, with the symbols τ , h and j denoting the deviatoric stress tensor, the heat flux vector and the scalar diffusion flux vector, respectively.Here, the latter are modelled by means of an eddy viscosity/diffusivity assumption with dissipation/diffusion coefficients proportional to the square root of the prognostic "turbulent kinetic energy" (TKE) (Schumann, 1991;Piotrowski et al., 2009).
The prognostic equations of the governing set (1) can be written in a generic conservation law form with ψ = ψ(t, x) denoting the transported mass-specific variable (i.e.ψ represents the three components of the velocity vector v or the potential temperature perturbation θ ) and R ψ summarises the corresponding rhs forcing terms.Following Smolarkiewicz and Margolin (1998), a compact description of the prognostic solution algorithm for Eq.(2) on the regular numerical grid (t n , x i ) can be given as where here A i symbolises the flux-form non-oscillatory forward-in-time (NFT) advection transport scheme MP-DATA (for multidimensional positive definite advection transport algorithm, Smolarkiewicz and Margolin, 1998), and ψ ≡ ψ n + 0.5 δt R ψ | n .The n, n + 1/2, and n + 1 superscripts in Eq. ( 3) denote the time level, the subscript i is the spatial grid vector index, and δt the time step increment.The scheme in Eq. ( 3) is implicit with respect to the dependent variables ψ resolved on the co-located grid (dissipation/diffusion terms are treated explicitly with first-order accuracy, see Smolarkiewicz and Margolin (1998)).Algebraic inversion and insertion into the mass continuity Eq. (1c) leads to an elliptic boundary value problem for the pressure variable π .The complete development and details about the applied preconditioned generalised conjugate residual iterative solver can be found in Prusa et al. (2008); Smolarkiewicz and Szmelter (2011) and Kühnlein et al. (2012).Finally, for second-order accuracy of the MPDATA module A in Eq. ( 3), the availability of an O(δt 2 ) estimate for the advective velocity v at the intermediate time level t n+1/2 is required, and various predictor schemes exist for that purpose (Smolarkiewicz and Margolin, 1998).For the wake vortex flows considered herein, a variant of a nonlinear predictor scheme has been found useful, which is highlighted in Appendix A. The applied nonlinear predictor scheme allows to run the model with a Courant-Friedrichs-Lewy (CFL) number up to 1 for the wake vortex flows (compared to a CFL number of ≈ 0.2 using a linear extrapolation scheme).Thus, the time step t could be substantially increased.
Recently, a microphysics module with Lagrangian ice particle tracking was coupled to EULAG and forms the model version EULAG-LCM (Sölch and Kärcher, 2010).For the present study, all microphysical processes are switched off and the simulation particles (representing the passive tracer) are simply advected with the grid point velocity plus an additional turbulent velocity.The turbulent velocity accounts for the sub-grid scale turbulent dispersion via the TKE equation.The EULAG model has been extensively used for aviationrelated numerical modelling in the past (Unterstrasser et al., 2008;Unterstrasser and Gierens, 2010b;Jeßberger et al., 2013).

NTMIX
NTMIX is a high-order solver developed for direct and large-eddy simulation of turbulent two-phase flows based on an Eulerian-Lagrangian approach.The fully compressible Navier-Stokes equations are solved in conservative form for the gaseous carrier phase while the dispersed phase is treated using a Lagrangian particle tracking method.As the details of the model and its derivation have been already presented elsewhere (Paoli et al., 2003(Paoli et al., , 2008) ) only a short summary is presented here.In the LES version of the model used in this study the transport equations for the velocity field v, total energy E ≡ c v T + (v • v)/ 2 (with c v the gas specific heat at constant volume), and scalar field C are filtered spatially so that any variable φ(x) is decomposed into a resolved part φ(x) and a sub-grid scale part φ (x), with φ(x) = φ(x)+φ (x).For compressible flows, Favre-filtered variables are used, defined as φ(x) = φ(x) + φ (x), with φ = ρφ/ρ.Using this approach, the Favre-averaged Navier-Stokes equations in Cartesian coordinates are (Paoli et al., 2003): where ∇ C is the molecular scalar diffusion flux vector, while the resolved pressure, temperature and density are related by the perfect gas law p = R ρ T (with R the specific gas constant).
The Prandtl and Schmidt numbers are set to Pr m = Sc m = 0.7.The turbulent (sub-grid scale) stress tensor, τ = −ρ v ⊗ u − v ⊗ ṽ , the turbulent heat flux vector, h = −ρc p T v − T v , and the turbulent scalar diffusion flux vector j = −ρ Cv − C v , are modelled through sub-gridscale eddy-viscosity concept (Moin et al., 1991;Lesieur and Comte, 2001): where the trace of the turbulent stress tensor tr (τ ) is negligible for weakly compressible flows as those studied here (Ng and Erlebacher, 1992;Lesieur and Comte, 2001), while turbulent Prandtl and Schmidt numbers are assumed to be constant with Pr t = Sc t = 0.419 as done in previous simulations of aircraft wakes flows (Gerz and Holzäpfel, 1999).The turbulent viscosity µ t is obtained through the Filtered Structure Function model (Métais and Lesieur, 1992;Ducros et al., 1996).Equations ( 4a)-(4d) are discretised on colocated meshes with non-uniform grid spacing using a threedimensional, finite differences Navier-Stokes solver.The spatial derivatives are computed using the sixth order compact scheme by Lele (1992) with modified coefficients to take into account the exact metrics of the mesh (Gamet et al., 1999).Compact schemes are known to provide spectral-like resolution, which is particularly helpful to represent accurately turbulent fluctuations in DNS and LES of turbulent flows (Lele, 1992).Finally, integration in time is performed by means of a low-storage third-order Runge-Kutta method with maximum CFL number fixed here to 0.5.The code has been used and validated for temporal simulations of corotating (Le Dizès and Laporte, 2002;Nybelen and Paoli, 2009) and counter-rotating (Laporte and Corjon, 2000) vortices, jet-vortex interactions (Paoli et al., 2003), and contrails (Paoli et al., 2004(Paoli et al., , 2013)).

Simulation setup
The domain has dimensions L x = 384 m, L y = 400 m and L z = 768 m.The transverse direction is along x, the flight direction (longitudinal) along y, and the vertical along z.The grid resolution is 1 m along the vertical and transverse direction and 2 m along flight direction.Periodic boundaries are applied in the horizontal and rigid boundaries in the vertical direction.The domain in the transverse direction L x is wide enough to neglect the influence of one vortex on the other across the period boundary.The size of L y is chosen to allow for the most unstable (Crow) mode of the elliptic instability of the counter-rotating vortex pair to develop in the domain (Leweke and Williamson, 1998;Laporte and Corjon, 2000).
In order to avoid interaction with the upper and lower boundaries, a relatively large value for L z is used.For the NTMIX simulations, the resolution along the flight direction is 4 m.In Sect.3.6, a EULAG sensitivity study with varied y shows that this parameter is of minor importance.Throughout the simulation, the domain size and grid resolution are kept constant.Thus, no mapping on a coarser grid or an extension to a larger domain is performed, which is advantageous for a consistent treatment of the ambient turbulence field.
The time step t depends on the setup (e.g.initial circulation 0 , vertical wind shear s) and is around 0.02-0.03s.At the later stages of a simulation, once the vortices are weaker, it is increased to 0.06 s.The simulations were run until T sim = 360 s and the number of time steps is of order 10 4 .In NTMIX, the time step t NTMIX = 0.002 s is considerably smaller for the limitation imposed by the acoustic waves and the simulations are more expensive.
A Lamb-Oseen vortex with core radius r c = 3 m and initial circulation 0 = 458 m 2 s −1 is prescribed.The cruise altitude is located 168 m below the upper boundary and identified with z = 0 m.The vortex separation distance is b 0 = 47.3 m and corresponds to an aircraft with wingspan b span = 60 m (type: Airbus A340 or Boeing B747).
The passive tracer is initialised within two discs (normal to the flight direction) of radius R init = 20 m which are centred around the vortex cores.This gives an initial plume cross section of 2.5 × 10 3 m 2 .The concentration C is uniformly C 0 = 1.0 m −3 and homogeneous along the flight direction.The tracer amount (per meter of flight path) is C tot = 2.5× 10 3 m −1 .In Sect.3.5, this rather simple initialisation is tested and compared to a more detailed initialisation from a simulation treating the interaction between the vortex and a model exhaust jet.
By default, two passive tracers are simulated in each run: the first tracer is initialised by a simple Eulerian variable and is advected using Eq.(1d), while the second is represented by means of Lagrangian particles and is advected using the LCM module.The two types of tracer do not interfere.The effect of the underlying advection scheme is analysed in Sect.3.6.Unless noted, solely tracers with Lagrangian treatment are presented in the following sections.Since the advection algorithm in NTMIX is not positive definite (only important for the Eulerian tracer), the initial concentration discontinuity at r = 20 m is smoothed with a standard tanhprofile in radial direction (e.g.Paoli et al., 2003).For EU-LAG, which uses the non-oscillatory MPDATA advection algorithm, and for the Lagrangian tracers, the initialisation as a step function is appropriate.
In the Lagrangian scheme the passive tracer is represented by N SIP = 150 simulations particles (SIP) per grid box.Overall, more than 75 million SIPs make up the plume.
For the reference case atmospheric conditions typical of the upper troposphere are chosen.The parameter values are listed in Table 1.
In order to introduce fluctuations in u, v, w, θ and p, a priori LES of decaying turbulence in a stably stratified atmosphere have been carried out as in Misaka et al. (2012) for all combinations of N BV and (background turbulence).
In addition, white noise with maximal u rms = 2 m s −1 is introduced along the core radius r c .The u rms -value decays exponentially with |r − r c | following Holzäpfel et al. (2001) (aircraft-induced turbulence).
No vertical shear of the cross-wind u is prescribed in the reference run.In a sensitivity study, a linear wind profile u(z) = s • (z − L z /2) was superimposed, where s is the wind shear.
For the analysis of wake vortex characteristics (vortex core position, circulation, topology), we use a 3-D vortex tracking algorithm developed by Hennemann and Holzäpfel (2011).During post-processing of data the algorithm first determines the three-dimensional path of the vortex core line and then computes piecewise curvature radii as a measure of vortex deformation.Hence, the coherent vortices can be tracked even after strong deformation and reconnection.This tool was found to be helpful for the interpretation of the simulation results.

Dilution of a passive tracer
In Sect.3.1, the wake vortex and tracer evolution is presented for the reference case (see Table 1 for the parameter settings).The enhanced tracer dilution due to wake dynamics is highlighted by comparing the reference case to simulations without wake vortices and aircraft-induced turbulence.
Thereafter, the individual subsections investigate the sensitivity of various parameters -thermal stratification (given by Brunt-Väisälä frequency N BV ), turbulence intensity (given by eddy dissipation rate ), vertical wind shear s, initial circulation 0 , initial spatial distribution and numerical parameters -on the extent and dilution of the plume.
In previous studies, the dispersion in the free atmosphere was computed for 2-D-Gaussian plumes (uniform along flight direction).Assuming that the plumes remain Gaussian, the evaluation of the variances of the tracer distribution (σ 2 v , σ 2 h , σ 2 s in the vertical, horizontal and shear directions, respectively) provide a notion on the half-widths of the plume.The increase in plume dimensions (determined by the σ -values) can be translated into turbulent diffusivities D v and D h of the turbulent field by analytical formulas (Konopka, 1995).As will be shown later, the plumes affected by wake dynamics become clearly non-Gaussian along the vertical direction (see e.g.Fig. 3).During the vortex phase, the dynamics are far more complex than in a simpler case with solely turbulent mixing in the free atmosphere.Thus, the plume dimensions and dilution cannot be sufficiently well described by σ and, consequently, D v and D h cannot be derived.We focus on discussing 1-D-profiles of the tracer amounts along all three dimensions: The tracer profiles are computed and Table 1.Numerical, atmospheric, aircraft and tracer parameters: Brunt-Väisälä frequency N BV , eddy dissipation rate , pressure p 0 and density of air ρ air at the lower boundary, vertical wind shear s, wingspan b span , vortex separation b 0 , initial circulation of each vortex 0 , reference timescale t 0 , initial descent speed w 0 , core radius r c , maximum rms of additional turbulent fluctuations around the core region u rms , radius of initial discs R init , tracer concentration C, number of simulation particles N SIP (in total and per grid box).

Numerical parameters
x, z i.e. the concentration distribution is integrated along the transverse and vertical direction and averaged along the longitudinal direction.For clarification, we remark that -except for Sect.4.2 -no profiles of average concentrations are discussed.
For the computation of the following quantities the Lagrangian tracer field was first smoothed with a running mean.Each SIP represents a constant concentration of C SIP = C 0 /N SIP .The smoothing explains concentrations C < C SIP .For the determination of the plume volume V , we consider all grid cells with concentration C > C low = 10 −3 m −3 : The (average) plume cross-sectional area is then given by A := V /L y .Similarly, for the computation of the median or mean concentration only grid boxes with C > C low are taken into account.The vertical plume centroid is computed as Similarly, Z c can be computed for each slice normal to the flight direction and is then denoted as Z c (y).

Reference simulation
The distribution of the (passive) tracer is controlled by the dynamics of the wake vortex evolution.This is exemplified for the reference case in Fig. 1 which displays (different) iso-surfaces of plume concentration and vorticity magnitude at three evolution stages.Since the vortex pair and the air within the oval-shaped streamline, which separates the air sinking with the pair from the environment (Lamb, 1879), descend in a stably stratified atmosphere, baroclinic forces develop between the two air masses inside and outside the oval (Holzäpfel et al., 2001).These forces create secondary vorticity structures (SVS) which wrap around the primary vortices in more or less regular longitudinal intervals as depicted in the figure after 1 min.The tracer is initially distributed in two areas around the vortex centres with radii of 20 m.That is, from the beginning a fraction of the exhaust is found close to the border of the oval where baroclinicity builds up and, thus, gets affected by the SVS which detrain some exhaust and transport it around the oval.At that time the primary vortices and also most of the tracer inside the oval are still organised in two straight tubes but the tracer concentration varies strongly along that longitudinal direction already at that early instant of time as can be imagined from Fig. 1 and will become apparent later when we discuss longitudinal variations of C l (y).
When the vortex pair continues to descend the baroclinicity further increases (it gets "accumulated" along the oval) and more exhaust material gets detrained by the SVS.These air and exhaust parcels are warmer than the local environment and rise adiabatically towards the original emission height at flight level.Hence, after t = 3 min, the detrained exhaust forms the so-called curtain between the current descent distance and the flight level.At that time the topology of the primary vortices shows the effect of the long-wave elliptic (Crow-) instability (Crow, 1970): The two vortices, which oscillate in inclined planes, are wider separated around y = 80 m and only very little separated around y = 280 m compared to the initial separation b 0 .The closer the vortices, the quicker they descend and the more rigorous the SVS produced by the accumulated baroclinicity become.Therefore, we see the thickest cloud of detrained exhaust in the segment around y = 280 m, where the primary vortices will soon link.When the two primary vortices eventually touch and connect (forming a vortex ring), the tracer fraction in that segment which is still stored inside the oval appears to "explode" and dilute rapidly.
At t = 6 min the primary vortices have already decayed leaving two "clouds" of incoherent vorticity, i.e. turbulence.Buoyancy effects lead to a re-accumulation of the tracer around flight altitude.The lower part of the plume displays the characteristic puffs due to vortex reconnection and ring formation (e.g.see photograph and LES results in Lewellen and Lewellen, 1996;Gerz et al., 1998).For a detailed analysis of the underlying mechanisms of vortex dynamics we  refer the reader to the papers of Holzäpfel et al. (2003) and Misaka et al. (2012).
Figure 2 shows tracer concentrations within the x − zplane.In order to demonstrate the longitudinal variations which arise solely due to dynamical processes, specific slices rather than a longitudinal average are displayed.Slices with minimum/maximal tracer amount C l (y) and smallest/largest vertical plume centroid Z c (y) at three selected times are illustrated in separate panels.Close to the segment of vortex reconnection (y ≈ 280 m) the slice with maximum tracer content and concurrently with minimal Z c appears at y ≈ 300 m.Especially, the area with high concentrations in the primary wake is larger than in the segments with minimal C l and maximal Z c .After the dissolution of the primary vortex pair to incoherent turbulence it is evident from the bottom row panel in Fig. 2 that on average the concentrations in the detrained and uplifted secondary wake are higher than in the plume left by the former primary wake for all four displayed slices.The slices with max C l and min Z c feature the puffs already known from Fig. 1, whereas for the two other cross-sections the former primary wake is more or less void of tracer.
Some aspects of the tracer evolution are more apparent in 1-D-profiles shown in Fig. 3. Initially, the vortex pair descends at a rate of w 0 = 0 /(2π b 0 ) = 1.55 m s −1 .Accordingly, the plume is centred at z = −90 m after 1 min.After 3 min the plume reaches a vertical extent of 300 m.Detrained and buoyant particles form the curtain between the current descent distance and the original emission height.After 6 min it even overshoots the flight level (z = 0).Moreover, the peak at flight level is larger and more material is contained in the secondary wake than in the primary wake.The final vertical plume structure may be approximated by a superposition of three Gaussian distributions (depicted with +-symbols) centred at the peaks of the primary wake, the original emission altitude and an intermediate level.The

Impact of stratification and turbulence
The decay of the wake vortex pair depends on stratification and ambient turbulence (Hennemann and Holzäpfel, 2011).
The selected values of Brunt-Väisälä frequency N BV and eddy dissipation rate (Table A1) cover the range typical of upper tropospheric conditions.The simulations reveal that for the chosen background values the variation of stratification has a much stronger effect on the plume evolution com-pared to the variation of the turbulence intensity.In particular, the vertical tracer distribution (first row of Fig. A3) is significantly affected by stratification.
In the strongly stable environment, the vortices break up faster, the tracer is at most 400 m below the emission altitude, and buoyant forces cause the tracer to rise and even protrude above the original emission level.In the case with weakly stable stratification, the vortex pair descends further down (nearly 500 m) and carries most of the tracer to much www.atmos-chem-phys.net/14/1/2014/Atmos.Chem.Phys., 14, 1-24, 2014  approximation can be computed by where N is a Gaussian distribution with mean z i and standard deviation σ i .The coefficients are given in Table 2.In Sect.3.2 we will see that the relative importance of the two peaks depends on stratification and no general approximation can be given.
Along the transverse direction, two peaks can be still identified after 1 min.Later, the two separate plumes merge into a single plume with concentration peaking in the middle (at least in an averaged sense).The profile tends to become Gaussian and a least-square fit yields σ h = 42 m for the 6 min-plume (see +-symbols).
Additionally, two NOVORTEX-simulations are performed where the velocity fields are initialised simply with the turbulent background fields (no wake vortices nor aircraft-induced turbulence are prescribed).The eddy dissipation rate is = 10 −5 m 2 s −3 .Hence, it is one hundred times larger than the reference value, which is at the lower end of typical upper  tropospheric and lower stratospheric values.Vertical shear is 0 or 0.006 s −1 .As expected the vertical extent of the NOVORTEX-plumes is much smaller since (natural) vertical mixing (especially in a stably stratified atmosphere) is lower than vortex-induced transport and mixing.Horizontal mixing is strong enough to partly merge the two initially separated jets.
In the first minutes, as long as the primary vortices are intact, they trap the exhaust and inhibit horizontal dispersion, as the NOVORTEX-plumes are broader than the reference plume at t = 3 min.The tracer concentration C l varies considerably (by a factor 2 to 4) along flight direction y.Initially, these variations of C l (and the detrainment pattern as in Fig. 1 at t = 1 min) occur randomly at small scales in y and are caused by the secondary vorticity structures.Consequently, these small-scale variations do not show up in the NOVORTEX-cases.At t = 3 and 6 min the variations are spatially smoother and are associated with the evolving longwave Crow instability and the formation of a vortex ring.This behaviour can often be observed in the sky when contrails have regular mammatus or smoke rings (see Fig. 1 at 6 min).We note that also for the NOVORTEX-cases the variations are eventually similar in amplitude.However, they ap-pear to be a bit more noisy, as turbulent structures are less coherent than the vortex-induced dynamical pattern.
Consistent with the stronger vertical transport and mixing, the reference run features a concentration PDF (Fig. 4 left) peaking at lower concentration value than the NOVORTEXsimulations.After three minutes the most frequent concentration value is about one order of magnitude smaller than the initial value C 0 = 1 m −3 .Similarly, the median of the concentration PDF drops down to 5 % of the initial value (Fig. 4 middle).Whereas the cross-sectional plume area barely doubled for the NOVORTEX-cases, enhanced mixing due to wake vortices leads to an increase by a factor of 17 after 6 min (Fig. 4 right).Interestingly, the impact of vertical wind shear on the turbulent mixing is low and no significant differences arise in the concentration PDFs, median and the plume area between the two NOVORTEX-cases within 6 min.

Impact of stratification and turbulence
The decay of the wake vortex pair depends on stratification and ambient turbulence (Hennemann and Holzäpfel, 2011).The selected values of Brunt-Väisälä frequency N BV and eddy dissipation rate (Table 1) cover the range typical of upper tropospheric conditions.The simulations reveal   that for the chosen background values the variation of stratification has a much stronger effect on the plume evolution compared to the variation of the turbulence intensity.In particular, the vertical tracer distribution (first row of Fig. 5) is significantly affected by stratification.
In the strongly stable environment, the vortices break up faster, the tracer is at most 400 m below the emission altitude, and buoyant forces cause the tracer to rise and even protrude above the original emission level.In the case with weakly stable stratification, the vortex pair descends further down (nearly 500 m) and carries most of the tracer to much lower altitudes compared to the strongly stratified case.Since buoyancy is weak, the plume top barely reaches the original emission level and the primary wake contains more material than the secondary wake.The final difference in the vertical plume centroid is more than 150 m between the two stratification scenarios (Fig. 6 right).Moreover, this figure nicely demonstrates the minor impact of the turbulence variations typical at cruise level.Here, the stronger turbulence causes a slightly earlier vortex break up.(Note that in a neutrally or convectively, i.e. unstably, stratified atmosphereas, for example, in the atmospheric boundary layer closer to the ground -vortex evolution would crucially depend on turbulence.) Along the transverse/longitudinal direction (second/third row of Fig. 5) the tracer distribution evolves similarly within the first minutes.After 3 to 4 min the plume becomes broader in the weakly stable case.This is connected to the formation of much wider vortex rings as the application of the vortex tracking algorithm by Hennemann and Holzäpfel (2011) reveals.More pronounced vortex rings also lead to a strongly inhomogeneous distribution along the flight direction with a strong tracer accumulation in segments of only 100 m length.
Although the plume shapes evolve fairly differently, no strong signature of the meteorological parameters can be found in the concentration PDFs (Fig. 7) nor in the temporal evolution of the median concentration and the plume area (Fig. 6).

Impact of the vertically sheared cross-wind
The linearly sheared cross-wind takes the values s = 0, 0.002, 0.004 or 0.006 s −1 , respectively.This resembles the typical range observed in the upper troposphere (Dürbeck and Gerz, 1996;Unterstrasser and Gierens, 2010a).
During the vortex descent the plume gets more and more tilted and consequently the width of the plume increases.The higher the wind shear, the stronger the broadening of the plume.After 6 min the plumes have widths of around 200 m, 300 m, 500 m and 600 m for the shear values listed above (not shown).
A linear profile of the background cross-wind s = dU/dz = const does not directly affect the primary vortices y , but it enhances the component ω x by tilting the component ω z of, e.g.background turbulent eddies or the evolving secondary vorticity structures which, in turn, act on the primary vortices.Since the timescale 1/s of the strongest shear rate case is about 3 min, thus 6 times larger than t 0 , the vortex pair has to descend by about 6 initial vortex spacings (280 m) before the shear has a sensible effect on the SVS and the mixing.This fact is corroborated by the evolution of the vertical tracer profiles which is shown in Fig. 8 for three different instants of times.Whereas the vertical tracer distribution is practically the same for all cases at t = 1 min and 3 min, it differs at t = 6 min.The higher the wind shear, the lower the location of the secondary wake.We speculate that shear-induced turbulence enhances mixing which dilutes the   buoyant plume more strongly and thus weakens the updraft.
The vertical plume centres diverge less than 30 m after 6 min and even the evolution of the plume areas shows only small differences, less than 10 % (not shown).

Impact of aircraft parameters
The vortex descent and decay is affected by the initial circulation 0 which is given by 0 For a fixed aircraft type (i.e.vortex separation b 0 unchanged), the initial circulation mainly depends on the aircraft mass M (payload, kerosene supply) since the density of air (ρ air = 0.4 kg m −3 ) and cruise speed (U = 0.78 Mach) do not change significantly in cruise conditions.The reference case with 0 = 458 m 2 s −1 corresponds to a mass at the lower end of the range for a A340/B747.Thus, the reference case is compared to cases with larger 0 -values of 520, 580 and 650 m 2 s −1 .The larger the initial circulation, the faster the initial descent (w 0 = 1.55, 1.75, 1.96 and 2.20 m s −1 , respectively) which is apparent in the different vertical displacements of the plumes in Fig. 9.The mutual velocity induction of the two vortices is largest for the high 0 -case and detrimental effects upon each other lead to a faster vortex decay (see e.g.Holzäpfel, 2006).The final displacement position of the primary wake is very similar (except for the 650 m 2 s −1 -case where vortices descended around 50 m further down) since the faster descent speed is compensated by a shorter vortex life time and a stronger rebound due to buoyancy.Qualitatively, the partitioning into a primary and secondary wake is similar for all four cases and the variation is much smaller than that caused by a change in stratification.
No qualitative changes occur for the distribution in longitudinal or transverse direction (not shown).
The cross-sectional area increases by a factor of 17-24 within 6 min and is larger for higher circulation (Fig. 10).In wake vortex research it is common to describe the wake vortex age in terms of a normalised time unit t 0 = (2π b 0 2 )/ 0 .This scaling is favourable when vortices of aircraft with differing b 0 and/or 0 are compared (Gerz et al., 2002).In our case, t 0 ranges from 21 s to 30 s (see Table 1).The curves lie much closer together, when the plume area evolution is plotted over t * = t/t 0 , which confirms the validity of the scaling approach.Interestingly, the trend is even reversed, i.e.A(t * ; 0 = 458 m 2 s −1 ) is slightly larger than A(t * ; 0 = 650 m 2 s −1 ).This is consistent with the N BVdependency of A(t * ) reported in Sect.3.2.In both sensitivity tests (N BV -and 0 -variation) A(t * ) is larger for higher normalised N * BV = N BV t 0 .

Variation of initial spatial distribution
Since the initial spatial distribution of the tracer is prescribed in a rather simple way in this study, we perform a sensitivity study varying the radius R init of the discs containing the initial tracer.The reference case with 20 m is supplemented by simulations changing R init by ±5 m.Two further simulations use a small disc with R init = 5 m and an annulus with inner radius R i = 10 m and outer radius R e = 15 m, respectively.Besides these idealised initialisations, we further employ a non-uniform concentration field.This initial tracer distribution is obtained from a 3-D large-eddy simulation of the nearfield wake that takes into account the interaction between the wake vortex and the exhaust jet (Paoli et al., 2013).Fig. 11 shows a 2-D contour plot of the initial scalar distribution generated by a large 4-engine aircraft (consistent with our choice of aircraft type).The main features of the interaction are the entrainment of the jet around the vortex core and the formation of vortex rings that are later dissipated by turbulent mixing and by the relaminarazing effect of vortex rotation (as shown for example by the snapshots of vorticity in Fig. 10 of Paoli et al., 2003).These mechanisms determine the vertical and horizontal extension of the plume and hence the wake structure at the beginning of the vortex phase.The initial vertical profiles of all initialisation variants are shown in the left panel of Fig. 12.In the non-reference simulations, the concentrations C were scaled such that the total tracer amount is equal to the reference value C tot = 2 π 20 2 m −1 .After two minutes an initialisation effect is still apparent.The stronger the tracer was initially concentrated around the vortex centre, the more material is contained in the primary wake.However, the differences reduce over time.The detrainment rate and the extent of the curtain are similar for all cases except for the very narrow initialisation (R init = 5 m).At the end of the vortex phase, the vertical distribution is qualitatively very similar.Again, it can be concluded that stratification has a much stronger impact on the relative importance of primary and secondary wake.Our interpretation is that turbulence tends to homogenise the tracer distribution    in such a way that the initial structure of the plume is 'forgotten' at the end of the vortex phase.This also reflects the linear character of the scalar transport Eq. ( 4d).(Note the extrapolation to non-passive scalars is not trivial as it can be biased by nonlinear effects such as chemical reactions or microphysical processes in the case of contrails).An exception to the above argument is represented by the case of initial radius R init = 5 m where, unlike the other cases, the tracer is initially concentrated inside the vortex core.This hampers the effect of turbulence diffusion with the consequence that a large fraction of the tracer remains trapped in the primary wake.
All idealistic initialisations assumed an initially uniform concentration across the specified cross-section.Only the jet phase model output has a non-uniform distribution whose PDF is depicted in Fig. 13.The PDFs after 2 and 4 min show similar patterns for both initialisations.This suggests that also the assumption of homogeneity is uncritical for our purposes.
We conclude that this sensitivity test, a posteriori, justifies the application of the rather simple default initial distribution used throughout this study.

Numerical issues
For the sake of completeness, this section reports several sensitivity studies testing miscellaneous numerical parameters and approaches.First, we present a (limited) exercise of model verification by comparing the Eulerian and Lagrangian approaches to represent scalar transport.Note this is an important issue in the numerical modelling with a mixed Eulerian/Lagrangian approach.Active simulation particles represent, e.g.contrail ice crystals, whereas water vapour is treated as a Eulerian variable.Since these quantities are coupled via microphysical processes, similar dispersion properties are desirable to exclude spurious interactions.It is known that for passive particles the quality of Lagrangian tracking depends essentially on the number of particles (samples) used to reconstruct Eulerian statistics.On the other hand in an Eulerian approach the transport properties rely mainly on the quality of the underlying mesh and the numerical scheme -ideally low dispersive and low dissipative.For the particular case analysed here the mesh resolution and the number of  2013).For the present analysis, the concentration fields of the non-reference simulations were scaled such that the total tracer amount equals the value of the reference run.This implies that the area enclosed by the various curves is equal for all simulations.Note the different axis ranges in both directions.2013).For the present analysis, the concentration fields of the non-reference simulations were scaled such that the total tracer amount equals the value of the reference run.This implies that the area enclosed by the various curves is equal for all simulations.Note the different axis ranges in both directions.numerical particles is sufficiently high that the results can be considered as "converged" with respect to these parameters.Indeed, Fig. 14 shows that the methods yields very similar results in terms of vertical profiles and PDF distributions of the tracer.As a default the simulations started with 150 SIPs per gridbox.Additional runs with 50 and 300 SIPs per grid box gave practically identical results and confirm the convergence of the Lagrangian approach.
Along flight direction the grid resolution is y = 2 m (for all EULAG simulations).NTMIX simulations presented in the subsequent section use y = 4 m for matters of cost effectiveness.A EULAG sensitivity study with y = 1 m and 4 m shows a minor sensitivity to y and indicates the NT-MIX resolution to be reasonable (also considering the higher precision of the numerical scheme).
Furthermore, we found that using a larger initial vortex core radius r c -though affecting the circulation evolution -  was of no relevance to the detrainment and distribution characteristics of the passive tracer.

Validation and comparison with observations
To increase confidence in the simulated detrainment and dispersion of the passive tracer, comparisons with an alternative LES model NTMIX and with aircraft in-situ measurements are presented in following subsections.

Model comparison EULAG -NTMIX
The comparison of EULAG and NTMIX is performed for three different environmental conditions, namely the reference case (N BV = 1.15 × 10 −2 s −1 , = 10 −7 m 2 s −3 ) and sensitivity runs with either reduced stability (N BV = 0.5 × 10 −2 s −1 ) or increased turbulence intensity ( = 10 −5 m 2 s −3 ).The corresponding EULAG runs were already discussed in detail in Sect.3.2, and showed a large sensitivity to thermal stratification and a weaker sensitivity to the turbulence intensity.The numerical characteristics of the two models are described in Sect. 2. Basically, they differ in terms of model equations (anelastic approximation to fully compressible), resulting numerical solution techniques and the consistency order of the employed difference schemes.A conclusive comparison is achieved by using identical setups in both models (vertical atmospheric profiles, spatial tracer initialisation, turbulence fields, vortex characteristics).Generally, we find an excellent qualitative agreement between the two models.We thus confine ourselves to simply presenting vertical profiles of the tracer distribution after 5 min in Fig. 15.For each of the three environmental cases, the final vertical displacement and the plume vertical extent are similar within a few tens of metres.Regarding the relative importance of the primary/secondary wake and the buoyancy-driven overshooting, the sensitivity to stratification and turbulence is similar in both models.This very convincing model comparison gives confidence in the general model capabilities and more specifically in the presented EULAG sensitivity studies (not covered by NTMIX simulations).

Comparison with in-situ NO-measurements
Here, we compare EULAG simulations to in-situ measurements of exhaust trace gas distributions during the vortex phase.We choose measurements of reactive nitrogen species (NO y ) behind an A-319 on 19 November 2008 obtained during the CONCERT campaign (Jurkat et al., 2011) onboard the DLR research aircraft Falcon.Primarily, the campaign aimed at probing contrails.In aircraft plumes of ages smaller than 3 min NO y is dominated by NO and NO 2 , and to a much lesser amount gas phase HONO and HNO 3 (Tremmel et al., 1998) from which only a small portion may be lost to the particulate phase at this early stage (Schäuble et al., 2009).For our purpose it may therefore be treated as passive tracer to a good approximation.The atmospheric condition and aircraft properties were recorded and are employed for the simulation with EULAG shown.They are listed in Table 1/"A319 case study" along with the vortex and plume parameters.The initial NO y mixing ratio enhancement of X NO y ,0 = 46.55 nmol mol −1 inside the two circular jet plumes after 20 s is calculated according to where m NO y is the molar mass of NO 2 also used in the calculated emission index EI NO x using the DLR fuel flow method described by Döpelheuer and Lecht (1999) with a fuel flow rate η = 0.9 t/h/engine at an airspeed of V = 224 m s −1 .The Falcon follows the A319 with slightly slower airspeed to probe several plume segments of steadily increasing age while trying to stay inside the plume at different vertical positions.In Fig. 16 the 1 s-average values of the normalised NO y in-situ measurements are indicated as diamonds for selected plume age intervals of 10 seconds centred around the indicated times.We mimic this sampling strategy in the simulation by defining a straight line segment in flight direction.This segment starts in a region with enhanced mixing ratio and we build an average value for mixing ratio along our domain size in flight direction, corresponding roughly to the 1 s averages in the measurements.We repeat this procedure for a large number of different flight segments to account for the spatial heterogeneity in the tracer distribution.Thus, we reproduce the pilots intent to stay inside the plume and the horizontal heterogeneity integrated in the 1 s-average values of a single measurement point.The comparison is limited by the fact that we only have one instance of a plume segment in the simulation in contrast to a variety of different plume segments for the same plume age interval in the measurement data.Nevertheless, the simulated and measured normalised mixing ratios agree favourably concerning magnitude and position.We especially want to highlight the evolution of the detrainment of the passive tracer.At 105 s the tracer is concentrated in a confined region of 50 m vertical extent with a pronounced peak in the primary wake and lower mixing ratios in the secondary wake slightly above, visible both in  measurements and simulated data.After 135 s the vertical extent of the plume increases with a growing amount of trace species present in the secondary plume.In the last panel after 185 s, measurements were only available in the upper part of the plume.These confirm the dilution of the mixing ratio and the detrainment rate in the simulation, indicating that a separation in primary and secondary wake has almost vanished and the trace species rebounds up to the initial flight level at z = 0.

Discussion
In this study the plume evolution has been analysed for a large variety of conditions, elucidating the sensitivity to atmospheric, aircraft and numerical modelling parameters.The turn-over (mixing) timescales for the simulated background turbulence range in the order of 30 min and 6 min for the reference and higher turbulent cases, respectively.The buoyancy timescales 2π/N have orders of magnitude of 10 min and 20 min for the reference and weaker stratified cases, respectively.The timescale of the background shear varies from infinity (REF) to 3 min (S6).In contrast, the timescale for the primary vortex to rotate air parcels around its axis, 2π/ω y , is in the order of a few seconds initially and still below one minute after break-up of the coherent vortex structures.Finally, the characteristic vortex-pair descent time t 0 is about 30 s. Hence and as expected, the dynamics of the trailing wake-vortex pair controls the tracer distribution to first order, followed by the effects of stratification, shear, and background turbulence.
Figure 17 now summarises the most important sensitivities for the plume characteristics in a concise form to facilitate discussion.It is apparent that -after the vertical transport of exhaust by the trailing vortex pair -stratification has the largest impact on the vertical distribution.Especially, the plume travels further down ("vertical displacement"), if the atmosphere is less stable.Also, the plume is lower in an averaged sense with fewer tracer in the secondary wake ("vertical centre").
Linear wind shear does not influence the unperturbed primary vortices but has an effect on the secondary vorticity structure and the detrained and updrafting air mass.However, for the simulated shear rates this effect is comparatively small within the first five to six minutes.The vertical extent seems to be generally a bit smaller compared to the reference case with no shear.The plume area increases with increasing wind shear.This is certainly a well-known effect from dispersion in the free atmosphere.However, here the situation is more intricate.If vertical wind shear had led to a substantial reduction of the vertical extent, the shear-induced spreading would have been less effective on the long run and later plume dimensions could have turned out smaller for larger wind shear.
The initial circulation 0 controls (linearly) the initial descent speed of the vortex system.Expectedly, the vertical displacement and extent are larger for stronger vortices.We find a less than linear relationship between w 0 and the plume quantities, as the effect of descent speed is partly balanced by a reduced vortex lifetime for larger 0 .Whereas the position of the primary wake is most strongly affected by a variation of 0 (compared to the other studied parameters), the tracer accumulation around cruise altitude is similar for all 0 unlike to a variation of stratification (cf.Figs. 5 and 9).
The initial spatial initialisation does not affect the vertical displacement and extent.Except for the very narrow case Fig. 117.Sensitivity analyses of various simulations (grouped in blocks "Stratification & turbulence", "Shear", "Vortex", "Init" and "Num").Depicted are plume characteristics after 5 minutes (from top to bottom): cross-sectional area, vertical centre relative to cruise altitude, the maximum vertical displacement below cruise altitude and the vertical extent.The meaning of the labels on the bottom xaxis for each specific run can be found in Table 11.Within block "NUM" the simulation "EUL" shows the Eulerian tracer of the reference simulation.The triangles show NTMIX results and the diamond symbol shows the Eulerian tracer of the respective 'R05' sensitivity run.For better orientation, the dashed horizontal lines show the value of the reference run.Fig. 17.Sensitivity analyses of various simulations (grouped in blocks "Stratification & turbulence", "Shear", "Vortex", "Init" and "Num").Depicted are plume characteristics after 5 min (from top to bottom): cross-sectional area, vertical centre relative to cruise altitude, the maximum vertical displacement below cruise altitude and the vertical extent.The meaning of the labels on the bottom x-axis for each specific run can be found in Table 1.Within block "NUM" the simulation "EUL" shows the Eulerian tracer of the reference simulation.The triangles show NTMIX results and the diamond symbol shows the Eulerian tracer of the respective "R05" sensitivity run.For better orientation, the dashed horizontal lines show the value of the reference run.
("R05") the detrainment rates are astonishingly similar, resulting in only slight differences in plume centre and area.Our results are consistent with earlier studies (Huebsch and Lewellen, 2006;Misaka et al., 2012) who also found a moderate and less than expected impact of the initial tracer placement on detrainment.
For the reference case, the Lagrangian and Eulerian tracer advection methods give similar results (run "EUL" in block "NUM").In the case of a very narrow spatial distribution (compare black with green symbols of run "R05"), the tracer distribution of the two numerical approaches differ.For loworder numerical schemes as those used in EULAG, the Eulerian approach suffers from numerical diffusion (at least for the chosen mesh resolution) and overestimates the detrainment from the primary wake and consequently the plume area.Higher order schemes as those used in NTMIX are designed to capture turbulent fluctuations that develop in the flow given their (almost) spectral resolution properties.However, if the initial conditions are not sufficiently smooth, artificial diffusion has to be added to stabilise the solution, which again may result in excessive diffusion of the scalar field at later times.The Lagrangian approach seems then to offer more versatility in the choice of the flow parameters.Note that the narrow case seems academic for the present wing engine geometrical configuration, as jet phase simulations of Paoli et al. (2013) show a broader distribution at the beginning of the vortex phase.For future aircraft designs the relative position between the vortex centre and the engine exit may change.If they are closer together making a narrower initial distribution more likely, numerical dilution studies should rely on Lagrangian approaches.
In box model studies, the plume expansion is usually analytically prescribed by some dilution/entrainment rate term.From the present study dilution rates during the vortex phase can be derived.We found roughly a linearly increasing plume area A(t) which implies a t −1 -dependence of the dilution rate ω(t)(:= Ȧ(t) / A(t)).We calculate the overall dilution during the vortex phase by := A init / A final .We found to vary most strongly within the "spatial init" sensitivity study.There, A init differs by as much as a factor of 6 (= (25 m / 10 m) 2 ) and A final turns out to be similar for all realistic spatial initialisations.In all other sensitivity studies, A init is constant by definition and A final varies at most by a factor of 1.5 implying smaller variations of .We recommend to prescribe a plume area of 4.0× 10 4 m 2 to 6× 10 4 m 2 after 6 min for the given aircraft type and compute the dilution rate after choosing the initial concentration/plume area within the specific application of a box model.
The following paragraphs discuss implications of our results for the contrail microphysical evolution.Considering persistent contrails in ice supersaturated ambient air, adiabatic warming in the downward moving vortex system can lead to subsaturation and a substantial loss of trapped ice crystals (Lewellen and Lewellen, 2001a;Unterstrasser and Sölch, 2010).Besides atmospheric parameters (temperature, relative humidity) the final vertical displacement of the primary wake (corresponding to a certain temperature increase) controls the extent of crystal loss.Contrary to the ice crystals in the primary wake, detrained ice crystals do not face such a strong adiabatic heating.These ice crystals are not prone to sublimation and their microphysical evolution is distinct from the "primary" ice crystals.Thus, the rate of detrainment and the relative tracer fractions in the primary and secondary wake affect the contrail microphysical evolution.
We found that the secondary vorticity structures, which evolve along the vortex-pair oval, and the primary vortex reconnection and ring formation control the heterogenization of the tracer along flight direction on small and large scales, respectively.In favourable atmospheric conditions with persistent contrails these initial "imbalances" might trigger turbulence due to differential latent or radiative heating and enhance the contrail spreading.This may affect the long-term evolution of contrails and its transition to cirrus.

Conclusions
A large eddy simulation (LES) model with Lagrangian particle tracking has been employed to study the dilution of aircraft emissions during the first minutes behind an aircraft.Aircraft emissions were treated as passive tracers in order to focus the analysis on dynamical effects.Chemical and microphysical processes were deliberately neglected in this study.
We confirmed and differentiated earlier findings that the exhaust plume evolution is strongly affected by the wake vortex dynamics.The descent of the vortex pair causes a substantial vertical plume spreading.This is an important process to be considered in larger-scale plume dispersion models, as the vertical extent of the plume strongly affects the later horizontal spreading by vertical wind shear.
For an A340/B747 aircraft at cruise altitude, we found the vertical plume extent to range from 400 m to 550 m and the plume cross-sectional area from 4.0 × 10 4 m 2 to 6 × 10 4 m 2 after 6 min (this is when aircraft-induced dynamical effects on plume spreading became negligible).From the provided final plume areas, dilution rates for box model applications can be estimated.
We studied the sensitivity of plume dilution to stratification, ambient turbulence, vertical wind shear, initial vortex circulation and initial spatial tracer distribution.We found the vertical displacement of the primary wake to depend most strongly on stratification and the initial circulation.Due to buoyancy, parts of the tracer are lifted back towards the initial emission altitude after vortex break-up.In the final stages of the vortex phase, two peaks (the primary wake and around cruise altitude) are apparent in the vertical tracer distribution.The relative fraction of tracer in the primary wake and around cruise altitude depends dominantly on ambient stratification.In a strongly stable atmosphere more tracer accumulates around cruise altitude than in the primary wake.This phenomenon could also be found in the presented in-situ measured plumes.The final vertical tracer distribution can be best approximated by a superposition of three Gaussian distributions (covering tracer around cruise altitude, the primary wake, and the curtain between).Integrating the tracer along flight and vertical direction, the distribution along transverse direction can be well approximated by a single Gaussian.
The model was initialised with a uniform concentration of a passive tracer.Although only dilution processes are regarded in this study, we end up with a fairly broad spectrum of concentrations.The concentration values range over one order of magnitude and these variations may be primarily attributed to small-scale turbulence effects.On the other hand, variations also occur on a larger scale.Due to vortex-induced instability phenomena, heterogenization of tracer amounts along flight direction is observed.
Moreover, we employed an alternative LES model for model verification.We found a high agreement between the two models EULAG and NTMIX which provides confidence in the presented results.

Fig. 11 .
Fig. 11.Exhaust plume and vorticity distributions from the reference run after 1, 3 and 6 minutes.Shown are iso-surfaces of concentration C and vorticity magnitude |ω| (in units s −1 ) with values given in the labels.

Fig. 1 .
Fig. 1.Exhaust plume and vorticity distributions from the reference run after 1, 3 and 6 min.Shown are iso-surfaces of concentration C and vorticity magnitude |ω| (in units s −1 ) with values given in the labels.
Fig. 12. 2D plume cross-sections from the reference run after 1, 3 and 6 minutes (from top to bottom).The four columns shows different slices at y = yp along flight direction, namely the slices with minimal/maximal tracer amount Cl(y) and lowest and highest vertical plume

Fig. 2 .
Fig. 2. 2-D plume cross-sections from the reference run after 1, 3 and 6 min (from top to bottom).The four columns shows different slices at y = y p along flight direction, namely the slices with minimal/maximal tracer amount C l (y) and lowest and highest vertical plume centre Z c (y).The label in each panel denotes y p , the normalised tracer amount C l /C tot and Z c .The fields are smoothed with a running mean using two adjacent points in each direction.The circles indicate the initial plume position.

Fig. 13 .
Fig. 13.Vertical, transverse and longitudinal profiles of accumulated tracer concentrations (from top to bottom).Plume ages are 1, 3 and 6 minutes as indicated on top.The results are shown for the reference simulation (black) and two NOVORTEX-runs with turbulent mixing only (no vortex dynamics).The NOVORTEX-runs use an eddy dissipation rate of = 10 −5 m 2 s −3 .The vertical wind shear is zero (red) or s = 0.006 s −1 (blue).The black +-symbols show Gaussian fits for the reference run.The fit parameters are given in Table12.The initial plume is centred around cruise altitude at z = 0 m.Note the different Cv-ranges on the abscissa in the top row.

Fig. 3 .
Fig. 3. Vertical, transverse and longitudinal profiles of accumulated tracer concentrations (from top to bottom).Plume ages are 1, 3 and 6 min as indicated on top.The results are shown for the reference simulation (black) and two NOVORTEX-runs with turbulent mixing only (no vortex dynamics).The NOVORTEX-runs use an eddy dissipation rate of = 10 −5 m 2 s −3 .The vertical wind shear is zero (red) or s = 0.006 s −1 (blue).The black +-symbols show Gaussian fits for the reference run.The fit parameters are given in TableA2.The initial plume is centred around cruise altitude at z = 0 m.Note the different C v -ranges on the abscissa in the top row.

Fig. 3 .
Fig. 3. Vertical, transverse and longitudinal profiles of accumulated tracer concentrations (from top to bottom).Plume ages are 1, 3 and 6 min as indicated on top.The results are shown for the reference simulation (black) and two NOVORTEX-runs with turbulent mixing only (no vortex dynamics).The NOVORTEX-runs use an eddy dissipation rate of = 10 −5 m 2 s −3 .The vertical wind shear is zero (red) or s = 0.006 s −1 (blue).The black +-symbols show Gaussian fits for the reference run.The fit parameters are given in Table 2.The initial plume is centred around cruise altitude at z = 0 m.Note the different C v -ranges on the abscissa in the top row.

Fig. 4 .
Fig. 4. Left: normalised frequencies of concentration values C. For the computation of the PDFs exponentially increasing bin sizes were used.Middle and right: Median of plume concentration and plume cross-sectional area (perpendicular to the flight direction) as a function of time.Colour coding as in Fig. 3.

Fig. 16 .
Fig. 16.Temporal evolution of the median of plume concentration C (left), the plume cross-sectional area (middle) and the centre of the plume in the vertical direction (right) for different stratification and turbulence intensities.Colour coding as in Fig. 15.

Fig. 6 .Fig. 17 .
Fig. 6.Temporal evolution of the median of plume concentration C (left), the plume cross-sectional area (middle) and the centre of the plume in the vertical direction (right) for different stratification and turbulence intensities.Colour coding as in Fig. 5. S. Unterstrasser et al.: Dimension of aircraft plumes: effect of wake vortices 21

Fig. 7 .
Fig. 7. Normalised frequencies of concentration values C at t = 1 min (solid) and t = 5 min (dotted) for different stratification and turbulence intensities.For the computation of the PDFs exponentially increasing bin sizes were used.Additionally the initial pdf (t = 0 min) is shown in grey (delta function for all simulations).Colour coding as in Fig. 5.
Fig. 110.Plume cross-sectional area for four different initial circulations.Colour coding as in Fig. 19.

Fig. 111 .
Fig. 111.2D contour plot of initial concentration as obtained from the NTMIX jetphase simulation.For better comparison, the radii of the various idealised spatial initialisation are added (R init = 5, 15, 20 and 25 m).The idealised initialisations use a uniform concentration C 0 = 1.0 m −3 within two discs of radius R init .

Fig. 11 .
Fig. 11.2-D contour plot of initial concentration as obtained from the NTMIX jetphase simulation.For better comparison, the radii of the various idealised spatial initialisation are added (R init = 5, 15, 20 and 25 m).The idealised initialisations use a uniform concentration C 0 = 1.0 m −3 within two discs of radius R init .

Fig. 112 .
Fig. 112.Vertical profiles of accumulated tracer concentrations at t = 0, 1, 2 and 4 min (from left to right) for various initial spatial distributions: Uniform concentration inside discs of radius Rinit = 20 m (black, reference run), 15 m (blue), 25 m (red) and 5 m (magenta) as well as inside an annulus with inner radius Ri = 10 m and outer radius Re = 15 m (green).A further simulation was initialised with a non-uniform concentration field (brown) provided by a jet phase simulation of Paoli et al. (2013).For the present analysis, the concentration fields of the non-reference simulations were scaled such that the total tracer amount equals the value of the reference run.This implies that the area enclosed by the various curves is equal for all simulations.Note the different axis ranges in both directions.

Fig. 12 .
Fig. 12. Vertical profiles of accumulated tracer concentrations at t = 0, 1, 2 and 4 min (from left to right) for various initial spatial distributions: Uniform concentration inside discs of radius R init = 20 m (black, reference run), 15 m (blue), 25 m (red) and 5 m (magenta) as well as inside an annulus with inner radius R i = 10 m and outer radius R e = 15 m (green).A further simulation was initialised with a non-uniform concentration field (brown) provided by a jet phase simulation of Paoli et al. (2013).For the present analysis, the concentration fields of the non-reference simulations were scaled such that the total tracer amount equals the value of the reference run.This implies that the area enclosed by the various curves is equal for all simulations.Note the different axis ranges in both directions.

SFig. 113 .
Fig. 113.Normalised frequencies of concentration values C for t = 0 (solid), 2 (dotted) and 4 min (dashed) for two initial spatial distributions: Default uniform concentration within discs of radius Rinit = 20 m (black, reference run) and non-uniform concentration field (brown) provided by a jet phase simulation (Paoli et al., 2013) (brown).The total tracer amount is the same in both simulations.

Fig. 13 .
Fig. 13.Normalised frequencies of concentration values C for t = 0 (solid), 2 (dotted) and 4 min (dashed) for two initial spatial distributions: Default uniform concentration within discs of radius R init = 20 m (black, reference run) and non-uniform concentration field (brown) provided by a jet phase simulation (Paoli et al., 2013) (brown).The total tracer amount is the same in both simulations.

Fig. 116 .
Fig. 116.Vertical distribution of in-situ measured (diamonds) and simulated (plus signs) NOy volume mixing ratios for plume ages of 105 s (top), 135 s (middle), and 185 s (bottom) normalised with the initial mixing ratio XNO y ,0.The sampling strategies are discussed in the text.Note that aircraft measurements are not always available throughout the whole plume extent at every time interval.

Fig. 16 .
Fig. 16.Vertical distribution of in-situ measured (diamonds) and simulated (plus signs) NO y volume mixing ratios for plume ages of 105 s (top), 135 s (middle), and 185 s (bottom) normalised with the initial mixing ratio X NO y ,0 .The sampling strategies are discussed in the text.Note that aircraft measurements are not always available throughout the whole plume extent at every time interval.

Table 2 .
Gaussian fit parameters of the tracer distribution in transverse and vertical direction for the reference run.Along the vertical direction the plume is a superposition of three Gaussian distributions.