Icarus 281 (2017) 55–72

Contents lists available at ScienceDirect

Icarusjournal homepage: www.elsevier.com/locate/icarus

Thermal structure of the upper atmosphere of Venus simulated by aground-to-thermosphere GCM
G. Gilli a,∗, S. Lebonnois a, F. González-Galindo b, M.A. López-Valverde b, A. Stolzenbach c,F. Lefèvre c, J.Y. Chaufray c, F. Lott aabcLaboratoire de Météorologie Dynamique (LMD)/IPSL, CNES, CNRS, Paris, FranceInstituto de Astrofísica de Andalucía (IAA), CSIC, Granada, Spain
Laboratoire Atmospheres, Milieux, Observations Spatiales (LATMOS), IPSL, Paris, Francea r t i c l e

i n f oArticle history:
Available online 15 September 2016a b s t r a c t
We present here the thermal structure of the upper atmosphere of Venus predicted by a full selfconsistent Venus General Circulation Model (VGCM) developed at Laboratoire de Météorologie Dynamique(LMD) and extended up to the thermosphere of the planet. Physical and photochemical processes relevantat those altitudes, plus a non-orographic GW parameterisation, have been added. All those improvements
make the LMD-VGCM the only existing ground-to-thermosphere 3D model for Venus: a unique tool toinvestigate the atmosphere of Venus and to support the exploration of the planet by remote sounding.
The aim of this paper is to present the model reference results, to describe the role of radiative, photochemical and dynamical effects in the observed thermal structure in the upper mesosphere/lower thermosphere of the planet. The predicted thermal structure shows a succession of warm and cold layers,
as recently observed. A cooling trend with increasing latitudes is found during daytime at all altitudes,
while at nighttime the trend is inverse above about 110 km, with an atmosphere up to 15 K warmer
towards the pole. The latitudinal variation is even smaller at the terminator, in agreement with observations. Below about 110 km, a nighttime warm layer whose intensity decreases with increasing latitudes
is predicted by our GCM. A comparison of model results with a selection of recent measurements shows
an overall good agreement in terms of trends and order of magnitude. Signiﬁcant data-model discrepancies may be also discerned. Among them, thermospheric temperatures are about 40–50 K colder and
up to 30 K warmer than measured at terminator and at nighttime, respectively. The altitude layer of thepredicted mesospheric local maximum (between 100 and 120 km) is also higher than observed. Possible interpretations are discussed and several sensitivity tests performed to understand the data-model
discrepancies and to propose future model improvements.
© 2016 Elsevier Inc. All rights reserved.

1. Introduction
The recent exploration by the Venus Express (VEx) orbiter andthe increasing number of ground based observations noticeably
improved the knowledge of our neighbouring planet. Prior to VEx,
satellite and ground-based observations of the upper mesosphereand thermosphere of the planet (90–180 km) were both limited in spatial and temporal coverage. Measurements by Pioneer
Venus (PV), both mass spectrometer and drag observations down
to 140 km (Niemann et al., 1980), together with accelerometer data below 130 km (Keating et al., 1985; Seiff, 1983), provided data on the composition and structure of the Venus ther-∗

Corresponding author.
E-mail address: ggilli@lmd.jussieu.fr (G. Gilli).

http://dx.doi.org/10.1016/j.icarus.2016.09.016
0019-1035/© 2016 Elsevier Inc. All rights reserved.

mosphere for the ﬁrst time. This enabled the construction of a 1-D
semi-empirical model of temperature and composition, the VTS3
by Hedin et al. (1983), mainly based on PV measurements above
140 km and on model extrapolations assuming hydrostatic equilibrium between 100–140 km.
After eight years in orbit (2006–2014), exceeding its planned
life, VEx collected a considerable amount of data from several instruments on board (Svedhem et al., 2015). Also, aerobraking manoeuvres carried the spacecraft down to 130 km altitude at theend of the mission (Svedhem and Müller-Wodarg, 2014), provided
direct measurements of local atmospheric densities and temperature of a region diﬃcult to sound (130–160 km). The short andlong-term density variability observed at those altitudes is surprisingly higher than expected and represents a challenge to current
theoretical models. All those new measurements revealed a dy-

56

G. Gilli et al. / Icarus 281 (2017) 55–72

namically active upper atmosphere, much more complex than predicted by VTS3. The interaction between lower and upper regions
and gravity wave propagation (Altieri et al., 2014; Zalucha et al.,2013) have to be taken into account to interpret the observed results.
Temperature values have been successfully retrieved by several
instruments on board VEx. The radio occultation measurement by
VeRa (Venus Radio Science experiment) provided high vertical resolution pressure and temperature proﬁles in the altitude range 40–
90 km (Tellmann et al., 2009). Strong IR absorption bands of CO2
at 4.3 μ m observed by VIRTIS (Visible and InfraRed Thermal Imaging Spectrometer) between 65 and 96 km altitudes have beenused to retrieve nighttime temperature and aerosol altitude proﬁles in the mesosphere (Grassi et al., 2010; 2014; Haus et al., 2013,2014; Migliorini et al., 2012). In addition, ground based observations provide higher spectral resolution and a better spatial coverage which complement satellite observations: sub-millimeter observations of CO lines from Earth constrain on the thermal structure between approximately 70 and 110 km (Clancy et al., 2012).
Kinetic temperature can also be calculated at about 110 km fromthe width of narrow CO2 emission lines using infrared heterodyne spectroscopy (Sonnabend et al., 2010). Those results showedthat the thermal structure of the Venus mesosphere varies strongly
with latitude: polar regions were surprisingly observed up to 20 K
warmer than equatorial latitudes at most pressure levels between
the cloud tops and 100 km. They also conﬁrmed the presence of acolder region (the so called “cold collar”) around 65°–70°S latitude,
at about 102 Pa (65 km altitude approximately).
The upper mesosphere/lower thermosphere was sounded by solar/stellar occultation experiments at the terminator and nighttime
by SOIR/VEx (Mahieux et al., 2015) and SPICAV/VEx (Piccialli et al.,2015), respectively. They show a complex thermal structure above
90 km with a large vertical and spatial variability, still not
fully explained by current models (Bougher et al., 2015). At theterminator, the observed temperature proﬁles are characterisedby a local maximum around 105 km, that was interpreted by
Bougher et al. (2015) to be produced by solar absorption in the
CO2 IR bands during daytime and then advected to the terminator. Above, a cold layer around 125–130 km is observed, which
Bougher et al. (2015) suggested to be caused by CO2 15 μm cooling.
Above 130 km the observed thermal structure shows features of the extreme EUV absorption, with increasing temperatures
above 140–150 km. At night, a warm region possibly resulting
from subsidence of day-to-night air circulation was observed by
SPICAV in the altitude range 90–100 km. On the dayside, temperatures from 100 to 150 km have been retrieved using infrared CO
non-LTE emission at 4.7 μm by VIRTIS-H (Gilli et al., 2015). A hintof local maximum (230 ± 40 K) around 115 km near the terminator at equatorial latitudes is observed, although with large uncertainties. However, at noon in the equatorial regions, those retrieval
results show lower temperatures than expected by previous models (Brecht et al., 2012) and in contrast to an upper mesospherein pure radiative balance. In addition to this, nadir measurementsof non-LTE CO2 emissions at 4.3 μm by VIRTIS-H on board VEx
have been used to retrieve daytime temperature in a wide altitude range, between 10 and 10−3 Pa (100–140 km, approximately)
by Peralta et al. (2016). They show a maximum of 195 ± 10 K
in the sub-solar region decreasing with latitude and local time, inagreement with Brecht et al. (2012). These nadir non-LTE measurement also show a “plateau” between 30°S and 30°N, not predictedby global circulation models. Daytime thermosphere temperatures
have been also inferred from a few CO2 UV dayglow observations
by SPICAV (Chaufray et al., 2012). These authors estimated the temperature near 145 km at 290 ± 60 K, comparable with noon temperatures predicted by the VTS3 model by Hedin et al. (1983).

In the last decade the scientiﬁc community has conducted
several studies to carry out a systematic and exhaustive intercomparison and validation of both satellite and ground-based telescope measurements after the arrival of VEx, taking also into account the geometry of observations and the different retrieval
techniques. Those efforts have been accomplished within the context of several international projects (e.g the scientiﬁc Team supported by the International Space Science Institute (ISSI), Bern,
Switzerland (Limaye et al., 2015)) and they contribute an unprecedented insight into the complex dynamical mechanisms ocurringin the atmosphere of our sister planet. Nevertheless, VEx results
also raised new challenging questions to be answered. Sophisticated climate models and synergies with the data analysis are fundamental to interpret the observed results and to help build a consistent picture of the spatial and temporal evolution of Venus atmosphere. Unfortunately, too many free parameters and theoreticaluncertainties (i.e limited gravity wave constraints for wave breaking, quenching rate coeﬃcient for the IR mesospheric cooling, etc)
prevent existing GCMs to use a unique set of wind ﬁelds, eddy
diffusion coeﬃcients and wave drag parameters to reproduce observed density, temperature and airglow variations (Bougher et al.,2006).
A 3D Thermospheric General Circulation Model (TGCM) for
Venus with coupled energetics, dynamics and composition described in Bougher et al. (1988) was recently updated to address
the key VEx observations in the altitude range between 70 and
200 km, in particular to reproduce the observations of the various airglow distributions (Brecht et al., 2011), dayside temperatures (Brecht and Bougher, 2012) and temperatures at the terminator (Bougher et al., 2015).
In this paper we present an improved version of the LMD
Venus GCM (Lebonnois et al., 2012, 2010), now extended up to
approximately 150 km. We aim here at exploring the mechanisms at work from 90 to 150 km altitudes, focusing mainly on
the thermal structure. The inclusion of the main processes contributing to the thermal balance of the upper atmosphere of Venus(i.e IR heating by CO2 , IR CO2 15 μm cooling, extreme UV heating, thermal conduction) and the coupling with a photochemical model (Stolzenbach et al., 2014), provide a fully self-consistentand quantitative description of the Venusian atmosphere. Those
improvements make the LMD-GCM the only existing ground-tothermosphere 3D model. In addition, a non-orographic gravity
waves parameterisation following Lott and Guez (2013) is includedin the LMD Venus GCM for the ﬁrst time.
Section 2 describes the model and the physical processes implemented in this version. The results of the reference simulation,in terms of thermal balance and thermal structure variation arepresented in the Section 3, with additional details on the composition given in Section 6. Comparison to observational data by Venus
Express mission and ground based telescopes is done in Section 4.
Section 5 discusses the sensitivity of the thermal structure to keymodel parameters, then conclusions are given in Section 7.
2. LMD-VGCM: model description and improvements
The Venus GCM developed at Laboratoire de Meteorologie
Dynamique (Lebonnois et al., 2010), hereinafter LMD-VGCM, isbased on the tools and experience gathered for the GCMs of the
Earth’s atmosphere (Hourdin et al., 2006) and Mars atmosphere(Forget et al., 1999). The dynamical core is a common part for
the three terrestrial planets, while the physical core is speciﬁc
for the conditions and characteristics of each planet. The model
solves the primitive equations of hydrodynamics in a sphere, using a ﬁnite difference discretisation scheme which conserves bothpotential enstrophy for barotropic non divergent ﬂows, and total
angular momentum for axisymmetric ﬂow. In the top four lay-

G. Gilli et al. / Icarus 281 (2017) 55–72
Table 1
List of chemical species included in the model. List of chemical
species included in the model. (l) stands for liquid phase, (g) for gas
phase.
1
2
3
4
5
6
7
8
9
10

CO2
CO
H2
H 2 O(g)H 2 O(l)O(1 D)OO2
O2 (1 g)
O3

11
12
13
14
15
16
17
18
19
20

H
OH
HO2
H2 O2
Cl
ClO
Cl2
HCl
HOCl
ClCO

21
22
23
24
25
26
27
28
29
30

ClCO3
COCl2
S
SO
SO2
SO3
S2 O2
OCS
HSO3
H2 SO4 (g)

31
32
33
34

H 2 O(l)
S2
ClSO2
OSCl

