Articles | Volume 10, issue 3
Research article
09 Jul 2019
Research article |  | 09 Jul 2019

A radiative-convective model based on constrained maximum entropy production

Vincent Labarre, Didier Paillard, and Bérengère Dubrulle

The representation of atmospheric convection induced by radiative forcing is a long-standing question mainly because turbulence plays a key role in the transport of energy as sensible heat, geopotential, and latent heat. Recent works have tried using the maximum entropy production (MEP) conjecture as a closure hypothesis in 1-D simple climate models to compute implicitly temperatures and the vertical energy flux. However, these models fail to reproduce realistic profiles. To solve the problem, we describe the energy fluxes as a product of a positive mass mixing coefficient with the corresponding energy gradient. This appears as a constraint which imposes the direction and/or limits the amplitude of the energy fluxes. It leads to a different MEP steady state which naturally depends on the considered energy terms in the model. Accounting for this additional constraint improves the results. Temperature and energy flux are closer to observations, and we reproduce stratification when we consider the geopotential. Variations in the atmospheric composition, such as a doubling of the carbon dioxide concentration, are also investigated.

1 Introduction

The climate system is complex and usually divided into different components: atmosphere, ocean, cryosphere, lithosphere, and biosphere (Peixoto and Oort1992). There are different approaches to climate modelling (Randall et al.2007). We can classify them in three main classes. Global climate models (GCMs) are the more sophisticated ones (see Dufresne et al.2013, for an example). They explicitly represent the circulation of the atmosphere and ocean. Earth models of intermediate complexity (EMICs) simulate the Earth system with more simplifications than GCMs (see Goosse et al.2010, for an example). These simplifications allow simulations over larger time periods, which is useful to study past climates. Simple climate models (SCMs) use only a few key processes to answer specific questions (see Paillard1998, for an example).

Both the complex and simple models have different strengths and weaknesses and are used for different applications. For example, GCMs are largely used to make climate projections for the next century. Since the numerical resolution of dynamical equations from the micro-scale (of the order 10-3 m for viscous dissipation) to the scale of interest (≃107 m for the typical size of the Earth) is still impossible, GCMs need to represent sub-grid processes such as small eddies, convection, or cloud formation. To do it, models usually express the intensity of fluxes due to unresolved phenomena as a function of the resolved variables. This approach is called a “turbulent closure” and usually requires the introduction of empirical parameters such as turbulent master length scale, turbulent velocity diffusion terms (Mellor and Yamada1982), or the use of different quantities like convective available potential energy (CAPE) or convective inhibition (CIN) to fix the convective intensity (Yano et al.2013). These parameterizations change from one model to another, resulting in different predictions (Stevens and Bony2013). They also require adjusting the numerous free parameters (“tuning”) in order to track observations (Hourdin et al.2017). SCMs appear as an interesting alternative when considering past or future climate over a larger period like glacial–interglacial cycles. Indeed, an SCM is a set of a reasonable number of equations and physical quantities, providing an easier assessment of the impact of the parameters of the models (such as the concentration of greenhouse gases for atmospheric models). Many SCMs are based on the idea that the computation of all the microscopic details may be unnecessary if we are interested in quantities at larger spatio-temporal scales.

Consequently, a lot of SCMs describe the Earth simply with energetic considerations (North et al.1981). Those models are called energy balance models (EBMs). The atmosphere is mainly driven by radiative forcing: solar radiation imparts energy to the Earth, which emits infra-red radiation to space. This heating is not homogeneous around the globe for different reasons such as geometry or variations in albedo with the nature of the ground. The insolation is more important for the tropics than the poles and leads to latitudinal heat transport. For the vertical axis, the amount of radiative energy absorbed by the Earth naturally depends on the atmospheric components. The ground usually receives more solar radiation because of the relative transparency of the atmosphere. As a result, the atmosphere is heated from below, which may lead to unstable situations where the temperature gradient exceeds the adiabatic gradient. Then, this causes atmospheric motion and vertical heat transport named convection. Since EBMs are usually based only on the energy budget, it is necessary to specify a relation between the energy fluxes and the temperature gradient, called a closure hypothesis. So EBMs mainly differ by the representation of the fluxes between the fluid layers of the Earth (atmosphere and ocean). For example, horizontal fluxes are sometimes represented by purely diffusion terms (North et al.1981). On the other hand, vertical energy transport has been modelled using different approaches such as the convective adjustment. The latter consists of computing the temperature profile at radiative equilibrium for stable regions and adjusting it where the critical gradient is exceeded (Manabe and Strickler1964). Representing both the horizontal and the vertical energy fluxes is an important issue for EBMs. This concerns the subject of the present paper.