ers of the GCM a sponge layer is included, with Rayleigh friction
to relax the horizontal winds and potential energy to the mean
zonal values. The Rayleigh friction is not used here as a generaldrag for the upper atmosphere. In the atmosphere of Venus theconservation of angular momentum is an important issue. In the
conﬁguration used in the LMD-VGCM the conservation of angular momentum is quite good, with numerical residuals around 10%of the amplitude of the surface physical forcings, as described in
Lebonnois et al. (2012). Compared to Lebonnois et al. (2010) (hereinafter L10), the main change in the physical parameterisation of
the lower atmosphere is the planetary boundary layer scheme.
The vertical eddy diffusion coeﬃcient and the time evolution of
the mixed variables are now computed following the Mellor and
Yamada (1982) scheme, used also for the other planetary GCMs
developed at LMD (Earth, Mars and Titan), and described in the
Appendix B of Hourdin et al. (2002). The time step is about 1.5
min, as in (Lebonnois et al., 2010). Concerning the radiative transfer below 100 km, the heating rates for the short-wave radiation
are taken from look-up tables based on the work by Crisp (1986).
Then the rates are interpolated for the correct solar zenith angle for each GCM grid point, to include the diurnal cycle. For the
long-wave radiation a radiative transfer code based on net exchange rate matrices is used, following the formalism developed by
Eymet et al. (2009). This version of the LMD Venus GCM was previously used to compute self-consistently the thermal structure of
the atmosphere of Venus and to interpret VEx observations aroundthe cloud tops (Grassi et al., 2014; Migliorini et al., 2012). However, the complex and dynamically active upper regions unveiled
by VEx demanded modelling efforts to include dominant radiative and chemical processes above about 100 km. Moreover, several modelling improvements have been done in the last 5 years
for the LMD-VGCM. For instance, new physical packages are under development: a microphysical model dedicated to study thecloud system by Määttänen et al. (2015), Guilbon et al. (2015) anda mesoscale module (Lefevre et al., 2016 (submitted)) (both not included in this version of the VGCM). In this section we describe
the model improvements focused on the vertical extension of the
LMD-VGCM up to the lower thermosphere (up to 150 km approximately).
2.1. Vertical extension
The ﬁrst step required for the vertical extension was to add
pressure levels to the previous version of the VGCM based on 50
hybrid coordinates (i.e p = a p + b p ps with ps the surface pressure
and ap and bb listed in Table 1, L10). A total of 28 pressure levels were added in order to reach thermospheric altitudes of about
150 km. They have been calculated through the hydrostatic equaz −ztion, P (z2 ) = P (z1 )exp(− 2 H 1 ), with z1 < z2 , where H is the scaleheight H = RT /gM, R the perfect gas constant (R = 8.32 J mol−1K−1 ), g the gravitational constant (g = 8.87 m s −1 ) and M themean molecular mass for Venus (M = 43.55 g/mol). Note that M is

57

explicitly calculated in the physics of the model. The vertical sampling is approximately 2–3 km between 100 and 150 km. The
ﬁnal grid consists of 78 vertical levels and a horizontal resolutionof 48 longitudes × 32 latitudes (7.5° × 5.625°).
The following step was to include into the GCM the physicalprocesses relevant to the upper atmosphere. At the low gas densities typical of upper atmospheric regions the impact of non-LTE effects is well known by a number of theoretical models of terrestrial
planet atmospheres (Fox and Bougher, 1991; Keating et al., 1980;
Lopez-Puertas and Taylor, 2001). In these non-LTE models there arenormally two separate terms of importance at thermospheric altitudes which require speciﬁc treatment, the IR cooling by CO2 at
15 μm and the solar heating by CO2 in its near-IR bands (1–5 μm).
However, a comprehensive on-line implementation of a complex
non-LTE model into a GCM is not straightforward and proper simpliﬁcations are required. The strategy adopted in this work was tofollow the parameterisations already developed for the LMD Mars
GCM, computationally compatible with 3D GCM simulations. The
Martian GCM is a valuable precursor for the Venus GCM since bothplanets share many radiative and chemical processes related to
CO2 atmospheres. The heating and cooling rates due to the radiative transfer in gaseous CO2 atmospheres at LTE and non-LTE conditions are calculated in both GCMs with separate parameterisations (see next sections for non-LTE parameterisation description).
To combine the radiative tendencies in the transition region between the LTE and non-LTE regimes, the proﬁles of heating/coolingrates are smoothly merged as it follows:

Q = α QNLT E + (1 − α )QLT E
pwith α = 1/(1 + pmerge)4 . Taking into account that on Venus thenon-LTE effects dominate above about 0.01 Pa (Gilli et al., 2009;
López-Valverde et al., 2011; Roldán et al., 20 0 0), this transition region is centered below that level. For the thermal cooling, pmerge
has been ﬁxed to 0.5 Pa. On another hand, for the solar heating
we ﬁrst veriﬁed a good match between the two radiative tendenQ
cies ( Q LT E = 1) below about 10 Pa, and then we ﬁxed pmerge =NLT E

102 Pa to get a soft transition.
Above pressure 0.01 Pa (about 120 km altitude), the tendencies calculated by the GCM are purely a sum of non-LTE contributions. The processes implemented in the vertical extension are
described in detail below and the heating/cooling rates calculatedby the model are shown in Figs. 1 and 2.
2.2. NIR solar heating rate
The absorption of solar radiation in the CO2 IR bands has a considerable impact on the thermal structure of Venus atmosphereabove 90–100 km. Detailed calculations performed by a full nonLTE model (Roldán et al., 20 0 0, hereinafter RO20 0 0) showed thatthe net solar radiation absorption depends mainly on the densityof the atmosphere, and to a smaller degree on the solar zenith
angle, thermal structure and atomic oxygen abundance. This assumption permits tabulation of correction factors to convert LTE
radiative transfer calculation rates into realistic non-LTE heating
rates. This tabulation was done for the Mars’ atmosphere by LopezValverde and Lopez-Puertas (1994), and it allowed to calculate
the heating rate in the Martian GCM using a 1-parameter formula (Eq. 1) that mimics the heating rates calculated by line-byline non-LTE simulations, as shown in Forget et al. (1999) and
González-Galindo et al. (2009). The essential aspects of the nonLTE 1D modelling are similar on Mars and Venus (Gilli et al., 2011;
López-Valverde et al., 2011): this analogy enables to implement asimilar non-LTE parameterisation on both planets. We assume theorbit of Venus circular r = r0 with a mean Venus-Sun distance r0= 0.72 AU, the axis perpendicular to the ecliptic plane, and no seasonal variation. At a pressure p, the solar heating is computed by

58

G. Gilli et al. / Icarus 281 (2017) 55–72

Fig. 1. Heating and cooling rates [K/s] of the reference simulation (after 10 Venus days) by the Venus LMD-GCM in the equatorial region (10°S–10°N): at midday, midnight,
morning and evening terminators, as indicated. Approximate altitudes axis is showed on the right panels, as reference.

Fig. 2. As in Fig. 1 but for mid-latitude regions (30°N–50°N).

the following function:


2

r
∂T
∂T
( p, r, μ ) =( p0 , r0 , 0 ) × 02
∂t
∂t
r




p −b

p0
μ˜ 1 + 1pp(1)

where at a pressure p0 = 1320 Pa and zero solar zenith angle μ =0 the heating rate per Venusian day (Vd) is ∂∂Tt ( p0 , r0 , 0 ) = 18.13
K/Vd (following precomputed tables for the solar ﬂux and heat-ing for Venus by Crisp (1986), p1 is the approximate pressure below which non-LTE effects are signiﬁcant ( p1 = 0.008 Pa ), andb = 1.362. The cosine of the solar zenith angle corrected for atmospheric refraction is μ˜ = [(1224μ2 + 1 )/1225 )]1/2 . Several tests
were performed to ﬁne-tune those parameters in order to reproduce the essentials of RO20 0 0. It had to be considered that therelative abundance of atomic oxygen is variable in the GCM, while

G. Gilli et al. / Icarus 281 (2017) 55–72

it was ﬁxed to VTS3 values in RO20 0 0. Also, in the ﬁne-tuning process particular attention was paid to get a good match in the lower
mesosphere (in LTE conditions) with the VGCM solar heating rates
based on Crisp (1986).
2.3. Thermal infrared
The thermal cooling by the CO2 bands at 15 μm plays a keyrole to maintain the observed temperature structure and to moderate solar EUV ﬂux changes during the solar cycle. Comprehensive
studies of LTE breakdown of the main CO2 levels helped explaining the long-standing problem of the cold lower thermosphere of
the planet observed by Pioneer Venus (Schubert et al., 1980) between 10−2 and 10−4 Pa and the small amplitude of its temperature variation during the 11-year solar cycle (Bougher et al., 1997;
Sharma and Wintersteiner, 1990). The collisions with atomic oxygen are especially effective in exciting the vibrational levels of CO2
molecules, thus enhancing the cooling at thermospheric heights
where non-LTE prevails. The effectiveness of this enhancement
process depends upon both the relative O density at a given pressure and the CO2 -O quenching rate coeﬃcient kv−t . On both Marsand Venus the collision of atomic oxygen with CO2 molecules produces vibrationally excited CO2 (O + CO2 → CO2 ∗ + O) which decays and radiates 15 μm photons through an optically thin path tospace, thus cooling the atmosphere. Both quantities ([O] and kv−t )are not very well known in the Venus atmosphere.
Regarding the [O] abundance, our knowledge is based only onnighttime mass spectrometer measurements down to ∼135 km by
ONMS/Pioneer Venus (for solar mean condition) (Kasprzak et al.,
1993; Niemann et al., 1980) and indirect observation of the O2night-glow by VIRTIS/Venus Express (Gérard et al., 2013; Soretet al., 2014). No daytime measurements are available in the lower
thermosphere of Venus, since PVO measurements were limited toaltitudes above ∼150 km. In this work, the atomic oxygen mixing ratio is calculated on-line by a coupled photochemical scheme(described in Section 2.7). Atomic oxygen is also vertically diffused
and transported by the general circulation.
Laboratory and atmospheric measurements on the Earth indicate values for kv−t between 1 × 10−12 cm3 s−1 and 6 × 10−12
cm3 s−1 . This uncertainty is still an open and important question for planetary atmosphere studies (Feoﬁlov and Kutepov, 2012;
Huestis et al., 2008; Sharma, 2014). The CO2 -O deactivation rate
(kv−t ) used here is a “median” value commonly used in terrestrial atmosphere models (3×10−12 cm3 s−1 ) (e.g. Bougher et al.,2015; González-Galindo et al., 2013; Brecht et al., 2011). We followed the same strategy developed for the LMD-MGCM (GonzálezGalindo et al., 2013). Our simpliﬁed non-LTE model takes into account the 5 strongest ro-vibrational bands of CO2 (instead of 92 forthe complete non-LTE model) and allows for a variable abundanceof atomic oxygen.
2.4. EUV heating
Extreme UV (EUV) absorption is also parameterised following
the scheme developed for the Mars GCM, described in GonzálezGalindo et al. (2013, 20 09, 20 05) (hereinafter FGG20 05). Fundamentally, the heating rates are obtained using pre-calculated
photo-absorption coeﬃcients. These coeﬃcients, integrated in 36
spectral subintervals of varying spectral width in the 1–800 nm
spectral range, were tabulated as a function of the column amount.
From this tabulation, the total photo-absorption coeﬃcients foreach species are calculated by interpolation to the actual column
amount and by a sum over the 36 spectral intervals. The abundance of the chemical species used in this parameterisation is selfconsistently calculated in the VGCM through the coupled photochemical model.

59

This implementation is suitable for the VGCM since both Martian and Venusian thermospheres are dominated by CO2 and its
dissociation products O and CO, as well as N2 and O2 . N2 is not
yet included in the photochemical scheme used here: the N2 volume mixing ratio is ﬁxed at 3.5 % below the homopause. We take
into account the absorption of solar radiation by the most abundant species (CO2 , CO, O and O2 , and N2 ). To include the variationof the UV solar ﬂux with the 11-year solar cycle, we also followed
the previous calculations by FGG2005 which reproduce the solar
cycle 22 and 23 by a sinusoidal ﬁt, and corrected it to Venus-Sun
distance. The ﬁrst year (2006) and the last year (2014) of the VEx
mission are close to the minimum and the maximum of the solar
cycle 24, respectively. In order to compare consistently our resultswith VEx measurements taken during those 8 years of the mission,
we chose a solar ﬂux value representative of “mean” solar conditions.
The EUV eﬃciency is another uncertain parameter in the modelling of terrestrial planets (Huestis et al., 2008). It is customary toadopt a nominal value and test for uncertainty range around thatvalue. The range used in this work is [16–28 %], in agreement withdetailed on-line calculations provided by Fox (1988), and compatible with acceptable 15 μm cooling rate kv−t parameters.
2.5. Thermal conduction and molecular viscosity
The thermal conduction is one the main cooling mechanisms
that, together with the CO2 15 μm cooling, balances the UV heating in the thermosphere of Venus. This process is governed by thefollowing equation:



∂T
1 ∂
∂T
=
k
∂t
ρ cp ∂ z ∂ z(2)where T is the temperature (K), ρ the density (kg/m3 ), cp the speciﬁc heat and k = AT 0.69 is the thermal conduction coeﬃcient (Jm−1 s−1 K−1 ), with a number density weighted average of the individual thermal conductivities.
A similar expression is used for the molecular viscosity, whose
effect is to smooth velocity gradients as given below:∂S 1 ∂∂Sμ=∂t ρ ∂ z∂z(3)

where S stands for the component of the horizontal wind (m s−1 )and μ is the coeﬃcient of molecular viscosity (kg m−1 s−1 ) that is
related to the thermal conduction coeﬃcient by k = 14 [9c p − 5(c p −R ) ] μ.
Here we followed the same implicit numerical scheme as described in Forget et al. (1999), in which both equations are discretised, allowing for longer and more stable time stepping.
2.6. Molecular diffusion
The thermosphere is also characterised by a transition region
called the homopause, below which atomic and molecular constituents are well mixed by winds and turbulence, and above
which the individual species begin to diffusively separate according to their individual masses or scale heights. Above that region
the molecular diffusion dominates the vertical distribution of different species, leading the lightest atomic species (CO, O, He, H) to
prevail up to the base of the exosphere.
In our model we use the algorithm developed by
Chaufray et al. (2015) which describes the behaviour of chemical
species by assuming that the diffusion velocity of each species is
independent of the other species’ diffusion velocity. The diffusive
equilibrium condition is assumed at the upper boundary. The
same procedure was previously implemented in the LMD-MGCM(see González-Galindo et al., 2013, for further details).

60

G. Gilli et al. / Icarus 281 (2017) 55–72

2.7. Photochemical model
The version of LMD-VGCM presented in this work includes a full photochemical model of the Venus atmosphere(Stolzenbach et al., 2015) which is based on the heritage of theMars photochemical model (Lefèvre et al., 20 08, 20 04). It provides a comprehensive description of the chemistry of carbon dioxide, oxygen, hydrogen, chlorine and sulfur compounds. The chemical module includes 34 species, listed in Table 1, which interact through 19 photodissociation and 102 chemical reactions. The
photodissociation rates are calculated with a modiﬁed version of
Tropospheric Ultraviolet and Visible (TUV) model (Madronich andFlocke, 1999) that takes into account the solar zenith angle, theoverhead column of CO2 and SO2 and the light scattering by Venusian aerosols to calculate the photodissociation rates. The optical depth from the clouds used is taken from Crisp (1986). The
chemical reaction rates are taken from the compilation of kinetics data by Sander et al. (2011). The chemical model includes aparametrised cloud representation based on a ﬁxed multimodal
droplet radius size distribution (Knollenberg and Hunten, 1980).
The cloud module also calculates at each physical time step thecomposition of the binary droplets of H2 SO4 /H2 O as well as thesedimentation rate of the condensed species, which is a functionof aerosol size.
This 3D photochemical model is able to reproduce key observations like the strong H2 O and SO2 decrease in the cloud region(Encrenaz et al., 2016). The implementation of a photochemical
model in the LMD-VGCM enables the study of the horizontal distribution and vertical transport of chemical species for the ﬁrst time
up to 150 km.

2.8. Gravity waves parameterisation: motivation and preliminary
implementation
The role of gravity wave (GW) propagation on the perturbation of density, temperature and cloud structure of the Venus atmosphere is commonly invoked in the literature to explain theobserved variations. Gravity wave activity in Venus middle atmosphere is conﬁrmed by several measurements, although theirsource is still unknown. Speciﬁcally, short time phase ﬂuctuations
were seen in radio occultation data by Pioneer Venus mission,
and were presumably caused by gravity waves (Woo et al., 1980).
More recently, VIRTIS on board VEx observed visible trains of oscillating cloud brightness in the UV for the upper cloud layer(∼70 km) on the dayside and around 50 km altitude from thethermal radiation on the nightside (Peralta et al., 2008). These authors detected wavelengths of 60–150 km propagating westward
with phase velocities similar to the zonal ﬂow, with no correlation between the location of the waves and surface topography,
latitude, local time or structure of the wind. On the contrary, astrong latitudinal dependence of the gravity wave activity was observed by VeRa Radio science experiment on board VEx, with anincreasing activity with latitude, probably associated with topography (Tellmann et al., 2012). Also, those results show evidenceof local time dependence in the equatorial region (stronger GW
activity during daytime), supporting the hypothesis of cloud-level
convection as a source of GW in this latitudinal range. Horizontal wavelengths observed by VeRa at 65 km are on the order of
300–450 km, with horizontal phase speeds of roughly 5–10 m/s.
In addition, Garcia et al. (2009) detected perturbations in the CO2
non-LTE emissions expected in the altitude range of 110–140 km
by the instrument VIRTIS/VEx. This study reported wave-like structure with horizontal wavelengths of 90–400 km, claimed to begenerated from the polar vortex. High latitude (60°-80°) waves ofdifferent type (long, medium, short and irregular) have also been