Since the 70s (Paltridge1975), maximum entropy production (MEP; Martyushev and Seleznev2006) is also used as a closure hypothesis in EBMs. This conjecture stipulates that the climatic system (or one of its component) optimizes its entropy production due to internal heat transfers. It allows for computing implicitly (i.e. without the computation of the dynamics) horizontal fluxes (O'brien and Stephens1995; Lorenz et al.2001) and vertical fluxes (Ozawa and Ohmura1997; Pujol and Fort2002) without the parameterizations required in more conventional models. Former MEP-based models (MEPMs) have been criticized for three main reasons. One is the absence of dynamics and the validity of MEP (Rodgers1976). The second criticism deals with the extra parameterizations or the assumptions used in MEPMs. Indeed, one may ask oneself if the successes of the models are really due to the MEP hypothesis or to tuning or other ingredients (Goody2007). The final criticism concerns the usually simplified description of the radiative forcing in these models. Recently, a MEPM overcoming the last two criticisms has been built in Herbert (2012). It includes a refined description of the radiative budget in the net exchange formalism, without extra assumptions. The only adjustable quantities concern the radiative budget, such as the albedo, and not the atmospheric or oceanic energy transport. The model provides a relatively good approximation for the temperature and horizontal heat fluxes (Herbert et al.2011b).

However, the vertical energy fluxes are still overestimated in such models when compared to observations or conventional radiative-convective models (RCMs) like Manabe and Strickler (1964). Furthermore, the energy fluxes are not always oriented against the energy gradient and it does not predict stratification in the upper atmosphere. This is not surprising because geopotential was not taken into account. Yet, we know from fluid mechanics that gravity plays a major role in natural convection (Rieutord2015). Gravity is also obviously responsible for stratification in the upper atmosphere. In this paper, we develop a MEPM that describes more properly the atmospheric convection. In the same spirit as the previous SCMs, we do not attempt to resolve the dynamical equations, but we add only some key features. Two ingredients are introduced to represent vertical heat fluxes more correctly. The first one is to describe energy transport as the product of a non-homogeneous mixing mass coefficient and the specific energy gradient. This brings a new constraint into a MEPM. The second one is to consider different energy terms: sensible heat, geopotential, and latent heat. We show that this simplified description of the energy transport, combined with the MEP closure hypothesis, can lead to relatively realistic results.

The outline of this paper is as follows. In the first part we describe our model, presenting the transport of heat by mixing. The formulation of the constrained MEP optimization problem is given (Sect. 2). Then, we compute the temperature, the specific energy, and the energy flux profiles. We give a physical interpretation of the effect of the constraint emerging from the positivity of the mass mixing coefficient. The impact of different expressions for energy is discussed (Sect. 3). A sensitivity test for the concentrations of O3 and CO2 is also performed. Finally, we discuss further works and objectives (Sect. 4). The computation of the geopotential is given in Appendix A and the resolution of the optimization problem is described in Appendix B.

2 Model

2.1 Vertical structure of the atmosphere

The atmosphere is divided into a column of N vertical layers. We work with prescribed pressure levels, so the elevation z depends on the temperature profile (see Appendix A). The CO2, O3, and water vapour profiles are fixed according to observations by McClatchey et al. (1972) and the ground is represented by a layer with a fixed surface albedo α. The atmosphere is supposed to be in hydrostatic equilibrium and is considered an ideal gas. The specific energy (energy per unit mass) in layer i – of mean elevation zi, temperature Ti, and mixing ratio qi (ratio between the mass of water vapour and total mass of the air for a given volume) – is the so-called moist static energy,

(1) e i = C p T i + g z i + L q i ,

where Cp=1005 J kg−1 K−1 is the heat capacity of the air, g=9.81 m s−2 is the terrestrial acceleration due to gravity, and L=2.5×106 J kg−1 is the latent heat of vaporization.

We note i, the net radiative energy input in layer i, by taking into account several effects: shortwave radiation, longwave radiation, reflection, and reabsorption. More explicitly,


where SWi is the downward radiative energy flux for short waves, SWi is the upward radiative energy flux for short waves, LWi is the downward radiative energy flux for long waves, and LWi is the upward radiative energy flux for long waves.

We use the code developed in Herbert et al. (2013) to compute the radiative budget. This model was developed to give a realistic description of the absorption properties of the more radiatively active constituents of the atmosphere while keeping a smooth dependence of the radiative flux with respect to the temperature profile. As suggested by the authors, this last requirement is important in the framework of a variational problem. The model is based on net exchange formalism (Dufresne et al.2005), where the basic variables are the net exchange rates between each pair of layers instead of radiative fluxes.

In the longwave domain, the code decomposes the spectrum into 22 narrow bands and, in each band, it accounts for absorption by water vapour and carbon dioxide only. The absorption coefficient is computed using the statistical model of Goody (1952) with the data from Rodgers and Walshaw (1966). For the spatial integration, the diffusive approximation is performed with the standard diffusion factor μ=1/1.66. Apart from the absorption data, given once and for all, the inputs of the model are the water vapour density, temperature profile, and carbon dioxide concentration. One may fix either absolute or relative humidity. In the shortwave domain, absorption by water vapour and ozone is accounted for by adapting the parameterization from Lacis and Hansen (1974). The input parameters for the model are the water vapour density and ozone density profiles, as well as surface albedo and solar constant. Clouds are not considered in the model. More details can be found in Herbert et al. (2013) and its supplementary material. The net radiative budget for the atmospheric layer i, i, is given by summing over all terms involving the layer in question. In particular, i is a function of all temperatures Tjj=0,,N in the profile. Then,

(2) R i T , q , O 3 , CO 2 , α = SW i q , O 3 , α + LW i T , q , CO 2 .

In the previous equation, and in the following, T, q, O3, and CO2 will refer to complete profiles (i.e. T={Ti}i=0,,N). Given that q (or h=q/qs(T) fixed relative humidity), O3, CO2, and α are fixed in our model, we will only indicate the T dependence.

The vertical energy flux is represented by mixing between adjacent layers. Then, the net upward energy flux between layers i and i−1 writes as

(3) F i = m i e i - 1 - e i ,

where mi is a mixing coefficient that represents a mass per unit time and per unit surface. It is not (necessarily) homogeneous in all the column. We notice that mi is typically the kind of coefficient that requires, at some point, a parameterization in usual climate models. Here m is obtained from the MEP procedure.

Taking into account the net radiative energy budget i, the energy balance at the stationary state for the layer i reads as

(4) F i - F i + 1 + R i = 0 .

Figure 1Discretization of an atmospheric column into N layers. The layer i – at temperature Ti, fixed pressure Pi, elevation zi (which depends on the temperature profile), and mixing ratio qi – has a specific energy ei=CpTi+gzi+Lqi; mi is the mass mixing coefficient between layers i−1 and i which leads to the net energy flux Fi=mi(ei-1-ei); i is the net radiative energy budget in the layer i. The ground is represented by layer 0 with fixed surface albedo α.


2.2 Maximum entropy production with constraint

The principle of our model is to determine the fluxes F and temperatures T with the maximization of the entropy production. In the thermodynamics of diffusive processes, the entropy production is expressed as the sum of products of the fluxes with their associated thermodynamic forces. But we aim here at representing convection, not diffusion. We need to represent how the air parcels are mixed but only after (turbulent) convective motions. We can usefully separate this process into two steps. First, we assume a pseudo-adiabatic motion of an air parcel from one layer to another due to convection. This step is purely mechanical, without entropy production since we neglect viscous dissipation. The energy of the air parcel is conserved but its temperature and composition may change during the motion. If the air parcel was initially in layer i and goes to the layer i+1, the temperature of the air parcel becomes, by conservation of energy,

(5) T i = T i + 1 C p g ( z i - z i + 1 ) + L ( q i - q i + 1 ) .

Here, we have assumed that the water vapour concentration changes pseudo-adiabatically during the convection and not due to the mixing. This is not fully consistent since we do not impose water conservation. Secondly, the air parcel is mixed by diffusion with the ambient air in layer i+1 and transfers an amount of sensible heat per unit mass CpTi. At the same time, air parcels leave the layer i+1 for the layer i with a sensible heat per unit mass CpTi+1. So the net flux of sensible heat due to this process is

(6) m i + 1 C p ( T i - T i + 1 ) = m i + 1 ( e i - e i + 1 ) = F i + 1 ,

where mi+1 is the mass mixing coefficient between layers i and i+1. So the entropy production that results due to the sensible heat exchange between layers i and i+1 is the product of the flux Fi+1 with the thermodynamic force associated with the sensible heat only (i.e. the gradient of inverse temperature)

(7) σ i + 1 = F i + 1 1 T i + 1 - 1 T i .

By summing over all layers, and using the fact that FN+1=0, we show that the total entropy production can be written as

(8) σ = i = 0 N F i - F i + 1 T i .

In thermodynamics, more terms may contribute to entropy production such as pressure-volume work, or mixing. As for other MEPMs (Kleidon2010), we only retain the sensible heat exchange term. The geopotential and latent heat terms appear in the entropy production only as a result of our representation of convective transport, which is supposed to occur as a mechanically induced mass transport without entropy production. We can easily express the entropy production with temperatures σ(T) by using the energy balance in stationary state, Fi-Fi+1+Ri(T)=0. Then the problem is usually solved in terms of temperature with a global constraint of energy conservation:

(9) max T 0 , , T N - i = 0 N R i ( T ) T i i = 0 N R i ( T ) = 0 .

Given the form of the energy transport (Eq. 3), we here need to have additional constraints. Namely, the mass mixing coefficients mi must be positive. It is then natural to solve the problem of flux by expressing the entropy production σ(F) and inequality constraints mi≥0 with energy fluxes. Assuming that the relation ℛ(T) is invertible (Appendix B), we can formally write

(10) F i + 1 - F i = R i ( T ) T i = R i - 1 ( F ) .

This results in the following optimization problem with inequality constraints:

(11) max F 1 , , F N i = 0 N F i - F i + 1 R i - 1 ( F ) m i 0 with F i = - m i e i - e i - 1 .

The constraint Fi=-miei-ei-1 with mi≥0 naturally depends on the specific energy e used in the model. The later simply imposes the energy fluxes to be opposed to the energy gradient. We point out that the energy conservation is implicit in the flux formulation of the variational problem and does not need to be imposed as a constraint here.

3 Results

We have computed temperature, specific energy (energy per unit mass), and energy flux profiles for different prescribed atmospheric compositions from McClatchey et al. (1972) corresponding to the tropical, mid-latitude summer, mid-latitude winter, sub-arctic summer, and sub-arctic winter conditions. We work with a fixed relative humidity profile (ratio of the partial pressure of water vapour to the equilibrium vapour pressure of water at a given temperature). Typical values of surface albedo are used: α=0.1 for the tropical and mid-latitude conditions and α=0.6 for sub-arctic ones. The atmosphere is discretized in N=20 vertical levels.

3.1 The effect of the constraint

We investigate the effect of the following energy terms on the constraint:

  • sensible heat CpT;

  • geopotential gz;

  • latent heat for a water-vapour-saturated air Lqs(T), where qs is the mixing ratio at the saturation point. Since we work in pressure coordinates, it depends only on local temperature. This is a first attempt to take into account the effect of humidity without an explicit derivation of the humidity profile and water cycle. However, the radiative budget is still computed using a fixed standard relative humidity profile.

For illustration purposes, we can consider the case with only two layers. We note F=m(e1-e2) the net energy flux from layer 1 to layer 2, where m is the mass mixing coefficient between layers (cf. Fig. 2). In this simple case, the entropy production is written as σ=F1/T2-1/T1 and is limited by the constraint m0F(e1-e2)0. We interpret the different energy terms on this constraint as follow.

Figure 2Energy exchanges between two layers of elevation z1 and z2 (z2z1), temperatures T1 and T2, and specific energy e1 and e2. We note F=m(e1-e2), the energy flux from layer 1 to layer 2, where m is the mass mixing coefficient between the two layers.


  • e=CpT: F≥0 if T1T2. The constraint simply imposes the energy transport from hot to cold regions.

  • e=CpT+gz: F≥0 if T1T2+g(z2-z1)/Cp. The geopotential gz limits the upward energy flux. We predict a warmer air at the bottom and a colder air at the top compared to the model with only sensible heat. Mechanically, we can see this as the expression of the fact that an air parcel from layer 1 may not have enough energy to move adiabatically to layer 2.

  • e=CpT+gz+Lqs: F≥0 if T1T2+g(z2-z1)+L(qs(T2)-qs(T1))/Cp. Since qs(T) is an increasing function, qs(T2)−qs(T1) has the same sign as T2T1. The temperature gradient is usually negative (i.e. T2T1). Adding the latent heat at saturation makes the upward transport of energy less constrained. Consequently, the atmospheric temperature gradient weakens. Mechanically, we can see this as the fact that moist convection is easier than dry convection because of the transport of latent energy from the bottom to the top of the atmosphere.

Figure 3Energy flux F, specific energy e, and temperature T for tropical atmospheric composition measured by McClatchey et al. (1972) and different expressions for energy on the constraint (e=CpT, e=CpT+gz, and e=CpT+gz+Lqs). The elevation is given in pressure level P. Results for the unconstrained model of Herbert et al. (2013) are represented. We also give the temperature profile corresponding to the measurements of McClatchey et al. (1972), labelled as “reference” for qualitative comparison.


3.2 General remarks

Various profiles are shown for tropical (Fig. 3) and sub-arctic winter (Fig. 4) conditions. The outputs of our constrained model are labelled by the energy terms taken into account in the constraint (e=CpT, e=CpT+gz, or e=CpT+gz+Lqs). The represented energy profiles are the energy corresponding to the constraint. For CpT, we represent the profile e=CpT, for CpT+gz we represent e=CpT+gz, and for CpT+gz+Lqs we represent e=CpT+gz+Lqs. For e=CpT, the specific energy is trivially more important for hot regions. For e=CpT+gz, the geopotential adds energy to upper layers. For e=CpT+gz+Lqs, the latent heat term adds energy to more humid layers. The results of the “unconstrained” model of Herbert et al. (2013) are represented with the associated thermal energy e=CpT (rigorously speaking, the model is constrained by the global conservation of energy, but we will refer to it as unconstrained since no constraint is imposed on fluxes). Temperature profiles measured by McClatchey et al. (1972) are also represented for qualitative comparison and labelled as “reference”.

Figure 4Energy flux F, specific energy e, and temperature T for sub-arctic winter atmospheric composition measured by McClatchey et al. (1972) and different expressions for energy on the constraint (e=CpT, e=CpT+gz, and e=CpT+gz+Lqs). The elevation is given in pressure level P. Results for the unconstrained model of Herbert et al. (2013) are represented. We also give the temperature profile corresponding to the measurements of McClatchey et al. (1972), labelled as “reference” for qualitative comparison.


For the unconstrained model, the energy flux is positive (i.e. upward) for the tropical (Fig. 3) and sub-arctic winter (Fig. 4) atmospheric compositions in all the column, despite the energy gradient inversion in the upper layers of the atmosphere. Therefore, the flux is in the same direction as the energy gradient in this region. This also corresponds to local negative entropy production and is not physically relevant. Consequently, the upward flux is overestimated, and the temperature gradient is weak. As discussed above (Sect. 2.1), the addition of the constraint mi≥0 imposes the energy flux to be opposed to the specific energy gradient everywhere. As a consequence, the local entropy production is also constrained to be positive as in Ozawa and Ohmura (1997). But with an explicit account for the mass transport as explained above, we can account not only for sensible heat but more generally for moist static energy transfers in a convective column. If an energy flux in this direction is not favourable in matters of entropy production, it vanishes and we have stratification. When the geopotential term is considered (i.e. e=CpT+gz or e=CpT+gz+Lqs), we observe the following.

  • An energy profile divided into three regions:

  1. An unstable surface layer with a decreasing energy profile;

  2. A neutral (slightly stable) mixed layer in the middle atmosphere with a vanishing energy gradient;

  3. An inversion layer at the top of the atmosphere where the energy is increasing with elevation.

  • A vanishing energy flux (stratification) in the upper part of the atmosphere (around P≃300 hPa for the tropics, Fig. 3, and P≃700 hPa for sub-arctic winter, Fig. 4).

We note that the thermal gradient is divided roughly by a factor of 2 when geopotential is considered, which gives a more realistic temperature profile.

3.3 Comparison between profiles

3.3.1 Model outputs for different climatic conditions

The constrained model is obviously sensitive to the water content. Considering e=CpT+gz or e=CpT+gz+Lqs in the constraint gives approximately the same results for sub-arctic winter conditions (Fig. 4) (since qs(T) is weak for low temperatures) while predictions differ for tropical conditions (Fig. 3) (where T and then qs(T) are more important). We have verified that the influence of the surface albedo explains a large part of the temperature modification when we compare different climatic conditions.

3.3.2 Model output vs. reference profile

Before discussing the differences between the outputs and observations, we point out that this conceptual model does not take into account some important physical processes:

  • Insolation is assumed to be constant, fixed at 1368∕4 W m−2 for all conditions. But in reality, it varies with respect to seasons, diurnal cycle, and latitude due to the Earth's geometry and obliquity. So, the model does not take into account the variation in the radiative budget because of these geometrical factors;

  • Horizontal energy fluxes are not considered in this 1-D description;

  • The effect of clouds, which plays an important role in the absorption and emission of radiation (Dufresne and Bony2008), is not implemented in the radiative code.

Therefore, the aim of this study is not to give realistic values of temperature profiles nor vertical energy fluxes but to give a qualitative evaluation of the model. However, we can make some remarks on the comparison of our results to the reference temperature profiles. We observe that our model with e=CpT+gz+Lqs provides better results for tropical conditions (Fig. 3), whereas the computed profiles are not so good for sub-arctic winter conditions (Fig. 4). Considering the previous remarks, one can explain the gap between our model and observations as follows. Constant insolation at the seasonal timescale is valid for the tropics, but it varies strongly for high latitudes. So our model is not adapted to represent a specific season at a high latitude like sub-arctic winter conditions.

Tropical regions are subject to strong vertical motion due to radiative heating. So horizontal energy fluxes are less important and the 1-D vertical description may be more adapted to this case. In contrast, radiative heating is less important for the Arctic (especially in winter), so convection is weaker. Then, the representation of horizontal energy fluxes is essential at high latitudes since they play a major role in heat transport from hot equatorial regions to cold poles. This probably explains why we underestimate temperature for high latitudes (Fig. 4).

3.4 Sensitivity to atmospheric composition

3.4.1 Ozone

When the influence of O3 is not considered in the radiative budget, we observe small changes for specific energy and large changes for temperature in the stratosphere (cf. Fig. 5). Indeed, O3 absorbs solar radiation at the top of the atmosphere, which induces heating of this region. It follows that the temperature in the high atmosphere is more important with ozone, and we even observe an inversion of the temperature gradient. It follows that less solar radiation heats the ground, resulting in smaller surface temperature. For our constrained model with e=CpT+gz+Lqs including ozone, we observe a downward convective energy flux at the top of the atmosphere (Fig. 5). Ozone is therefore associated with heating from the top with an inversion of the temperature gradient and downward energy fluxes in the high atmosphere. This last effect only appears when both geopotential energy and O3 are taken into account. When geopotential is not considered, the upward energy flux is so overestimated that the effect is undetectable (Figs. 3 and 4).

Figure 5Energy flux F, specific energy e, and temperature T, with and without ozone, for the constrained model with e=CpT+gz+Lqs and for the unconstrained model and for tropical atmospheric composition by McClatchey et al. (1972). Oc: with O3, constrained. Nc: without O3, constrained. Ou: with O3, unconstrained. Nu: without O3, unconstrained.


3.4.2 Carbon dioxide

We also have performed the classic experiment of doubling CO2 concentration (Randall et al.2007). The climate sensitivity, defined as the surface temperature differences between computations with [CO2]=560 ppm and [CO2]=280 ppm, is reported for the different atmospheric compositions in Table 1. Conventional models usually represent various processes like water vapour, ice-albedo, lapse rate, and cloud feedbacks. They play an important role in amplifying the climate sensitivity (Forster and Gregory2006). When comparing our values with the literature, we must keep in mind that our model does not represent all those feedbacks. The lapse rate feedback is taken into account. Water vapour feedback is partially represented in a crude way by fixing relative humidity (changes in temperature have an impact on water content and change the radiative budget) but there is no explicit representation of the hydrological cycle. It is technically possible to include the ice-albedo feedback in a MEPM (Herbert et al.2011a), but this is not the case here. Clouds are not represented in the model so we will focus on the comparison between the constrained and unconstrained models using the same radiative scheme. Nevertheless, typical values of climate sensitivity for multi-model averages with only relevant feedbacks are given (Dufresne and Bony2008) for qualitative comparison. For fixed absolute water profile, our values are compared to Dufresne and Bony (2008) accounting only for the lapse rate feedback; for the fixed relative humidity, we compare them to Dufresne and Bony (2008) accounting for both lapse rate and water vapour feedbacks.

The sensitivity values computed here for the unconstrained model differ from Herbert et al. (2013). We have checked that it is only due to the fact that we use N=20 atmospheric layers here instead of N=9 in Herbert et al. (2013) (see Appendix B for the convergence of the algorithm with N). The climate sensitivity is higher for the constrained model than the unconstrained one (despite one exception with fixed absolute moisture for mid-latitude winter). This result may be interpreted as follow. If we start with a radiative forcing induced by a CO2 doubling, it is the same for the two models since we fix identical atmospheric compositions and surface albedos. However, upward energy fluxes are limited for the constrained case which induces more important warming of the lower part of the atmosphere. The induced ground temperature increase is, therefore, more important for the constrained model. This effect is observable when we look at the perturbation of energy flux, specific energy, and temperature (Fig. 6). The constrained model provides more realistic values of sensitivity for fixed relative humidity, particularly for the tropics. Indeed, the sensitivity of 1.60K computed in this case is closer to the literature value of 2.1±0.2K (Dufresne and Bony2008) compared to the unconstrained model.

(Dufresne and Bony2008)(McClatchey et al.1972)

Table 1Climate sensitivity (warming, in K, of the surface due to a doubling of CO2 concentration) of the constrained model with e=CpT+gz+Lqs, unconstrained model of Herbert et al. (2013), and literature (Dufresne and Bony2008). We give the values for different atmospheric compositions and for fixed absolute or relative water vapour profiles.

Download Print Version | Download XLSX

Figure 6Differences in convective energy flux F, specific energy e, and temperature T between [CO2]=560 ppm and [CO2]=280 ppm for tropical atmospheric composition by McClatchey et al. (1972) represented for constrained model with e=CpT, e=CpT+gz, and e=CpT+gz+Lqs, and for the unconstrained model by Herbert et al. (2013).


4 Discussion

MEPMs are different from the usual GCMs or EMICs. Generally, atmospheric models are based on the following.

  1. Kinematics: equations describing how the fluid moves.

  2. Dynamics: equations describing why the fluid moves. They are based on Navier–Stokes equations linking the fluid acceleration to the forces.

  3. Thermodynamics: energy budget equation involving dissipation, radiation, phases changes, etc.

  4. An equation of state such as the ideal gas relation and approximations such as hydrostatic equilibrium to simplify the problem.

  5. Closure hypothesis and/or parameterizations to represent sub-grid processes.

In usual climate models, energy transport is obtained after the computation of the velocity and other fields, whereas in MEPMs energy fluxes are computed implicitly without consideration of the dynamics. According to Dewar (2003) and Dewar (2009), MEP might be viewed as a statistical inference (such as the information theory interpretation of statistical physics by Jaynes1957) rather than a physical law. From this point of view, MEP allows drawing predictions (heat fluxes, specific energy, and temperature profiles) from the partial knowledge of the radiative budget and the energy content, while the effect of the sub-grid-scale processes is unknown. There are recent attempts to link MEP to other variational principles for dynamical systems or non-equilibrium statistical physics like the maximum of Kolmogorov–Sinaï entropy (Mihelich2015). Similar methods like the maximization of dynamical entropy or “Maximum Caliber” (Monthus2011; Dixit et al.2018) may also be relevant to understand where MEP arises from a more formal point of view.

The gap between the model and observations is easily understood and explains why our 1-D description with constant insolation is more adapted for the tropics than arctic conditions. Further improvements are needed to solve these problems. Firstly, a more general 2- or 3-dimensional mass scheme transport is required. Secondly, it is formally possible to compute winds (Karkar and Paillard2015), moisture profiles, and humidity fluxes using MEP. Finally, we need to include a time dependence in our model for seasonal or diurnal cycles. A long-term objective might be to construct an SCM, with a limited number of adjustable parameters.

5 Conclusions

We have investigated the possibility of computing the vertical energy fluxes and temperature in the atmosphere using the MEP closure hypothesis into a simple climate model. The fluxes are then computed in an implicit way, which avoids tuning parameters. Contrary to some authors (Ozawa and Ohmura1997; Pujol and Fort2002; Herbert et al.2013; Pascale et al.2012), we have given a description of how the energy is transported. This paper provides the first attempt, to the best of the authors' knowledge, to introduce such a representation in a MEPM. Different energy terms can be considered: sensible heat, geopotential, and latent heat for a saturated profile. We have shown that this better energetic description of convection allows for obtaining more physically relevant temperature, specific energy, and energy flux profiles, still without any adjustable parameter for the dynamics. In particular, considering geopotential leads to stratification in the upper atmosphere and allows us to reproduce a temperature gradient closer to the observed one. We have investigated the sensitivity of the model when the atmospheric composition is modified. The results were compared to previous MEPMs and the literature. Our model is more sensitive to CO2 than in Herbert et al. (2013) because the geopotential limits upward energy fluxes.

We hope that the present model may be helpful to construct SCMs with a reduced number of adjustable parameters.

Code availability

A Python code, based on the module scipy.optimize, that reproduce the results presented in this paper can be found at (Labarre et al.2019).

Appendix A: Computing geopotential

We show here how to compute the geopotential and the dry static energy eid=CpTi+gzi as a function of temperature. One first writes

(A1) z i = z i - z i - 1 2 + j = 1 i - 1 Δ z j ,

where Δzj=zj+12-zj-12 is the height of the layer j. So if ρ is the density of the air, R is the specific air constant and we assume the atmosphere is an ideal gas at hydrostatic equilibrium

(A2) g Δ z j = z j - 1 2 z j + 1 2 g d z = - p j - 1 2 p j + 1 2 d p ρ = - p j - 1 2 p j + 1 2 R T d p p .

Then, we can compute the mean elevation of a layer with two possible prescriptions:

  • Isothermal layers (T=Tj in the integrand).

    (A3) g Δ z j = R T j ln p j - 1 2 p j + 1 2 .

    So the geopotential reads

    (A4) g z i = R T i ln p i - 1 2 p i + j = 1 i - 1 T j ln p j - 1 2 p j + 1 2 .
  • Dry isentropic layers (T=TjppjRCp in the integrand):

    (A5) g Δ z j = C p T j p j - 1 2 p j R C p - p j + 1 2 p j R C p .

    So the geopotential reads

    (A6) g z i = C p T i p i - 1 2 p i R C p - 1 + j = 1 i - 1 T j p j - 1 2 p j R C p - p j + 1 2 p j R C p .

In both cases, for imposed pressure levels, we obtain the following expression of the specific energy,

(A7) e i d j = 0 N C p δ i j + G i j T j j = 0 N E i j d T j ,

where δij is the Kronecker symbol and Gij and Eijd are the coefficients of constant matrices G and Ed.

Appendix B: Resolution

In order to solve the optimization problem (Eq. 11), we express it in Lagrangian formalism, assuming strong duality holds (Boyd and Vandenberghe2004). We therefore search the critical points of the Lagrangian associated with this problem

(B1) L = σ - i = 1 N μ i m i with m i 0 , μ i 0 and μ i m i = 0 for i = 1 , , N ,

where μ1,,μN are Lagrange multipliers associated with the constraint (mass flux positivity). In order to formulate the problem in terms of energy fluxes F, we must express the inverse temperature X=1/T and the mass mixing m with F.

We use an iterative method to solve this non-linear optimization problem:

  1. We linearize the radiative budget and specific energy around a given temperature profile.

  2. The entropy production and constraints are then quadratic forms of energy fluxes that can be solved numerically.

  3. We reiterate step 1 by linearizing around the temperature profile obtained in step 2 until convergence.

This is a rather standard procedure for optimization though there is no guarantee of finding the global solution in case of multiple local maxima. To overcome this issue, we start with various random initial temperature profiles that may lead to different local maxima. In the end, we retain only the best maximum, which is assumed to be the maximum of entropy production.

B1 Flux–temperature relation

The energy balance equation in stationary state can be written as follows:

(B2) R i ( X ) + F i - F i + 1 = 0 .

At each iteration, we linearize the radiative budget i(X) around a reference temperature profile X0:

(B3) R i ( X ) R i ( X 0 ) + j = 0 N R i j ( X j - X j 0 ) ,

where R is a square matrix of size N and ℛ(X0) is the radiative budget for the profile X0. If we assume R to be invertible, the energy flux can be computed with

(B4) X i = X i 0 - j = 0 N R i j - 1 F j - F j + 1 + R j ( X 0 ) .

B2 Flux–mass relation

We first consider the dry static energy ed=CpT+gz. Considering the atmosphere is an ideal gas at hydrostatic equilibrium, and for prescribed pressure levels, layer volume depends only on temperature. Therefore, the elevation of a layer is a function of temperatures of layers below only and we can express the energy of a layer as (see Appendix A)

(B5) e i d = k = 0 N E i k d T k ,

where Ed is a square, triangular matrix of size N. If we linearize around X0, one obtains

(B6) e i d k = 0 N E i k d X k 0 2 ( 2 X k 0 - X k ) .

Using Eq. (B4) gives the expression of energy as a function of flux

(B7) e i d ( F ) k = 0 N E i k d X k 0 2 X k 0 + j = 0 N R k j - 1 F j - F j + 1 + R j ( X 0 ) .

We also can take into account the latent heat for a water-vapour-saturated atmosphere. The mixing ratio at the saturation point, qs, depends only on temperature T (in K) and pressure p (in Pa). It is given by the Bolton equation (Bolton1980):

(B8) q s ( T , p ) = 622.0 h s ( T ) p - h s ( T ) with h s ( T ) = 6.112 exp 17.62 ( T - 273.15 ) T - 30.03 ,

where hs is saturation vapour pressure of the mixture (in Pa). At fixed pressure, the moist static energy at saturation es=CpT+gz+Lqs of layers is only a function of temperatures. If we linearize around the profile X0,

(B9) e i s = e i d + L q s 1 X i e i d + L q s 1 X i 0 - L q s T 1 X i 0 X i - X i 0 X i 0 2 .

Then, we can use the same reasoning as for the dry static heat and replace the matrix Ed by Es to consider the effect of latent energy for a saturated moisture profile. However, the radiative budget is still computed with reference water vapour profiles. In the following, e can represent ed or es.

B3 Constraint

By multiplying both sides of

(B10) F i = - m i e i - e i - 1

by ei-ei-1, we obtain

(B11) F i e i - e i - 1 = - m i e i - e i - 1 2 .

So the constraint mi≥0 is equivalent to

(B12) α i ( F ) - F i e i ( F ) - e i - 1 ( F ) 0 .

B4 Associated Lagrangian in flux space

Using the linearized energy budget (Eq. B4) and the constraint (Eq. B12), the problem (Eq. 11) is supposed to be equivalent to the search of critical points of the following Lagrangian

(B13) L ( F , μ ) = σ ( F ) - i = 1 N μ i α i ( F ) = i = 0 N X i F i - F i + 1 - i = 1 N μ i α i

(B14) i = 0 N X i 0 - j = 0 N R i j - 1 F j - F j + 1 + R j ( X 0 ) F i - F i + 1 - i = 1 N μ i α i ( F ) ,

while respecting the Karush–Kuhn–Tucker (KKT) conditions

(B15) L F i = 0 with α i ( F ) 0 , μ i 0 , and μ i α i ( F ) = 0 for i = 1 , , N .

The problem is solved numerically by using an interior point method (Boyd and Vandenberghe2004).

B5 Convergence of the algorithm

In practice, the algorithm may fail to find the global optimum for large N (≥50) but is robust for N≤40. As N increases, the algorithm converges rapidly to a solution (Fig. B1). The choice N=20 is a good compromise between computation time and resolution.

Figure B1Energy flux F, specific energy e, and temperature T computed by our constrained model with e=CpT+gz+Lqs for imposed tropical atmospheric composition measured by McClatchey et al. (1972), with various resolutions N.


Author contributions

DP developed the model code. VL performed the simulations and prepared the manuscript under the supervision of DP and BD.

Competing interests

The authors declare that they have no conflict of interest.

Special issue statement

This article is part of the special issue “Thermodynamics and optimality in the Earth system and its subsystems (ESD/HESS inter-journal SI)”. It is not associated with a conference.


We thank Stan Schymanski and the anonymous reviewer for their useful comments that helped us to improve the article.

Review statement

This paper was edited by Michel Crucifix and reviewed by Stan Schymanski and one anonymous referee.


Bolton, D.: The Computation of Equivalent Potential Temperature, Mon. Weather Rev., 108, 1046–1053,<1046:TCOEPT>2.0.CO;2, 1980. a

Boyd, S. and Vandenberghe, L.: Convex Optimization, Cambridge University Press, New York, NY, USA, 2004. a, b

Dewar, R.: Information theory explanation of the fluctuation theorem, maximum entropy production and self-organized criticality in non-equilibrium stationary states, J. Phys. A-Math. Gen., 36, 631–641,, 2003. a

Dewar, R. C.: Maximum Entropy Production as an Inference Algorithm that Translates Physical Assumptions into Macroscopic Predictions: Don’t Shoot the Messenger, Entropy, 11, 931–944,, 2009. a

Dixit, P. D., Wagoner, J., Weistuch, C., Pressé, S., Ghosh, K., and Dill, K. A.: Perspective: Maximum caliber is a general variational principle for dynamical systems, J. Chem. Phys., 148, 010901,, 2018. a

Dufresne, J.-L. and Bony, S.: An Assessment of the Primary Sources of Spread of Global Warming Estimates from Coupled Atmosphere–Ocean Models, J. Climate, 21, 5135–5144,, 2008. a, b, c, d, e, f, g

Dufresne, J.-L., Fournier, R., Hourdin, C., and Hourdin, F.: Net Exchange Reformulation of Radiative Transfer in the CO2 15-µm Band on Mars, J. Atmos. Sci., 62, 3303–3319,, 2005. a

Dufresne, J.-L., Foujols, M.-A., Denvil, S., Caubel, A., Marti, O., Aumont, O., Balkanski, Y., Bekki, S., Bellenger, H., Benshila, R., Bony, S., Bopp, L., Braconnot, P., Brockmann, P., Cadule, P., Cheruy, F., Codron, F., Cozic, A., Cugnet, D., de Noblet, N., Duvel, J.-P., Ethé, C., Fairhead, L., Fichefet, T., Flavoni, S., Friedlingstein, P., Grandpeix, J.-Y., Guez, L., Guilyardi, E., Hauglustaine, D., Hourdin, F., Idelkadi, A., Ghattas, J., Joussaume, S., Kageyama, M., Krinner, G., Labetoulle, S., Lahellec, A., Lefebvre, M.-P., Lefevre, F., Levy, C., Li, Z. X., Lloyd, J., Lott, F., Madec, G., Mancip, M., Marchand, M., Masson, S., Meurdesoif, Y., Mignot, J., Musat, I., Parouty, S., Polcher, J., Rio, C., Schulz, M., Swingedouw, D., Szopa, S., Talandier, C., Terray, P., Viovy, N., and Vuichard, N.: Climate change projections using the IPSL-CM5 Earth System Model: from CMIP3 to CMIP5, Clim. Dynam., 40, 2123–2165,, 2013. a

Forster, P. M. F. and Gregory, J. M.: The Climate Sensitivity and Its Components Diagnosed from Earth Radiation Budget Data, J. Climate, 19, 39–52,, 2006. a

Goody, R.: Maximum Entropy Production in Climate Theory, J. Atmos. Sci., 64, 2735–2739,, 2007. a

Goody, R. M.: A statistical model for water-vapour absorption, Q. J. Roy. Meteor. Soc., 78, 165–169,, 1952. a

Goosse, H., Brovkin, V., Fichefet, T., Haarsma, R., Huybrechts, P., Jongma, J., Mouchet, A., Selten, F., Barriat, P.-Y., Campin, J.-M., Deleersnijder, E., Driesschaert, E., Goelzer, H., Janssens, I., Loutre, M.-F., Morales Maqueda, M. A., Opsteegh, T., Mathieu, P.-P., Munhoven, G., Pettersson, E. J., Renssen, H., Roche, D. M., Schaeffer, M., Tartinville, B., Timmermann, A., and Weber, S. L.: Description of the Earth system model of intermediate complexity LOVECLIM version 1.2, Geosci. Model Dev., 3, 603–633,, 2010. a

Herbert, C.: Applications de la mécanique statistique à la modélisation du climat: thermodynamique et dynamique de l'atmosphère, PhD thesis, available at: (last access: 20 October 2018), thèse de doctorat dirigée par Paillard, Didier et Dubrulle, Bérengère Géophysique Paris 6, 2012. a

Herbert, C., Paillard, D., and Dubrulle, B.: Entropy production and multiple equilibria: the case of the ice-albedo feedback, Earth Syst. Dynam., 2, 13–23,, 2011a. a

Herbert, C., Paillard, D., Kageyama, M., and Dubrulle, B.: Present and Last Glacial Maximum climates as states of maximum entropy production, Q. J. Roy. Meteor. Soc., 137, 1059–1069,, 2011b. a

Herbert, C., Paillard, D., and Dubrulle, B.: Vertical Temperature Profiles at Maximum Entropy Production with a Net Exchange Radiative Formulation, J. Climate, 26, 8545–8555,, 2013. a, b, c, d, e, f, g, h, i, j, k

Hourdin, F., Mauritsen, T., Gettelman, A., Golaz, J.-C., Balaji, V., Duan, Q., Folini, D., Ji, D., Klocke, D., Qian, Y., Rauser, F., Rio, C., Tomassini, L., Watanabe, M., and Williamson, D.: The Art and Science of Climate Model Tuning, B. Am. Meteorol. Soc., 98, 589–602,, 2017. a

Jaynes, E. T.: Information Theory and Statistical Mechanics, Phys. Rev., 106, 620–630,, 1957. a

Karkar, S. and Paillard, D.: Inferring global wind energetics from a simple Earth system model based on the principle of maximum entropy production, Earth Syst. Dynam. Discuss., 6, 407–433,, 2015. a

Kleidon, A.: A basic introduction to the thermodynamics of the Earth system far from equilibrium and maximum entropy production, Philos. T. Roy. Soc. B, 365, 1303–1315,, 2010. a

Lacis, A. A. and Hansen, J.: A Parameterization for the Absorption of Solar Radiation in the Earth's Atmosphere, J. Atmos. Sci., 31, 118–133,<0118:APFTAO>2.0.CO;2, 1974. a

Labarre, V., Paillard, D., and Dubrulle, B.: Supplementary Materials to: A Radiative Convective Model based on constrained Maximum Entropy Production, Zenodo, 18 March 2019,, 2019. a

Lorenz, R. D., Lunine, J. I., Withers, P. G., and McKay, C. P.: Titan, Mars and Earth: Entropy production by latitudinal heat transport, Geophys. Res. Lett., 28, 415–418,, 2001. a

Manabe, S. and Strickler, R. F.: Thermal Equilibrium of the Atmosphere with a Convective Adjustment, J. Atmos. Sci., 21, 361–385,<0361:TEOTAW>2.0.CO;2, 1964. a, b

Martyushev, L. and Seleznev, V.: Maximum entropy production principle in physics, chemistry and biology, Physics Reports, 426, 1–45,, 2006. a

McClatchey, R. A., Fenn, R. W., Selby, J. E. A., Volz, F. E., and Garing, J. S.: Optical properties of the atmosphere, 3rd Edn., Tech. Rep. AFCRL-72–0497, 108 pp., Air Force Geophys. Lab., Hanscom AFB, Mass., 1972. a, b, c, d, e, f, g, h, i, j, k

Mellor, G. L. and Yamada, T.: Development of a turbulence closure model for geophysical fluid problems, Rev. Geophys., 20, 851–875,, 1982. a

Mihelich, M.: Vers une compréhension du principe de maximisation de production d'entropie, PhD thesis, available at: (last access: 20 October 2018), thèse de doctorat dirigée par Dubrulle, Bérengère Physique Paris Saclay 2015, 2015. a

Monthus, C.: Non-equilibrium steady states: maximization of the Shannon entropy associated with the distribution of dynamical trajectories in the presence of constraints, J. Stat. Mech.-Theory E., 2011, P03008,, 2011. a

North, G. R., Cahalan, R. F., and Coakley Jr., J. A.: Energy balance climate models, Rev. Geophys., 19, 91–121,, 1981. a, b

O'brien, D. M. and Stephens, G. L.: Entropy and climate. II: Simple models, Q. J. Roy. Meteor. Soc., 121, 1773–1796,, 1995.  a

Ozawa, H. and Ohmura, A.: Thermodynamics of a Global-Mean State of the Atmosphere – A State of Maximum Entropy Increase, J. Climate, 10, 441–445,<0441:TOAGMS>2.0.CO;2, 1997. a, b, c

Paillard, D.: The timing of Pleistocene glaciations from a simple multiple-state climate model, Nature, 391, 378–381, 1998. a

Paltridge, G. W.: Global dynamics and climate – a system of minimum entropy exchange, Q. J. Roy. Meteor. Soc., 101, 475–484,, 1975. a

Pascale, S., Gregory, J. M., Ambaum, M. H. P., Tailleux, R., and Lucarini, V.: Vertical and horizontal processes in the global atmosphere and the maximum entropy production conjecture, Earth Syst. Dynam., 3, 19–32,, 2012. a

Peixoto, J. P. and Oort, A. H.: Physics of climate, American Institute of Physics, New York, 1992. a

Pujol, T. and Fort, J.: States of maximum entropy production in a onedimensional vertical model with convective adjustment, Tellus A, 54, 363–369,, 2002. a, b

Randall, D. A., Wood, R. A., Bony, S., Colman, R., Fichefet, T., Fyfe, J., Kattsov, V., Pitman, A., Shukla, J., Srinivasan, J., Stouffer, R. J., Sumi, A., and Taylor, K. E.: Climate models and their evaluation, in: Climate Change 2007: The Physical Science Basis. Contribution of Working Group I to the fourth assessment report of the Intergovernmental Panel on Climate Change, edited by: Solomon, S., Qin, D., Manning, M., Chen, Z., Marquiz, M., Averyt, K. B., Tignor, M., and Miller, H. L., Cambridge University Press, Cambridge, United Kingdom/New York, NY, USA, 2007. a, b

Rieutord, M.: Fluid dynamics: an introduction, Graduate Texts in Physics, Springer, ISBN 978-3-319-09351-2, 2015. a

Rodgers, C. D.: Comments on paltridge's “minimum entropy exchange” principle, Q. J. Roy. Meteor. Soc., 102, 455–458,, 1976. a

Rodgers, C. D. and Walshaw, C. D.: The computation of infra-red cooling rate in planetary atmospheres, Q. J. Roy. Meteor. Soc., 92, 67–92,, 1966. a

Stevens, B. and Bony, S.: What Are Climate Models Missing?, Science, 340, 1053–1054,, 2013. a

Yano, J.-I., Bister, M., Fuchs, Ž., Gerard, L., Phillips, V. T. J., Barkidija, S., and Piriou, J.-M.: Phenomenology of convection-parameterization closure, Atmos. Chem. Phys., 13, 4111–4131,, 2013. a

Short summary
We tried to represent atmospheric convection induced by radiative forcing with a simple climate model based on maximum entropy production. Contrary to previous models, we give a minimal description of energy transport in the atmosphere. It allows us to give better results in terms of temperature and vertical energy flux profiles.
Final-revised paper