Table 2
Wave characteristics used in the reference simulation. Values in the bracket indicate
the extreme of the probability distribution.
Phase Velocity
[m/s]
[1–61

Horizontal wavel.
amplitude [km][50–500]

EP-ﬂuxkg m−1 s−2[0–0.005]

Saturation
0.85

Diffusion
kg m−1 s−1
0.1

identiﬁed by the VMC camera on board VEx and assumed to be
gravity waves (Piccialli et al., 2014).
Several general circulation models (Bougher et al., 1988;
Hoshino et al., 2013; Zalucha et al., 2013; Zhang et al., 1996)
also have investigated the effect of gravity waves on the variation of wind velocity in the Venusian lower thermosphere.
Bougher et al. (1988) assumed Rayleigh friction and suggested
that wave drag would decelerate the SS-AS winds in the thermosphere. In Zalucha et al. (2013), the authors suggested that gravity
waves did not propagate above an altitude of about 115 km at
the terminator, because of total internal reﬂection. In the work by
Hoshino et al. (2013) the source of gravity waves causing perturbations in the thermosphere is supposed to be at 80 km (near the
lower boundary layer of their GCM). They concluded that the westward momentum which drives the retrograde superrotation zonal(RSZ) winds in the nightside around 110 km is provided by gravity
waves breaking between 115–130 km altitude.
Encouraged by observational evidence of GW and their potential impact on the upper mesosphere/lower thermosphere ofVenus, we implemented a preliminary parameterisation of non stationary (i.e non-orographic) gravity waves for the ﬁrst time in
the LMD-VGCM. At this stage of the work, the main purpose isto reduce the strong jets obtained before the morning terminator
at equatorial and low-latitudes and perform more stable runs. In
practice, the model is unstable and it stops running when numerical instabilities at some points of the grid produce unrealistic temperatures/densities.
Non-orographic gravity waves on Venus are thought to be emitted above the thick convective layer of the middle cloud region (50–60 km), propagate upwards and break in the thermosphere providing a signiﬁcant source of momentum and energy (Alexander et al., 1992). In this work we followed the formalism developed for the Earth GCM and described in Lott andGuez (2013) based on a stochastic approach, where a large ensemble of monochromatic GW is generated at a ﬁxed launching altitude by emitting a few waves at each time step, and by adding theeffect of these waves to that of the waves launched before. The
source of the gravity waves is here chosen uniform, without latitudinal variation, and ﬁxed at roughly 55 km, near the top of theconvective layer. The routine is now fully implemented in the code
and a detailed sensitivity study to tuneable values used in this GW
parameterization will be done in a future dedicated paper. Here
we provide the baseline parameters used in the current version of
VGCM. Following Lott and Guez (2013) the wave characteristics arechosen randomly, with an arbitrarily ﬁxed probability distribution.
As listed in Table 2, the absolute phase speed is between 1 and60 m/s, the horizontal wave number is between 50 and 500 km(within observational range). The maximum value of the EP-ﬂux atthe launching altitude is ﬁxed at 0.005 kg m−1 s−2 . The parameters controlling the breaking (saturation parameter) and the decayof EP-ﬂux (diffusion parameter) are 0.85 and 0.1 kg m−1 s−1 , respectively (see Equations (11) and (12) in Lott and Guez, 2013 fordetails).
3. Results: reference simulations
In this section we describe the results obtained after 10 Venus
days (about 3 terrestrial years). All the physical processes andG. Gilli et al. / Icarus 281 (2017) 55–72

61

Fig. 3. Reference averaged temperature proﬁles predicted by the LMD-VGCM at 103 –10−5 Pa (about 65–150 km approximately), at equatorial regions (10°S–10°N) (left panel),
mid-latitudes (30°N–50°N) (middle panel) and at high latitudes (70°N–80°N) (right panel) at four different local time conditions: midday (blue), midnight (red), morningterminator (green) and evening terminator (black). Approximate altitudes are also indicated in the three panels, as reference . (For interpretation of the references to colourin this ﬁgure legend, the reader is referred to the web version of this article.)

the non-orographic GW parameterisation described in the previous sections are included. The initial state is taken from previous
simulations on the 0–100 km domain, performed for nearly 10 0 0Venus days. Above 100 km, the initial horizontal wind is ﬁxed tozero and the initial temperature is set to the 100 km values andheld constant until the top of the model.
For the reference simulation we used: 3×10−12 cm3 s−1 for theCO2 -O deactivation rate (kv−t ), 23 % for the EUV eﬃciency, close
to the upper limit of standard eﬃciencies calculated by Fox (1988),
and solar cycle “mean” ﬂux values (to compare our results with
Venus Express data taken between 2006 and 2014). Sensitivity
tests showing the impact of the variation of those key parameterson the thermal structure are shown in Section 5. We describe here
the thermal structure predicted by the LMD-VGCM above 100 km,
in terms of vertical, diurnal and horizontal variation.
3.1. Thermal balance
The thermal structure predicted by the LMD-VGCM reﬂects the
combined effects of radiative, chemical and dynamical processes.
Figs. 1 and Fig. 2 show the components of the thermal balance
averaged in equatorial (10°S-10°N) and mid-latitude (30°N–50°N)
regions, respectively, at four different local times: midday, midnight, morning terminator (MT) and evening terminator (ET). As
discussed previously, during daytime, and above 10 Pa (approximately 95 km) the heating occurs via solar absorption by CO2
near-infrared bands (2.7 μm, 4.7 μm and 1–2 μm), and above
10−3 Pa (130 km approximately) EUV absorption dominates. Thermal cooling takes place via CO2 non-LTE transitions around 15 μm
which, together with thermal conduction, competes with the heating terms to control the temperature in the lower thermosphere.
The net radiative rate is balanced by the heating/cooling terms obtained from the dynamical core (adiabatic ﬂuid advection, including both horizontal and vertical transport).
At midday and at equatorial regions (Fig. 1), the mean heating rates peak around 0.027 K/s (about 2300 K/Earth day) and
0.058 K/s (about 50 0 0 K/Earth day), for the NIR and EUV rates,
respectively. Smaller values are found at high-latitudes (not shownhere), related to the decrease of solar radiation with increasing
latitudes. The net radiative heating at midday indicates a warm

layer between 1 and 0.1 Pa (about 105 and 110 km) at all latitudes. Radiative balance largely maintains the midday structure,
while the dynamical terms are minimal compared to these radiative terms. This behaviour is also simulated by the VTGCM (e.gBrecht and Bougher, 2012). At other times of the day the net radiative heating is negative. At the terminators the thermal balance suggests a colder thermosphere above 10−3 Pa (approximately
130 km), since the solar ﬂux is weaker and hence the EUV heating
is not strong enough to balance the thermal conduction. In contrast to daytime, the dynamical heating at the terminator dominates above the other radiative terms. At nighttime, the dynamical rate is slightly stronger at the equator than at higher latitudes,
as shown in Figs. 1 and 2. This suggests a link with the vertical
downwelling of subsolar-to-antisolar (SS-AS) circulation air at theanti-solar (AS) point. The subsidence, together with the transport
from the warm dayside (advection), heats the atmosphere around
0.5 Pa (about 105 km). In addition, radiative cooling is weak atnighttime at those layers.
Those NIR heating/cooling rates are compatible with the results by the VTGCM of Brecht et al. (2011) (see their ﬁgure 5) andBrecht and Bougher (2012) for similar latitude-local time conditions. However, the EUV heating rates predicted at daytime by ourmodel are larger by about a factor 2. This is due to the different
solar cycle conditions used in the models. These authors used solar minimum ﬂuxes (F10.7=70) to compare their simulations withearly VEx data, while our simulations assume mean solar conditions. It should also be noted that both models follow RO20 0 0 for
the 4.3 μm heating and 15 μm cooling rates, although the parameterisation of these non-LTE processes is different in the models.
Moreover, the VTGCM described in Brecht et al. (2011) does notinclude the coupling between deep and upper regions of the atmosphere, as it is the case in our model. This may play an important role on the whole thermal structure and could partially explain the differences between both models found at the altitudes
considered in this paper. In spite of these differences, the coolingand solar heating terms are similar. Nevertheless, a proper modelinter-comparison has to be carried out to better explain those discrepancies/similarities between GCMs, as it was done for the Martian GCMs in (González-Galindo et al., 2009).

62

G. Gilli et al. / Icarus 281 (2017) 55–72

Fig. 4. Local time-pressure cross sections of temperature ﬁelds [K] of the reference
VGCM simulation after 10 Venus days, averaged at equatorial region 30°S–30°N (upper panel), and high northern latitudes 70°N–80°N (lower panel). Approximate altitudes are also indicated in the two panels, as reference.

3.2. Vertical structure and diurnal variation
The temperature structure of the upper mesosphere/lower thermosphere of Venus is strongly affected by the net day-to-nightvariation of solar ﬂux, and the vertical structure is therefore strictly
linked to the diurnal variation. Fig. 3 shows vertical temperature
proﬁles of the reference simulation at four selected local times (0,6, 12, 18 LT) in the pressure range 103 –10−5 Pa (corresponding to
65–150 km, approximately). The three panels of the ﬁgure represent the averaged proﬁles at equatorial regions (10°S–10°N), at
middle (30°N–50°N) and at high (70°N–80°N) northern latitudes.
Similar results are obtained at southern latitudes (not shown).
A local warm region is found at dayside mesospheric layers between 1 and 10−1 Pa (approximately 105 and 110 km) due to theimpact of near IR heating rates during daytime. This warm layer
corresponds to the net radiative heating peak shown in Figs. 1 and
2. On the nightside, a warm region around 0.5 Pa (approximately
105 km) whose intensity decreases with increasing latitudes, ispredicted by the LMD-VGCM. It is interpreted here as a combination of advection from the dayside and subsidence at the AS-point(see Section 3.1).

Above this warm layer a colder region is produced by 15 μm
non-LTE CO2 bands. In the thermosphere the temperature increases
due to EUV solar absorption. Also, a cold nightside thermospherictemperature of less than 150 K is found, consistent with the nightside “cryosphere” above 140–150 km observed by Pioneer Venus(Keating et al., 1980; Schubert et al., 1980).
In general, the predicted vertical proﬁle of temperature is inagreement with recent observations by VEx (Gilli et al., 2015;
Mahieux et al., 2015; Piccialli et al., 2015) as shown with more details in Section 4, and with other modelling efforts in the upper
atmosphere of Venus (Bougher et al., 2015; Brecht and Bougher,
2012; Brecht et al., 2011).
A further insight into the vertical and diurnal variations of thereference LMD-VGCM simulation is gained by analysing the local
time maps of latitudinally averaged temperatures versus pressurelevels, shown in Fig. 4. The two panels correspond to the averagetemperatures between 30°S and 30°N and at high northern latitudes, between 70°N and 80°N. The nightside warmer layer between 1 Pa and 10−1 Pa (approximately 105–110 km) at equatorial
regions is about 30–40 K warmer than at higher latitudes, as expected by strong subsidence at the AS point. This is also clearlyvisible in the upper panel of Fig. 5a, in which temperature values
at 1 Pa are plotted versus local time, for 6 selected latitudes in thenorthern hemisphere (0°, 16°, 33°, 50°, 67°and 84°). At those levels (∼105–110 km), nighttime temperature values at the equatorreach about 240 K, which are consistent with the permanent warm
mesospheric layer observed by SPICAV/VEx at lower altitudes (90–110 km) (Piccialli et al., 2015). A more detailed discussion is provided in Section 4.
The lower panel of Fig. 5a shows the variation of the ratio [O/CO2 ] with local time to investigate possible impact on thethermal structure. Two minima of temperature (about 190 K and
195 K) predicted at 4–5 LT and at 19–20 LT, just before the morning terminator and after the evening terminator, respectively, are
not linked to this ratio but to upward transport and strong zonal
winds, as it will be detailed in Section 3.5. Fig. 5a also shows astrong enhancement of atomic oxygen at night, peaking at midnight. This is produced by the eﬃcient day-to-night transport andit corresponds to the subsiding region (i.e local warm region).
Similarly, Fig. 5b shows the diurnal evolution of the temperatureand the ratio [O/CO2 ] near the top of the model (approximately
140 km). The LMD-VGCM predicts a nightside equator-to-polewarming (up to 15 K) at those altitudes. Figs. 1 and 2 suggest that
this trend is due to dynamical heating in the lower thermosphere(in the absence of other heating processes at nighttime). This emphasises the large role of dynamical processes in the evolution of
the thermal structure and composition in the upper atmosphere ofVenus, notably throughout nighttime and terminators, while during daytime radiative processes dominate at all altitudes. This importance of dynamical and radiative processes is also simulated by
the VTGCM (Brecht and Bougher, 2012).

3.3. Latitudinal variation
Fig. 6 shows the mean temperature ﬁelds as function of pressure and latitude, for different local time conditions. Temperature
is averaged for each Local Time and pressure level (i.e in a Solarﬁxed reference frame) and the associated standard deviation (notshown here) is in general less than 1K, except in the upper mesosphere (∼100–120 km) before dawn, where it is up to 2K. A further insight into the latitudinal variation is given in Fig. 7, showing the local time-latitude cross sections at six reference pressurelevels (10−4 Pa, 10−3 Pa, 10−1 Pa, 1 Pa, 5 Pa, 10 Pa), covering
the upper layers of the atmosphere of Venus studied in this work(from 95 to 140 km altitudes, approximately).

G. Gilli et al. / Icarus 281 (2017) 55–72

63

Fig. 5. Panels (a) Diurnal variation of reference temperature (upper panel) and ratio [O/CO2 ] (lower panel) at a constant pressure (1 Pa) corresponding to 105 km approximately, at six selected latitudes in the northern hemisphere: equator (green), 16°(blue), 33°(black), 50°(yellow), 67°(pink), 84°(purple). Panels (b) Same as in (a) but at theconstant pressure level 10−4 Pa (140 km approximately). (For interpretation of the references to colour in this ﬁgure legend, the reader is referred to the web version of thisarticle.)During daytime (7–17 LT), in general the model predicts a decreasing temperature from the SS point to higher latitudes at all
altitudes in the upper mesosphere/lower thermosphere. This trend
is in good agreement with VIRTIS/VEx daytime nadir observations
described in Peralta et al. (2016) and with other model results by
Brecht and Bougher (2012), although the intensity of the maximum around the SS point differs, being up to 40 K larger in ourmodel than in Brecht and Bougher (2012). However, this value isin agreement with the observational range obtained by groundbased telescopes as reported in Sonnabend et al. (2010) and
Krause et al. (2016). Our predicted local maximum between 10.1 Pa at daytime peaks around 280 K at the SS point and it decreases by about 50 K near the poles.
At terminators (left panels of Fig. 6) the latitudinal variation of
the temperature above 10−1 Pa (110 km approximately) is rather
small, in overall agreement with SOIR/VEx observations and with
the VTGCM results (Bougher et al., 2015; Mahieux et al., 2015).
Maps in Fig. 7 show that in the mesosphere, around 10 Pa (panelf), a colder region is found at mid-low latitudes (40°S–40°N) near
the terminator. This structure evolves and propagate upwards ac-

cordingly to the strength of the zonal winds at those altitudes, as
discussed in Section 3.4.
At nighttime, the warming predicted at mesospheric altitudes(see previous sections) is stronger at the equator than at higher
latitudes, conﬁrming a strong subsidence at the anti-solar point.
Above 10−1 Pa the latitudinal trend is inverse: the atmosphere is
up to 15 K warmer at polar regions than at equator latitudes (seealso panels a, b and c in Fig. 7). Unfortunately, this variation cannot be validated yet with available measurements: nighttime observations by SPICAV/VEx are scarcer at polar regions compared toequatorial measurements.
3.4. Impact of gravity waves parameterisation on the thermal
structure: a preliminary study
Fig. 8 shows the impact of our parameterisation of nonorographic gravity waves (see Section 2.8) on the thermal structure
and winds predicted by the LMD-VGCM around 105 km (1 Pa).
Simulations performed without non-orographic GW are comparedwith those obtained including this parameterisation (after 3Vd and
10Vd). The tuneable GW parameters used here are listed in Table 2.

64

G. Gilli et al. / Icarus 281 (2017) 55–72

Fig. 6. Latitude-pressure cross sections of reference simulated temperatures [K] at morning terminator (6 LT), midday (12 LT), evening terminator (18 LT) and midnight (0LT), as indicated. The pressure levels (102 –10−5 Pa) correspond to 80–150 km altitude, approximately.

The cold region visible between 2–5 LT in the temperature ﬁeld in
the left column (no GW) is produced by strong equatorial jets andvertical ascending ﬂux from the lower to the upper mesosphere.
These jets are noticeably reduced after the non-orographic GW
scheme implementation (second and third column), notably comparing simulations after the same day (3 Vd). Both zonal and vertical winds are smoothed at mid-high latitudes by the GW propagation, and consequently the formation of the cold upwelling pool isdelayed and weakened (see results after 10 Vd). Its origin and evolution remains to be explored in a dedicated paper but it is either
linked to angular momentum transport by the thermal tides, or toplanetary scale waves propagation from the lower mesosphere. The
current GW parameterisation, even if preliminary, helps reducing
the temperature differences between MT and ET, which are not observed in the data.
A number of sensitivity tests have been performed to evaluate
the impact of GW tuneable parameters on the thermal structure, inparticular those linked to the waves characteristics, still not wellconstrained by observations. For example, we found that GW aremostly absorbed where zonal winds are stronger (i.e near the terminator). Gravity waves slightly accelerate the zonal wind ﬁeld in
the region of its peak (around 103 Pa) and decelerate it at higher
layers. The higher the relative phase velocity, the larger this effect.
Moreover, larger horizontal amplitude GW deposit less momentum than smaller waves. GW with larger initial vertical momentum ﬂux are rapidly absorbed near the launching altitude and theyinject more kinetic momentum compared to waves of the samecharacteristics (phase velocity amplitude, horizontal wavenumber,launching altitude). As a consequence, they decelerate the zonal
winds in the whole atmosphere.
In conclusion, the variation of the main parameters used in the
GW routine (i.e., wave characteristics, saturation, wave dissipation

parameters, vertical momentum ﬂux) has a signiﬁcant impact inthe stability of the model, since it reduces (or ampliﬁes) the evolution of the dynamical structure discussed in previous section andits effect on the atmospheric layers above 90 km. The criteria toselect the waves characteristics values in this work was consistency
with the limited observational constrains (altitude of the source,
waves amplitude, phase velocity) and a better agreement with observed temperature values at the terminator.
3.5. The transition region: key to understand the variability of theupper atmosphere
Recent VEx observations revealed a surprisingly more variable
upper atmosphere (Gérard et al., 2014; Soret et al., 2014) notably
between 90 and 120 km, were RSZ ﬂow and day-to-night circulation overlap. Our model simulations can help to investigate the exchanges between the deep atmosphere and the upper layers, andto study the speciﬁc processes responsible for driving the variations in this puzzling “transition” region. The latitude-local timemaps at a constant pressure (z = 1 Pa) in Fig. 8 exhibit the strong
correlation between temperature, zonal winds and vertical transport in that region. This is particularly discernible in the left column of the ﬁgure, showing results without the GW implementation. In that case, strong RSZ winds (up to 250 m/s) shift the cold
air westward and consequently the vertical upward ﬂux does notmatch perfectly the cold temperature region.
We suggest a link between the angular momentum transport
by the thermal tides and the upwelling before the morning terminator. A detailed wave analysis by the LMD-VGCM presented in
Lebonnois et al. (2016) conﬁrmed the role of diurnal and semidiurnal tides in vertical angular momentum transport in the upper
clouds, as shown before in L10. Above 102 Pa (about 80 km) at
high latitudes diurnal tides dominate, while between 104 and 102

G. Gilli et al. / Icarus 281 (2017) 55–72

65

Fig. 7. Local time-latitude cross sections of longitudinally averaged temperature ﬁelds [K], after 10 Vd at six different pressure levels: as indicated. They correspond to
140 km (a), 130 km (b), 110 km (c), 105 km (d), 100 km (e) and 95 km (f) altitudes, approximately.

(60–80 km, approximately) at low-middle latitudes a combinationof semi-diurnal tides and vertical transport prevails. However, it
must be kept in mind that in the simulations described here, thezonal wind below the clouds is a factor 4–5 weaker than observed.
This certainly affects the angular momentum budget in the cloud
region, but it may also have an impact on the circulation in thelayers above, in particular in the transition region discussed here.
4. Comparison with VEx and ground-based temperatures above
100 km
In this section we compare our reference simulation, performed
after 10 Venus days, with a selection of Venus Express and ground-based retrieved temperatures. Data and simulations in Figs. 9–11
are averaged in mid-latitudes bins (30°N and 50°), at different local
time conditions. The comparison is performed in altitude insteadof pressure coordinates, since the selection measurements usedhere (i.e solar/stellar occultation, limb measurements) are taken infunction of altitudes.
For the dayside (7-17LT) the measurements are very scarce anduncertain. Fig. 9 shows available daytime retrieved temperatures
above 100 km altitude from CO IR non-LTE emissions around
4.7 μm observed by VIRTIS/VEx (Gilli et al., 2015), and groundbased observations of heterodyne CO2 emissions at 10 μm, around
110 km by THIS/HIPWAC (Wischnewski et al., 2015). As represented

66

G. Gilli et al. / Icarus 281 (2017) 55–72

Fig. 8. Latitude-local time cross section of longitudinally averaged temperature [K], zonal wind and vertical ﬁelds, without (left column) and with the inclusion of a GW
parameterisation (middle and right columns), after 3 Vd and 10 Vd, as indicated. Pressure level is 1 Pa (105 km approximately) in all panels. Positive zonal wind values
are westward, and positive vertical wind values correspond to downward motion. Temperature ﬁeld contours (black lines) are over-plotted in all panels and wind vectors (inm/s) only in the upper panels.

in the ﬁgure, VIRTIS errors bars (30–40 K) and the variability observed by ground-based telescopes are very large. Thus, it is diﬃcult to draw strong conclusions, but the model results are withinerror bars at almost all altitudes. Compared to the VTS3 reference model, both LMD-VGCM and the observed temperatures arewarmer below 125 km altitude, and colder above.
At nighttime (19-5LT), SPICAV/VEx averaged proﬁles in the latitude boxes 30°N–50°N and 30°S–50°S (updated from Piccialli et al.,2015) are shown in Fig. 10, with standard deviation bars. Also,
the VTS3 model at nighttime is plotted as reference. LMD-VGCM
temperatures are generally larger than observations and VTS3 values everywhere above 100 km. Moreover, the vertical variationof temperature in the VGCM is different from SPICAV data below
about 120 km, although the intensity of the peak is similar in
the model and in the observations (around 210 K). One explanation is that the apparent warm layer found by SPICAV between
90–100 km is instead predicted by our model at higher altitudes.
Another possibility is that our model produces too strong heat-ing around 110 km. Above 120 km, at thermospheric altitudes,
LMD-VGCM nighttime simulations are about 30 K warmer thanobserved, probably caused by a too strong dynamical heating in
the model. Gravity-wave induced drag might produce the deceleration of the day-to-night winds necessary to isolate the dayside andnightside suﬃciently to yield the observed cryosphere (Alexander
et al., 1992; Bougher et al., 1997). Therefore, it is possible that theGW drag on the SS-AS component is not strong enough in our simulation to reduce the temperature on the nightside. A detailed sensitivity study of the model to GW parameters is foreseen in future
works.
Regarding the terminator, SOIR averaged proﬁles (6 and 18 LT)in the selected latitude boxes from Mahieux et al. (2015) andground-based observations from THIS/HIPWAC (Krause et al., 2015)
around 110 km are shown in Fig. 11, together with model predictions. Similar to daytime conditions, a succession of warm and coldlayers is found in both observed and predicted proﬁles, reﬂecting
the combined inﬂuence of radiative features and global circulation.

G. Gilli et al. / Icarus 281 (2017) 55–72

Fig. 9. Temperature vertical averaged proﬁles as function of altitude, predictedby the LMD-VGCM in the latitudinal range 30°N–50°N (solid black line) at midday (10-12 LT), with standard deviation (dashed black line), together with VIRTIS/VEx temperature retrieval results (blue circles) and error bars at the same
latitude and local time boxes. (extracted from Fig. 14 in Gilli et al., 2015).
Temperature results at approximately 110
km from ground based observations in Wischnewski et al. (2015) are shown in purple. The model VTS3 by
Hedin et al. (1983) at daytime is also plotted (red dashed line), as reference. (Forinterpretation of the references to colour in this ﬁgure legend, the reader is referredto the web version of this article.)

67

Fig. 11. Temperature vertical mean proﬁles as function of altitude, predicted by the
LMD-VGCM in the latitudinal range 30°N–50°N (solid black line) at the terminator,
with standard deviation of the mean (dashed black line), together with SOIR/VEx
temperature retrieval results (pink dashed line for the morning and solid for theevening terminator) at the same latitudes bin. SOIR standard deviations are also
represented with shaded pink area (after Mahieux et al., 2015).Temperature results
at 110 km from ground based observations (Krause et al., 2015) are also shownin purple. (For interpretation of the references to colour in this ﬁgure legend, thereader is referred to the web version of this article.)

VGCM of the uncertainty of the main model input. Figs. 12 and 13
show the sensitivity to three of the main radiative input parameters: the EUV solar ﬂux, the CO2 -O deactivation rate, and the nearIR non-LTE heating rate. Results obtained at the equatorial region(10°S–10°N) at midday (12 LT) are shown in Fig. 12, for both thermal balance and temperature, as representative of daytime conditions. Fig. 13 shows the impact of those parameters on the temperature proﬁles, at midnight (0 LT) and evening terminator (18LT), at equatorial regions. The variation of the ﬁrst two parameters(EUV solar ﬂux and quenching rate coeﬃcient) affects only the layers above approximately 125 km (10−2 Pa) for all local times. In
contrast, changes in the SHR parameterisation play a role mainly inthe mesosphere during daytime, but also in the thermosphere atnighttime and at terminator through heat transport from the dayside.
5.1. EUV solar ﬂux
Fig. 10. Temperature vertical averaged proﬁles as function of altitude, predicted by
the LMD-VGCM in the latitudinal range 30°N–50°N (solid black line) at nighttime(19-5 LT), with standard deviation (dashed black line), together with SPICAV/VEx
temperature retrieval results (in green): 30°N–50°N (dashed line) and 30°S–50°S(solid line). Standard deviations are also represented with the shaded green area
(after Piccialli et al., 2015). The model VTS3 by Hedin et al. (1983) at nighttime isalso plotted (red dashed line), as reference. (For interpretation of the references tocolour in this ﬁgure legend, the reader is referred to the web version of this article.)

In spite of the large variability exhibited by SOIR measurements(20–30 K), signiﬁcant data-model differences may be discerned.
First, our simulations are warmer (up to 30 K) between about 110
and 125 km altitude, while thermospheric temperatures predictedabove 135 km are 40–50 K colder than observed. Second, as also
noticed at nighttime, our model predicts a mesospheric peak athigher altitudes compared to observations.
Those discrepancies are investigated by means of sensitivity
studies in the next section.
5. Sensitivity tests to input parameters
Several sensitivity tests have been performed in order to evaluate the impact on the thermal structure predicted by the LMD-

Compared to the reference temperature proﬁles described in
Section 3, calculated for mean solar cycle conditions, the solar cycle maximum/minimum conditions produce an increase/decreaseof temperatures up to 25 K above 135 km during daytime, andless then 10 K at the terminator. This indicates that the solar activity alone cannot explain the 40–50 K difference between dataand our model at the terminator in the thermosphere. In addition, Fig. 13 shows that the variation of the EUV solar ﬂux has no
impact on the temperature proﬁle at nighttime. Another parameter whose uncertainty may affect the temperature in the thermosphere is the EUV eﬃciency (not shown here): a variation of 5 %produce changes of temperatures of about 10–12 K during daytime(smaller at other local times).
5.2. Collisional rate coeﬃcient kv−t
The uncertainty in the knowledge of the deactivation/collisional
rate coeﬃcient kv−t used in the 15 μm cooling could partially explain the discrepancies between VEx temperatures and model simulations at thermospheric layers. The uncertainty brackets values
from 1 × 10−12 cm3 s−1 to 6 ×10−12 cm3 s−1 . As shown in the ﬁg-

68

G. Gilli et al. / Icarus 281 (2017) 55–72

Fig. 12. Sensitivity tests to some model input parameters. Thermal balance (left panels) and vertical temperature proﬁles (right panels) at midday (12 LT), averaged in
the latitude box 10°S–10°N are shown as function of pressure (corresponding to 65–150 km altitudes, approximately). The reference simulation is plotted in solid lines.
Upper panels: variation EUV solar cycle conditions (“minimum” dotted-dashed line, “maximum” dashed line). Middle panels: variation of the quenching rate coeﬃcient kv−t(“maximum” dashed line, “minimum” dotted-dashed line). Lower panels: variation of the Solar Heating rate (“doubled” dotted-dashed line, “halved” dashed line).

ures, the smaller the rate coeﬃcient kv−t , the weaker the coolingabove 125–130 km. The increase in temperature with respect tothe reference simulation values is up to 30 K at all local times. This
test indicates that a lower value of kv−t would reduce the datamodel discrepancy and improve the comparison above 130 km at
daytime and at terminator. However, at night the radiative coolingat 15 μm competes with the dynamical and the conduction terms,
and a smaller kv−t would also produce larger temperature, therefore amplifying the data-model difference.
5.3. Non-LTE solar heating
The last sensitivity tests shown in the Figs. 12–13 (bottom panels) concerns the eﬃciency of the non-LTE solar heating by CO2
bands in the near-IR. The solar heating rate (SHR) has been multiplied/divided by a factor 2 and the effect on the thermal structure
is to enhance/reduce the intensity of the local temperature maximum mainly at mesospheric altitudes. The changes of the temperature proﬁles around the mesospheric peak are up to 70 K at daytime, and up to 40 K and 30 K, at the terminator and nighttime,

respectively. The increase/decrease in SHR only affects the magnitude of the peak of temperature, not the altitude. In other words,
a reduction of SHR with respect to our “reference” rates (based on
Roldán et al., 20 0 0) would produce a better agreement with data
for all local times below about 125 km. However, the altitude oftemperature peak at mesospheric layers does not change since it
is strictly linked to the altitude of the NIR heating peak (i.e related
to the choice of the parameters described in Section 2.2). Further
modelling efforts are required to improve this non-LTE parameterisation by a ﬁne-tuning process which fulﬁls both a good matchwith the LTE conditions and the stability of the model.
These tests also show that radiative processes play only a partial role in the thermal balance of the Venus upper atmosphere,
and that the contribution of dynamics should be better investigated to improve the comparison with the observations. On-goingwork on the impact of the GW parameterisation on winds andtemperatures is very promising and a dedicated study is foreseen
to evaluate the impact of tuneable GW parameters on the winds
and to reach a better agreement with observations.

G. Gilli et al. / Icarus 281 (2017) 55–72

69

Fig. 13. Same as in Fig. 12 but at nighttime (left column) and at evening terminator (right column).

6. CO and O density results
Although this work is mainly focused on the VGCM predicted thermal structure, other valuable insights can be obtained
from the model. Noticeably, a better understanding of modelled atmospheric variables, as winds and composition, helps
to constrain temperature ﬁelds. As also shown by GonzálezGalindo et al. (2013) in the case of Mars GCM, temperatures in theupper mesosphere/thermosphere are more consistent with observations when temporally and spatially variable atomic oxygen concentrations are used.
Our aim here is to check the composition of the upper atmosphere calculated by the LMD-VGCM. In particular, the atomic
oxygen thermosphere content is an important parameter of the
15 μm routines. For this purpose we have compared our CO andO abundances with available VEx observations. These two species
are both produced in the Venus upper layers during daytime by
CO2 photodissociation and are thought to be rapidly transferred to
the nightside by subsolar-to-anti-solar circulation (Bougher et al.,2006; Clancy and Muhleman, 1985; Lellouch et al., 1997). Theyboth serve as sensitive tracers of transport in the lower thermosphere, in particular they can be used to monitor the large

scale circulation patterns between 90–120 km (the atomic oxygen
through O2 infrared airglows analysis) (Soret et al., 2012). Simulated density proﬁles shown in Fig. 14 represent averages in theequatorial region (10°S-10°N) at three LT conditions (0, 12, 18 LT).
A selection of CO measurements obtained at similar latitudelocal time bins are also plotted. Daytime CO densities are extracted
from the work of Gilli et al. (2015) where VIRTIS/VEx non-LTE COemissions around 4.7 μm were used, while at the terminator an
average proﬁle at equatorial regions (0°-30°N) from SPICAV/VEx retrievals (Vandaele et al., 2016) is considered. Typical measurement
error bars (not shown) are a factor 10 lower than the densities, butthey increase above 130 km. Indirect atomic oxygen density measurements can be deduced from O2 night-glow emissions observed
by VIRTIS/VEX. The O proﬁle in Fig. 14 (bottom panel) is extracted
from Fig. 7 in Soret et al. (2012), and it represents -24°S latitude.
Measured O density proﬁle averaged at latitudes 10°S-10°N are not
yet available in the literature.
Atomic oxygen concentrations predicted by the LMD-VGCM
peak at nighttime around 100
km. This is consistent withthe proﬁle by Soret et al. (2012), with VTGCM results by
Brecht et al. (2012) (see their Table 2) and with the photochemical
model by Krasnopolsky (2012). The agreement with the measured

70

G. Gilli et al. / Icarus 281 (2017) 55–72

Fig. 14. Density proﬁles of CO (top panel) and O (bottom panel) predicted by
the LMD-VGCM, averaged in equatorial regions (10°S–10°N) for different local time
conditions: LT=12 (solid lines), LT= 18 (dashed lines) and LT= 0 (dotted-dashedline). CO abundances retrieved from Venus Express observations are shown forcomparison: VIRTIS/VEx measurements at midday, in the latitude range 10°S–10°N(Gilli et al., 2015) and SOIR/VEx densities at the terminator, averaged in the latitude range 0–30°(Vandaele et al., 2016). O density proﬁle (00:44 LT, −24°S) fromSoret et al. (2012) is also showed with triangles.

range 1–5× 1011 cm−3 by VIRTIS/VEx (Gérard et al., 2009) is generally good except at the anti-solar point, where our densities arelarger than observed (about a factor 2).
Focusing on CO data set below 130 km (with smaller errorbars), daytime values are up to a factor 2 larger than at terminator,
with a discernible change of slope of the vertical gradient between
110 and 120 km, as also observed at similar layers (Vandaele et al.,2016). Again, a good agreement is found with results by a photochemical model (Krasnopolsky, 2012). Nighttime CO abundances
are up to a factor 7 larger than daytime above 100 km, and anighttime mesospheric bulge is predicted by our model around
95 km, in agreement with observations by ground based telescopes (Clancy et al., 2012; Gurwell et al., 1995). In the altitude
region between 100 and 120 km model proﬁles are smaller thanobserved values, both at daytime and terminator.
This is partially an effect of thermal structure: smaller temperatures (∼40 K) than predicted lead to larger densities (about 25%) at
those altitudes. Another possible reason to explain those discrepancies is that the SS-AS component of the wind is too strong in
the model at about 100–120 km, and CO is too rapidly transported
from the dayside to nightside. This interpretation also explains thelarger LMD-VGCM nighttime temperatures at the same altitudes. A
too strong SS-AS circulation in the model also leads to an overestimation of nighttime atomic oxygen. Larger abundances of O above
100 km would cool the temperature by CO2 collisional quenching
above about 130 km. But as demonstrated in Section 5, the increasing of the CO2 -O deactivation rate (i.e increasing the O abundances by the same factor) should be much larger than a factor 2
to have a considerable impact on the temperature in the thermosphere at nighttime.
7. Summary and conclusions
Current GCMs simulations indicate that both radiative and dynamical effects play a crucial role in determining the thermal
structure of the upper atmosphere of Venus, and that previous 1Dmodels (e.g VTS3) are unable to describe the high variability observed by VEx instruments at those altitudes. The GCM used here

is an improved and extended version of the model described in
Lebonnois et al. (2012); 2010) and it is the only existing groundto-thermosphere Venus GCM. It includes the main processes contributing to the thermal balance of the atmosphere of Venus from90 to 150 km (i.e IR heating by CO2 , IR CO2 15 μm cooling, extreme UV heating, thermal conduction, molecular diffusion) as well
as a photochemical model (Stolzenbach et al., 2014) and a gravity
wave parameterisation Lott and Guez (2013).
The main radiative and dynamical features predicted by ourmodel are consistent with current literature. Focusing on the thermal structure, a succession of warm and cold layers is found,in agreement with VEx observed trends and simulations performed with the VTGCM by Bougher et al. (2015); Brecht andBougher (2012); Brecht et al. (2011). High temperatures produced
by EUV absorption above 140–150 km sit on top a cold layer
around 125–130 km, where the 15 μm cooling also plays a rolein the thermal balance. Below that region, a local maximum produced by solar absorption by CO2 IR bands during daytime is thenadvected to the terminator and to the nighttime. The net radiative heating and the dynamical heating indicate a persistent warm
layer between 100 and 120 km at all local times and latitudes.
Decreasing temperatures with increasing latitudes are predicted
during daytime at all altitudes above 100 km. Above 110 km(10−1 Pa) this latitudinal trend is inverse at nighttime, with anatmosphere up to 15 K warmer towards the pole. However, thisvariation is smaller than the errors in the available measurements
and it cannot be conﬁrmed by data comparison. At the terminator the variation is small, as also observed by SOIR/VEx. Below
110 km, the nightside warm layer predicted by our GCM around
105 km (about 0.5 Pa) is interpreted as a combination of advection from the dayside and subsidence at the AS-point, similar to that described by the VTGCM (Brecht and Bougher, 2012).
The cold thermosphere (less than 150 K) detected at night by
Alexander et al. (1992) is reproduced by our model. Daytime temperatures are in reasonable agreement with the measurements inthe thermosphere, taking into account the large measured error
bars. CO and O density results by the LMD-VGCM have been also
compared to a few observations, and an overall good agreement isfound both in order of magnitude and shape of proﬁles.
Signiﬁcant data-model discrepancies may be also discerned.
First, thermospheric temperatures (above 135 km) at the terminators are about 40–50 K colder than measured by SOIR/VEx, whilenighttime temperatures at similar altitudes are up to 30 K warmer
than observed by SPICAV/VEx. Second, the temperature maximum
in the mesosphere (between 100 and 120 km) is higher than observed by VEx (between 90 and 110 km), although in good agreement with 3D model results by Brecht and Bougher (2012).
CO density proﬁles predicted by our model are smaller than selected measurements (VIRTIS/VEx and SOIR/VEx) between 100 and
120 km, both at daytime and terminator. Simulated O density atnighttime are larger than observed (about a factor 2) above about
100 km. This may indicate that the SS-AS component of the winds
is too strong in the model at those altitudes and the CO (and O) is
too rapidly transported from the dayside to nightside. This interpretation might also explain the larger nightside temperatures inthe thermosphere.
Sensitivity tests to input radiative parameters show that a reduction of radiative cooling rate at 15 μm (an increase of temperatures) would improve the data-model comparison at daytime andterminator above about 130 km, but not at nighttime. In addition,
the solar activity alone is not able to explain the 40–50 K differences at the terminator in the thermosphere, and the variation of
the solar heating rate only affects the magnitude of the peak, but
not the altitude. Further modelling efforts are required to improve
the non-LTE parameterisation currently implemented in the LMDVGCM. In addition, these tests indicate that dynamical effects have

G. Gilli et al. / Icarus 281 (2017) 55–72

to be further studied to improve the comparison with the observations and that radiative processes have a major role in the thermalbalance of the Venus upper atmosphere only during daytime.
A ﬁrst effort to improve data-model comparison was done hereby implementing into the Venus GCM a non-orographic GW parameterisation. The comparison with measured temperatures isnoticeably improved at terminators. However, the GW drag on
the SS-AS component is still weak. A full sensitivity study to GW
parameters with the VGCM is the next step. This work is very
promising since it will help understanding the large observed variability in the upper atmosphere and its interaction with the layers
below. In addition, future work is ongoing at LMD to develop amesoscale model for the VGCM, based on the same physics as the
GCM. This will also allow a better constraint of the GW parameterisation and to take into account the variation of the distribution,
amplitude and intensity of the gravity waves, which are currently
largely uncertain.
Acknowledgements
G. Gilli has been supported by CNES postdoc grant. F. GonzálezGalindo and M.A. López-Valverde thank the Spanish MINECO for
funding support through the CONSOLIDER program “ASTROMOL”
CSD20 09-0 0 038 and through project AYA2011-23552/ESP. We also
thank E. Millour for his technical support with the model, and two
anonymous referees for constructive reviews.
References
Alexander, M.J., Stewart, A.I.F., Bougher, S.W., 1992. Local-time asymmetries in thevenus thermosphere. LPI Contribut. 789, 1–2.
Altieri, F., Migliorini, A., Zasova, L., Shakun, A., Piccioni, G., Bellucci, G., 2014. Modeling VIRTIS/VEX o2 nightglow proﬁles affected by the propagation of gravity
waves in the venus upper mesosphere. J. Geophys. Res. (Planets) 119, 2300–
2316. doi:10.10 02/2013JE0 04585.
Bougher, S.W., Alexander, M.J., Mayr, H.G., 1997. Upper atmosphere dynamics: global
circulation and gravity waves. In: Bougher, S.W., Hunten, D.M., Phillips, R.J.
(Eds.), Venus II: Geology, Geophysics, Atmosphere, and Solar Wind Environment,
p. 259.
Bougher, S.W., Brecht, A.S., Schulte, R., Fischer, J., Parkinson, C.D., Mahieux, A., Wilquet, V., Vandaele, A., 2015. Upper atmosphere temperature structure at thevenusian terminators: a comparison of SOIR and VTGCM results. Planet Sp. Sci.
113, 336–346. doi:10.1016/j.pss.2015.01.012.
Bougher, S.W., Rafkin, S., Drossart, P., 2006. Dynamics of the venus upper atmosphere: outstanding problems and new constraints expected from venus express. Planet Sp. Sci. 54, 1371–1380. doi:10.1016/j.pss.2006.04.023.
Bougher, S.W., Roble, R.G.E., Dickinson, R.E., Ridley, E.C., 1988. Venus mesosphereand thermosphere. III - three-dimensional general circulation with coupled
dynamics and composition. Icarus 73, 545–573. doi:10.1016/0019-1035(88)90064-4.
Brecht, A.S., Bougher, S.W., 2012. Dayside thermal structure of venus’ upper atmosphere characterized by a global model. J. Geophys. Res. (Planets) 117, 8002.
doi:10.1029/2012JE004079.
Brecht, A.S., Bougher, S.W., Gérard, J.C., Parkinson, C.D., Rafkin, S., Foster, B., 2011.
Understanding the variability of nightside temperatures, NO UV and o2 IR nightglow emissions in the venus upper atmosphere. J. Geophys. Res. (Planets) 116,8004. doi:10.1029/2010JE003770.
Brecht, A.S., Bougher, S.W., Gérard, J.C., Soret, L., 2012. Atomic oxygen distributionsin the venus thermosphere: comparisons between venus express observations
and global model simulations. Icarus 217, 759–766. doi:10.1016/j.icarus.2011.06.
033.
Chaufray, J.Y., Bertaux, J.L., Leblanc, F., 2012. First observation of the venus UV
dayglow at limb from SPICAV/VEX. Geophys. Res. Lett. 39, 20201. doi:10.1029/
2012GL053626.
Chaufray, J.Y., Gonzalez-Galindo, F., Forget, F., Lopez-Valverde, M.A., Leblanc, F., Modolo, R., Hess, S., 2015. Variability of the hydrogen in the martian upper atmosphere as simulated by a 3d atmosphere-exosphere coupling. Icarus 245, 282–
294. doi:10.1016/j.icarus.2014.08.038.
Clancy, R.T., Muhleman, D.O., 1985. Diurnal CO variations in the venus mesospherefrom CO microwave spectra. Icarus 64, 157–182. doi:10.1016/0019-1035(85)90084-3.
Clancy, R.T., Sandor, B.J., Moriarty-Schieven, G., 2012. Thermal structure and CO distribution for the venus mesosphere/lower thermosphere: 20 01–20 09 inferior
conjunction sub-millimeter CO absorption line observations. Icarus 217, 779–
793. doi:10.1016/j.icarus.2011.05.032.
Crisp, D., 1986. Radiative forcing of the venus mesosphere. i - solar ﬂuxes and heating rates. Icarus 67, 484–514. doi:10.1016/0019-1035(86)90126-0.

71

Encrenaz, T., Greathouse, T., Richter, M., DeWitt, C., Widemann, T., Bézard, B.,
Fouchet, T., Atreya, S., Sagawa, H., 2016. Imaging spectroscopy of venus in thethermal infrared: variability of HDO and SO2, and daily variations of the thermal structure. In: EGU General Assembly Conference Abstracts, p. 2484.
Eymet, V., Fournier, R., Dufresne, J.L., Lebonnois, S., Hourdin, F., Bullock, M.A., 2009.
Net exchange parameterization of thermal infrared radiative transfer in venus’
atmosphere. J. Geophys. Res. (Planets) 114, 11008. doi:10.1029/20 08JE0 03276.
Feoﬁlov, A.G., Kutepov, A.A., 2012. Infrared radiation in the mesosphere and lower
thermosphere: energetic effects and remote sensing. Surv. Geophys. 33, 1231–
1280. doi:10.1007/s10712-012-9204-0.
Forget, F., Hourdin, F., Fournier, R., Hourdin, C., Talagrand, O., Collins, M., Lewis, S.R.,
Read, P.L., Huot., J.P., 1999. Improved general circulation models of the Martian atmosphere from the surface to above 80 km. J. Geophys. Res. 104,
24,155–24,176.
Fox, J.L., 1988. Heating eﬃciencies in the thermosphere of Venus reconsidered.
Planet Space Sci 36, 37–46. doi:10.1016/0032-0633(88)90144-4.
Fox, J.L., Bougher, S.W., 1991. Structure, luminosity, and dynamics of the venus thermosphere. Sp. Sci. Rev. 55, 357–489. doi:10.10 07/BF0 0177141.
Garcia, R., Drossard, P., Piccioni, G., Lopez-Valverde, M., Occhipinti, G., 2009. Gravity waves in the upper atmosphere of venus revealed by CO2 nonlocal thermodynamic equilibrium emissions. J. Geophys. Res. 114, E00B32. doi:10.1029/
20 08JE0 03073.
Gérard, J.C., Saglam, A., Piccioni, G., Drossart, P., Montmessin, F., Bertaux, J.L., 2009.
Atomic oxygen distribution in the venus mesosphere from observations of o2
infrared airglow by VIRTIS-venus express. Icarus 199, 264–272. doi:10.1016/j.
icarus.2008.09.016.
Gérard, J.C., Soret, L., Migliorini, A., Piccioni, G., 2013. Oxygen nightglow emissionsof venus: vertical distribution and collisional quenching. Icarus 223, 602–608.
doi:10.1016/j.icarus.2012.11.019.
Gérard, J.C., Soret, L., Piccioni, G., Drossart, P., 2014. Latitudinal structure of thevenus o2 infrared airglow: a signature of small-scale dynamical processes in theupper atmosphere. Icarus 236, 92–103. doi:10.1016/j.icarus.2014.03.028.
Gilli, G., López-Valverde, M.A., Drossart, P., Piccioni, G., Erard, S., Cardesín
Moinelo, A., 2009. Limb observations of CO2 and CO non-LTE emissions in thevenus atmosphere by VIRTIS/venus express. J. Geophys. Res. (Planets) 114, 0–+.
doi:10.1029/20 08JE0 03112.
Gilli, G., López-Valverde, M.A., Funke, B., López-Puertas, M., Drossart, P., Piccioni, G.,
Formisano, V., 2011. Non-LTE CO limb emission at 4.7μm in the upper atmosphere of venus, mars and earth: observations and modeling. Planet Sp. Sci. 59,
1010–1018. doi:10.1016/j.pss.2010.07.023.
Gilli, G., López-Valverde, M.A., Peralta, J., Bougher, S., Brecht, A., Drossart, P., Piccioni, G., 2015. Carbon monoxide and temperature in the upper atmosphere ofvenus from VIRTIS/venus express non-LTE limb measurements. Icarus 248, 478–
498. doi:10.1016/j.icarus.2014.10.047.
González-Galindo, F., Chaufray, J.Y., López-Valverde, M.A., Gilli, G., Forget, F.,
Leblanc, F., Modolo, R., Hess, S., Yagi, M., 2013. Three-dimensional martian ionosphere model: I. the photochemical ionosphere below 180 km. J. Geophys. Res.
(Planets) 118, 2105–2123. doi:10.1002/jgre.20150.
González-Galindo, F., Forget, F., López-Valverde, M.A., Angelats i Coll, M., 2009. A
ground-to-exosphere martian general circulation model: 2. atmosphere duringsolstice conditions. thermospheric polar warming. J. Geophys. Res. (Planets) 114,
8004–+. doi:10.1029/2008JE003277.
González-Galindo, F., López-Valverde, M.A., Angelat i Coll, M., Forget, F., 2005. Extension of a martian general circulation model to thermospheric altitudes: Uv
heating and photochemical models. J. Geophys. Res. (Planets) 110, 9008–+.
Grassi, D., Migliorini, A., Montabone, L., Lebonnois, S., Cardesìn-Moinelo, A., Piccioni, G., Drossart, P., Zasova, L.V., 2010. Thermal structure of venusian nighttime
mesosphere as observed by VIRTIS-venus express. J. Geophys. Res. (Planets) 115,
9007. doi:10.1029/2009JE003553.
Grassi, D., Politi, R., Ignatiev, N.I., Plainaki, C., Lebonnois, S., Wolkenberg, P.,
Montabone, L., Migliorini, A., Piccioni, G., Drossart, P., 2014. The venus nighttime atmosphere as observed by the VIRTIS-m instrument. average ﬁelds fromthe complete infrared data set. J. Geophys. Res. (Planets) 119, 837–849. doi:10.
10 02/2013JE0 04586.
Guilbon, S., Määttänen, A., Montmessin, F., 2015. Modeling venus’ clouds with themoment method: paving the way for 3d GCM simulations. European Planetary
Science Congress 2015, Held 27 September - 2 October, 2015 in Nantes, France,
EPSC2015–614, 10.
Gurwell, M.A., Muhleman, D.O., Shah, K.P., Berge, G.L., Rudy, D.J., Grossman, A.W.,
1995. Observations of the CO bulge on venus and implications for mesosphericwinds. Icarus 115, 141–158. doi:10.1006/icar.1995.1085.
Haus, R., Kappel, D., Arnold, G., 2013. Self-consistent retrieval of temperature proﬁles and cloud structure in the northern hemisphere of venus using VIRTIS/VEX and PMV/VENERA-15 radiation measurements. Planet Sp. Sci. 89, 77–
101. doi:10.1016/j.pss.2013.09.020.
Haus, R., Kappel, D., Arnold, G., 2014. Atmospheric thermal structure and cloud features in the southern hemisphere of venus as retrieved from VIRTIS/VEX radiation measurements. Icarus 232, 232–248. doi:10.1016/j.icarus.2014.01.020.
Hedin, A.E., Niemann, H.B., Kasprzak, W.T., Seiff, A., 1983. Global empiricalmodel of the venus thermosphere. J. Geophys. Res. 88, 73–83. doi:10.1029/
JA088iA01p0 0 073.
Hoshino, N., Fujiwara, H., Takagi, M., Kasaba, Y., 2013. Effects of gravity waves on
the day-night difference of the general circulation in the venusian lower thermosphere. J. Geophys. Res. (Planets) 118, 2004–2015. doi:10.1002/jgre.20154.

72

G. Gilli et al. / Icarus 281 (2017) 55–72

Hourdin, F., Couvreux, F., Menut, L., 2002. Parameterization of the dry convective
boundary layer based on a mass ﬂux representation of thermals.. J. Atmosph.
Sci. 59, 1105–1123. doi:10.1175/1520-0469(2002)0591105:POTDCB2.0.CO;2.
Hourdin, F., Musat, I., Bony, S., Braconnot, P., Codron, F., Dufresne, J.L., Fairhead, L.,
Filiberti, M.A., Friedlingstein, P., Grandpeix, J.Y., Krinner, G., Levan, P., Li, Z.X.,
Lott, F., 2006. The LMDZ4 general circulation model: climate performance andsensitivity to parametrized physics with emphasis on tropical convection. Climate Dyn. 27, 787–813. doi:10.10 07/s0 0382-0 06-0158-0.
Huestis, D.L., Bougher, S.W., Fox, J.L., Galand, M., Johnson, R.E., Moses, J.I., Pickering, J.C., 2008. Cross sections and reaction rates for comparative planetary
aeronomy. Sp. Sci. Rev. 139, 63–105. doi:10.10 07/s11214-0 08-9383-7.
Kasprzak, W.T., Niemann, H.B., Hedin, A.E., Bougher, S.W., Hunten, D.M., 1993. Neutral composition measurements by the pioneer venus neutral mass spectrometer during orbiter re-entry. Geophys. Res. Lett. 20, 2747–2750. doi:10.1029/
93GL02241.
Keating, G.M., Bertaux, J.L., Bougher, S.W., Dickinson, R.E., Cravens, T.E., Hedin, A.E.,
1985. Models of venus neutral upper atmosphere - structure and composition.
Adv. Sp. Res. 5, 117–171. doi:10.1016/0273- 1177(85)90200- 5.
Keating, G.M., Nicholson, J.Y., Lake, L.R., 1980. Venus upper atmosphere structure. J.
Geophys. Res. 85, 7941–7956. doi:10.1029/JA085iA13p07941.
Knollenberg, R., Hunten, D., 1980. The microphysics of the clouds of Venus: resultsof the pioneer venus particle size spectrometer experiment.. J. Geophys. Res. 85,
8039–8058.
Krasnopolsky, V.A., 2012. A photochemical model for the venus atmosphere at 47–
112 km. Icarus 218, 230–246. doi:10.1016/j.icarus.2011.11.012.
Krause, P., Sornig, M., Wischnewski, C., Herrmann, M., Stangier, T., Sonnabend, G.,
Kostiuk, T., Livengood, T., 2015. Long-term variation of temperature and dynamic
at the morning terminator in venus upper atmosphere from ground-based infrared heterodyne spectroscopy. In: European Planetary Science Congress 2015,
held 27 September - 2 October, 2015 in Nantes, France, 10, pp. EPSC2015–544.
Krause, P., Wischnewski, C., Sornig, M., Stangier, T., Sonnabend, G., Herrmann, M.,
Wiegand, M., Kostiuk, T., Livengood, T., 2016. Temperature measurements invenus upper atmosphere between 2007 and 2015 from ground-based infrared heterodyne spectroscopy. In: EGU General Assembly Conference Abstracts, p. 16468.
Lebonnois, S., Covey, C., Grossman, A., Parish, H., Schubert, G., Walterscheid, R., Lauritzen, P., Jablonowski, C., 2012. Angular momentum budget in general circulation models of superrotating atmospheres: a critical diagnostic. J. Geophys. Res.
(Planets) 117, 12004. doi:10.1029/2012JE004223.
Lebonnois, S., Hourdin, F., Eymet, V., Crespin, A., Fournier, R., Forget, F., 2010. Superrotation of venus’ atmosphere analyzed with a full general circulation model. J.
Geophys. Res. (Planets) 115, 6006. doi:10.1029/20 09JE0 03458.
Lebonnois, S., Sugimoto, N., Gilli, G., 2016. Wave analysis in the atmosphere of venus
below 100-km altitude, simulated by the LMD venus GCM. Icarus 278, 38–51.
doi:10.1016/j.icarus.2016.06.004.
Lefèvre, F., Bertaux, J.L., Clancy, R.T., Encrenaz, T., Fast, K., Forget, F., Lebonnois, S.,
Montmessin, F., Perrier, S., 2008. Heterogeneous chemistry in the atmosphereof mars. Nature 454, 971–975. doi:10.1038/nature07116.
Lefèvre, F., Lebonnois, S., Montmessin, F., Forget, F., 2004. Three-dimensional modeling of ozone on mars. J. Geophys. Res. 109, E07004. doi:10.1029/20 04JE0 02268.
Lefevre, M., Spiga, A., Lebonnois, S., 2016. Three-dimensional mesoscale modeling of
the venusian cloud layer and associated gravity waves. Submitted to JGR.
Lellouch, E., Clancy, T., Crisp, D., Kliore, A.J., Titov, D., Bougher, S.W., 1997. Monitoring of mesospheric structure and dynamics. In: Bougher, S.W., Hunten, D.M.,
Phillips, R.J. (Eds.), Venus II: Geology, Geophysics, Atmosphere, and Solar Wind
Environment, p. 295.
Limaye, S.S., Bougher, S., Chamberlain, S., Clancy, R.T., Gilli, G., Grassi, D., Haus, R.,
Herrmann, M., Imamura, T., Kohler, E., Krause, P., Lebonnois, S., Mahieux, A.,
Sandor, A., Sornig, M., Svedhem, H., Tellmann, S., Vandaele, A.C., Widemann, T.,
Wilson, C., Mueller-Wodarg, I., Zasova, L., 2015. Improved knowlegde of venus
atmospheric structure. In: European Planetary Science Congress 2015, held 27
September - 2 October, 2015 in Nantes, France, 10, pp. EPSC2015–758.
Lopez-Puertas, M., Taylor, F.W., 2001. Non-LTE Radiative Transfer in the Atmosphere.
World Scientiﬁc Pub., Singapore, Series on Atmospheric Oceanic and Planetary
Physics, vol. 3 ISBN 9810245661..
Lopez-Valverde, M.A., Lopez-Puertas, M., 1994. A non-local thermodynamic equilibrium radiative transfer model for infrared emission in the atmosphere of mars.
2: daytime populations of vibrational levels. J. Geophys. Res. 99, 13117–13132.
doi:10.1029/94JE01091.
López-Valverde, M.A., López-Puertas, M., Funke, B., Gilli, G., García-Comas, M.,
Drossart, P., Piccioni, G., Formisano, V., 2011. Modelling the atmospheric limb
emission of CO2 at 4.3 μm in the terrestrial planets. Planet. Space Sci. 59, 988–
998. doi:10.1016/j.pss.2010.02.001.
Lott, F., Guez, L., 2013. A stochastic parameterization of the gravity waves due toconvection and its impact on the equatorial stratosphere. J. Geophys. Res. (Atmosph.) 118, 8897–8909. doi:10.1002/jgrd.50705.
Määttänen, A., Guilbon, S., Stolzenbach, A., Bekki, S., Montmessin, F., 2015. VenLA:
the LATMOS venus cloud model. In: European Planetary Science Congress 2015,
held 27 September - 2 October, 2015 in Nantes, France, 10, pp. EPSC2015–634.
Madronich, S., Flocke, S., 1999. The role of solar radiation in atmospheric chemistry.
In: Environmental Photochemistry. Springer, pp. 1–26.
Mahieux, A., Vandaele, A.C., Bougher, S.W., Drummond, R., Robert, S., Wilquet, V., Chamberlain, S., Piccialli, A., Montmessin, F., Tellmann, S., Pätzold, M.,
Häusler, B., Bertaux, J.L., 2015. Update of the venus density and temperature
proﬁles at high altitude measured by SOIR on board venus express. Planet Sp.
Sci. 113, 309–320. doi:10.1016/j.pss.2015.02.002.

Mellor, G.L., Yamada, T., 1982. Development of a turbulent closure model for geophysical ﬂuid problems. Rev. Geophys. Space Phys. 20, 851–875.
Migliorini, A., Grassi, D., Montabone, L., Lebonnois, S., Drossart, P., Piccioni, G., 2012.
Investigation of air temperature on the nightside of venus derived from VIRTISh on board venus-express. Icarus 217, 640–647. doi:10.1016/j.icarus.2011.07.013.
Niemann, H.B., Kasprzak, W.T., Hedin, A.E., Hunten, D.M., Spencer, N.W., 1980. Mass
spectrometric measurements of the neutral gas composition of the thermosphere and exosphere of venus. J. Geophys. Res. 85, 7817–7827. doi:10.1029/
JA085iA13p07817.
Peralta, J., Hueso, R., Sánchez-Lavega, A., Piccioni, G., Lanciano, O., Drossart, P., 2008.
Characterization of mesoscale gravity waves in the upper and lower clouds ofvenus from VEX-VIRTIS images. J. Geophys. Res. (Planets) 113, 0. doi:10.1029/
20 08JE0 03185.
Peralta, J., López-Valverde, M.A., Gilli, G., Piccialli, A., 2016. Dayside temperatures inthe venus upper atmosphere from venus express/VIRTIS nadir measurements at
4.3 μm. A & A 585, A53. doi:10.1051/0 0 04-6361/201527191.
Piccialli, A., Montmessin, F., Belyaev, D., Mahieux, A., Fedorova, A., Marcq, E.,
Bertaux, J.L., Tellmann, S., Vandaele, A.C., Korablev, O., 2015. Thermal structureof venus nightside upper atmosphere measured by stellar occultations with SPICAV/venus express. Planet Sp. Sci. 113, 321–335. doi:10.1016/j.pss.2014.12.009.
Piccialli, A., Titov, D.V., Sanchez-Lavega, A., Peralta, J., Shalygina, O., Markiewicz, W.J.,
Svedhem, H., 2014. High latitude gravity waves at the venus cloud tops as observed by the venus monitoring camera on board venus express. Icarus 227,
94–111. doi:10.1016/j.icarus.2013.09.012.
Roldán, C., López-Valverde, M.A., López-Puertas, M., Edwards, D.P., 20 0 0. Non-LTE
infrared emissions of CO2 in the atmosphere of venus. Icarus 147, 11–25. doi:10.
10 06/icar.20 0 0.6432.
Sander, S.P., Abbatt, J., Barker, J., Burkholder, J., Friedl, R., Golden, D., 2011. Chemical
Kinetics and Photochemical Data for Use in Atmospheric Studies, Evaluation No.
17. Technical Report. JPL.
Schubert, G., Covey, C., del Genio, A., Elson, L.S., Keating, G., Seiff, A., Young, R.E.,
Apt, J., Counselman, C.C., Kliore, A.J., Limaye, S.S., Revercomb, H.E., Sromovsky, L.A., Suomi, V.E., Taylor, F., Woo, R., von Zahn, U., 1980. Structure andcirculation of the venus atmosphere. J. Geophys. Res. 85, 8007–8025. doi:10.
1029/JA085iA13p08007.
Seiff, A., 1983. Thermal structure of the atmosphere of venus 215–279.
Sharma, R.D., 2014. Technical note: a new mechanism of 15 μm emission in
the mesosphere-lower thermosphere (MLT). Atmosph. Chem. Phys. Discuss. 14,25083–25093. doi:10.5194/acpd- 14- 25083- 2014.
Sharma, R.D., Wintersteiner, P.P., 1990. Role of carbon dioxide in cooling planetary
thermospheres. Geophys. Res. Lett. 17, 2201–2204.
Sonnabend, G., Kroetz, P., Sornig, M., Stupar, D., 2010. Direct observations of venus
upper mesospheric temperatures from ground based spectroscopy of CO2 . Geophys Res. Lett. 37, L11102. doi:10.1029/2010GL043335.
Soret, L., Gérard, J.C., Montmessin, F., Piccioni, G., Drossart, P., Bertaux, J.L., 2012.
Atomic oxygen on the venus nightside: global distribution deduced from airglow mapping. Icarus 217, 849–855. doi:10.1016/j.icarus.2011.03.034.
Soret, L., Gérard, J.C., Piccioni, G., Drossart, P., 2014. Time variations of o2 (a1 )nightglow spots on the venus nightside and dynamics of the upper mesosphere.
Icarus 237, 306–314. doi:10.1016/j.icarus.2014.03.034.
Stolzenbach, A., Lefèvre, F., Lebonnois, S., Määttänen, A., Bekki, S., 2014. Three-dimensional modelling of venus photochemistry. In: EGU General Assembly Conference Abstracts, p. 5315.
Stolzenbach, A., Lefèvre, F., Lebonnois, S., Maattanen, A.E., Bekki, S., 2015. Three-dimensional modelling of venus photochemistry. AGU Fall Meeting Abs..
Svedhem, H., Müller-Wodarg, I., 2014. Results of the venus express aerobraking campaign. In: AAS/Division for Planetary Sciences Meeting Abstracts, p. 302.09.
Svedhem, H., Wilson, C., Piccioni, G., 2015. Introduction to the special issue on venus
exploration. Planet Sp. Sci. 113, 1. doi:10.1016/j.pss.2015.06.003.
Tellmann, S., Häusler, B., Hinson, D.P., Tyler, G.L., Andert, T.P., Bird, M.K., Imamura, T.,
Pätzold, M., Remus, S., 2012. Small-scale temperature ﬂuctuations seen by thevera radio science experiment on venus express. Icarus 221, 471–480. doi:10.
1016/j.icarus.2012.08.023.
Tellmann, S., Pätzold, M., Häusler, B., Bird, M.K., Tyler, G.L., 2009. Structure of
the venus neutral atmosphere as observed by the radio science experimentvera on venus express. J. Geophys. Res. (Planets) 114, E00B36. doi:10.1029/
20 08JE0 03204.
Vandaele, A.C., Chamberlain, S., Mahieux, A., Ristic, B., Robert, S., Thomas, I.,
Trompet, L., Wilquet, V., Belyaev, D., Fedorova, A., Korablev, O., Bertaux, J.L.,
2016. Contribution from SOIR/VEX to the updated venus international reference
atmosphere (VIRA). Adv. Sp. Res. 57, 443–458. doi:10.1016/j.asr.2015.08.012.
Wischnewski, C., Sornig, M., Sonnabend, G., Stangier, T., Krause, P., Herrmann, M.,
Kostiuk, T., Livengood, T., 2015. Thermal structure of venus’ dayside in 110 km
altitude based on ground-based heterodyne observations between 2007 and2014. In: European Planetary Science Congress 2015, held 27 September - 2 October, 2015 in Nantes, France, 10, pp. EPSC2015–546.
Woo, R., Armstrong, J.W., Ishimaru, A., 1980. Radio occultation measurements of turbulence in the venus atmosphere by pioneer venus. J. Geophys. Res. 85, 8031–
8038. doi:10.1029/JA085iA13p08031.
Zalucha, A.M., Brecht, A.S., Rafkin, S., Bougher, S.W., Alexander, M.J., 2013. Incorporation of a gravity wave momentum deposition parameterization into the venus
thermosphere general circulation model (VTGCM). J. Geophys. Res. (Planets) 118,
147–160. doi:10.1029/2012JE004168.
Zhang, S., Bougher, S.W., Alexander, M.J., 1996. The impact of gravity waves on thevenus thermosphere and o2 IR nightglow. J. Geophys. Res. 101, 23195–23206.
doi:10.1029/96JE02035.