Temperatures from energy balance models: the effective heat capacity matters

Energy balance models (EBMs) are highly simplified models of the climate system, providing admissible conceptual tools for understanding climate changes. The global temperature is calculated by the radiation budget through the incoming energy from the Sun and the outgoing energy from the Earth. The argument that the temperature can be calculated by this simple radiation budget is revisited. The underlying assumption for a realistic temperature distribution is explored: one has to assume a moderate diurnal cycle due to the large heat capacity and the fast rotation of the Earth. Interestingly, the global mean in the revised EBM is very close to the originally proposed value. The main point is that the effective heat capacity and its temporal variation over the daily and seasonal cycle needs to be taken into account when estimating surface temperature from the energy budget. Furthermore, the time-dependent EBM predicts a flat meridional temperature gradient for large heat capacities, reducing the seasonal cycle and the outgoing radiation and increasing global temperature. Motivated by this finding, a sensitivity experiment with a complex model is performed where the vertical diffusion in the ocean has been increased. The resulting temperature gradient, reduced seasonal cycle, and global warming is also found in climate reconstructions, providing a possible mechanism for past climate changes prior to 3 million years ago.

1 Introduction

Energy balance models (EBMs) are among the simplest climate models. They were introduced almost simultaneously by Budyko (1969) and Sellers (1969). Because of their simplicity, these models are easy to understand and facilitate both analytical and numerical studies of climate sensitivity (Peixoto and Oort, 1992; Hartmann, 1994; Saltzmann, 2001; Ruddiman, 2001; Pierrehumbert, 2010). A key feature of these models is that they eliminate the climate's dependence on the wind field, ocean currents and the Earth's rotation, and thus have only one dependent variable: the Earth's near-surface air temperature T.

With the development of computer capacities, simpler models have not disappeared; on the contrary, a stronger emphasis has been given to the concept of a hierarchy of models as the only way to provide a linkage between theoretical understanding and the complexity of realistic models (von Storch et al. 1999; Claussen et al. 2002). In contrast, many important scientific debates in recent years have had their origin in the use of conceptually simple models (Le Treut et al., 2007; Stocker, 2011) and also as a way to analyze data (Köhler et al., 2010) or complex models (Knorr et al., 2011).

Figure 1 Schematic view of the energy absorbed and emitted by the Earth following Eq. (1). Modified following Goose (2015).

Pioneering work has been done by North (North, 1975a, b, 1981, 1983) and these models were applied subsequently (e.g., Ghil, 1976; Su and Hsieh, 1976; Fraedrich, 1979; Ghil and Childress, 1987; Short et al., 1991; Stocker et al., 1992; North and Kim, 2017). Later the EBMs were equipped by the hydrological cycle (Chen et al., 1995; Lohmann et al., 1996; Fanning and Weaver, 1996; Lohmann and Gerdes, 1998) to study the feedbacks in the atmosphere–ocean–sea ice system. One of the most useful examples of a simple but powerful model is the one- or zero-dimensional energy balance model. As a starting point, a zero-dimensional model of the radiative equilibrium of the Earth is introduced (Fig. 1)

$$\begin{array}{}\text{(1)}& (\mathrm{1}-\mathit{\alpha })S\mathit{\pi }{R}^{\mathrm{2}}=\mathrm{4}\mathit{\pi }{R}^{\mathrm{2}}\mathit{ϵ}\mathit{\sigma }{T}^{\mathrm{4}},\end{array}$$

where the left-hand side represents the incoming energy from the Sun (size of the disk is equal to the shadow area πR²), while the right-hand side represents the outgoing energy from the Earth (Fig. 1). T is calculated from the Stefan–Boltzmann law assuming a constant radiative temperature; S is the solar constant (the incoming solar radiation per unit area), about 1367 Wm⁻²; and α is the Earth's average planetary albedo, measured to be 0.3. R is Earth's radius =  6.371×10⁶ m, σ is the Stefan–Boltzmann constant $=\mathrm{5.67}\times {\mathrm{10}}^{-\mathrm{8}}\mathrm{J}{\mathrm{K}}^{-\mathrm{4}}{\mathrm{m}}^{-\mathrm{2}}{\mathrm{s}}^{-\mathrm{1}}$, and ϵ is the effective emissivity of Earth (about 0.612) (e.g., Archer, 2009). The geometrical constant πR² can be factored out, giving

$$\begin{array}{}\text{(2)}& (\mathrm{1}-\mathit{\alpha })S=\mathrm{4}\mathit{ϵ}\mathit{\sigma }{T}^{\mathrm{4}}.\end{array}$$

Solving for the temperature,

$$\begin{array}{}\text{(3)}& T=\sqrt[\mathrm{4}]{{\displaystyle \frac{(\mathrm{1}-\mathit{\alpha })S}{\mathrm{4}\mathit{ϵ}\mathit{\sigma }}}}.\end{array}$$

Since the use of the effective emissivity ϵ in (1) already accounts for the greenhouse effect, we gain an average Earth temperature of 288 K (15 ^∘C), very close to the global temperature observations and reconstructions (Hansen et al., 2010) at 14 ^∘C for 1951–1980. The implicit assumption in Eqs. (2) and (3) is that we have a single temperature on the Earth, although we know that the Equator-to-pole surface temperature gradient is on the order of 50 K and that the incoming solar radiation at the Equator is about twice that at the poles (Peixoto and Oort, 1992).

Furthermore, Eq. (3) does not contain parameters like the heat capacity of the planet. We will explore the idea that this is essential for the temperature of the Earth's climate system. We will evaluate the effect of the effective heat capacity in the climate system. Schwartz (2007) stressed that the effective heat capacity is not an intrinsic property of the climate system but is reflective of the rate of penetration of heat energy into the ocean in response to the particular pattern of forcing and background state.

Wang et al. (2019) showed that a pronounced low Equator-to-pole gradient in the annual mean sea surface temperatures is found in a numerical experiment conducted with a coupled model consisting of an atmospheric general circulation model coupled to a slab ocean model in which the mixed-layer thickness is reduced. In the present paper, it is shown that the heat capacity is linked to the question of a low Equator-to-pole gradients during the Paleogene–Neogene climate (Markwick, 1994; Wolfe, 1994; Sloan and Rea, 1996; Huber et al., 2000; Shellito et al., 2003; Tripati et al., 2003; Mosbrugger et al., 2005). These published temperature patterns resemble the high-latitude warming (with moderate low-latitude warming) and reduced seasonality.

2 Results

2.1 The spatial distribution of the energy balance

Let us have a closer look onto Eq. (1) and consider the local radiative equilibrium of the Earth at each point. Figure 2 shows the latitude–longitude dependence of the incoming shortwave radiation. The global mean temperatures are not affected by the tilt (Berger and Loutre 1991, 1997; Laepple and Lohmann 2009). We assume an idealized geometry of the Earth, no obliquity and no precession, which makes an analytical calculation possible. Furthermore, ϵ and α are assumed to be constants.

Figure 2 Latitudinal (a) and longitudinal (b) dependence of the incoming shortwave radiation. On the right-hand side, the insolation as a function of latitude φ and longitude Θ with maximum insolation (1−α)S is shown. See the text for the details.

The incoming radiation goes with the cosine of latitude φ and longitude Θ, and there is only sunshine during the day. Figure 2a shows the latitudinal dependence. As we assume no tilt (this assumption is later relaxed), the latitudinal dependence is a function of latitude only: cos φ. On the right-hand side, the function is shown. Figure 2b shows that the latitudinal dependence is a function of longitude: cos Θ for the sun-exposed side of the Earth, and for the dark side of the Earth it is zero. For simplicity, we can define the angle Θ anti-clockwise for the sun-exposed side between $-\mathit{\pi }/\mathrm{2}$ and π∕2. We always define the maximal insolation at Θ=0, which is moving in time. In the panel, the Earth's rotation is schematically sketched as the red arrow, and we see the time dependence on the right-hand side. It is noted that the geographical longitude can be calculated by mod($\mathrm{\Theta }-\mathrm{2}\mathit{\pi }\cdot t/\mathrm{24},\mathrm{2}\mathit{\pi }$) where t is measured in hours and mod is the modulo operation. Summarizing our geometrical considerations, we can now write the local energy balance as follows:

$$\begin{array}{}\text{(4)}& \begin{array}{lll}\mathit{ϵ}\mathit{\sigma }{T}^{\mathrm{4}}& =(\mathrm{1}-\mathit{\alpha })S\cdot \cos \mathit{\phi }\cdot \cos \mathrm{\Theta }& \text{for }-\mathit{\pi }/\mathrm{2}<\mathrm{\Theta }<\mathit{\pi }/\mathrm{2}.\\ & =\mathrm{0}& \text{elsewhere }\end{array}\end{array}$$

Temperatures based on the local energy balance without a heat capacity would vary between T_min=0 K and $T_{\max }=\sqrt[\mathrm{4}]{\frac{(\mathrm{1}-\mathit{\alpha })S}{\mathit{ϵ}\mathit{\sigma }}}=\sqrt{\mathrm{2}}\cdot \sqrt[\mathrm{4}]{\frac{(\mathrm{1}-\mathit{\alpha })S}{\mathrm{4}\mathit{ϵ}\mathit{\sigma }}}=\sqrt{\mathrm{2}}\cdot \mathrm{288}\mathrm{K}=\mathrm{407}$ K.

Integration of Eq. (4) over the Earth surface is as follows:

$$\begin{array}{}\text{(5)}& \begin{array}{rl}& \underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}\left(\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}\mathit{ϵ}\mathit{\sigma }{T}^{\mathrm{4}}R\cos \mathit{\phi }d\mathrm{\Theta }\right)Rd\mathit{\phi }=\\ & (\mathrm{1}-\mathit{\alpha })S\underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}R{\cos }^{\mathrm{2}}\mathit{\phi }d\mathit{\phi }\cdot \underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}R\cos \mathrm{\Theta }d\mathrm{\Theta }\\ & \mathit{ϵ}\mathit{\sigma }{R}^{\mathrm{2}}{\displaystyle \frac{\mathrm{4}\mathit{\pi }}{\mathrm{4}\mathit{\pi }}}\underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}\left(\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}{T}^{\mathrm{4}}\cos \mathit{\phi }d\mathrm{\Theta }\right)d\mathit{\phi }=\\ & (\mathrm{1}-\mathit{\alpha })S{R}^{\mathrm{2}}\underset{\frac{\mathit{\pi }}{\mathrm{2}}}{\underbrace{\underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}{\cos }^{\mathrm{2}}\mathit{\phi }d\mathit{\phi }}}\cdot \underset{\mathrm{2}}{\underbrace{\underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}\cos \mathrm{\Theta }d\mathrm{\Theta }}}\\ & \mathit{ϵ}\mathit{\sigma }\mathrm{4}\mathit{\pi }\stackrel{\mathrm{‾}}{T^{\mathrm{4}}}=(\mathrm{1}-\mathit{\alpha })S\mathit{\pi }\end{array},\end{array}$$

giving a similar formula as that in Eq. (3) with the definition for the average

$$\stackrel{\mathrm{‾}}{T^{\mathrm{4}}}={\displaystyle \frac{\mathrm{1}}{\mathrm{4}\mathit{\pi }}}\underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}d\mathit{\phi }\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}\cos \mathit{\phi }d\mathrm{\Theta }{T}^{\mathrm{4}}.$$

What we really want is the mean of the temperature $\stackrel{\mathrm{‾}}{T}$. Therefore, we take the fourth root of Eq. (4):

$$\begin{array}{}\text{(6)}& \begin{array}{rl}T& =\sqrt[\mathrm{4}]{{\displaystyle \frac{(\mathrm{1}-\mathit{\alpha })S\cos \mathit{\phi }\cos \mathrm{\Theta }}{\mathit{ϵ}\mathit{\sigma }}}}\text{for }-\mathit{\pi }/\mathrm{2}<\mathrm{\Theta }<\mathit{\pi }/\mathrm{2},\\ & \text{and zero elsewhere.}\end{array}\end{array}$$

If we calculate the zonal mean of Eq. (6) by integration at the latitudinal cycles we have

$$\begin{array}{}\text{(7)}& \begin{array}{rl}T\left(\mathit{\phi }\right)& ={\displaystyle \frac{\sqrt{\mathrm{2}}}{\mathrm{2}\mathit{\pi }}}\underset{\sqrt{\mathit{\pi }}\mathrm{\Gamma }(\mathrm{5}/\mathrm{8})/\mathrm{\Gamma }(\mathrm{9}/\mathrm{8})}{\underbrace{\underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}(\cos \mathrm{\Theta }{)}^{\mathrm{1}/\mathrm{4}}d\mathrm{\Theta }}}\sqrt[\mathrm{4}]{{\displaystyle \frac{(\mathrm{1}-\mathit{\alpha })S}{\mathrm{4}\mathit{ϵ}\mathit{\sigma }}}}(\cos \mathit{\phi }{)}^{\mathrm{1}/\mathrm{4}}\\ & =\underset{\approx \mathrm{0.608}}{\underbrace{{\displaystyle \frac{\mathrm{1}}{\sqrt{\mathrm{2}\mathit{\pi }}}}{\displaystyle \frac{\mathrm{\Gamma }(\mathrm{5}/\mathrm{8})}{\mathrm{\Gamma }(\mathrm{9}/\mathrm{8})}}}}\cdot \sqrt[\mathrm{4}]{{\displaystyle \frac{(\mathrm{1}-\mathit{\alpha })S}{\mathrm{4}\mathit{ϵ}\mathit{\sigma }}}}(\cos \mathit{\phi }{)}^{\mathrm{1}/\mathrm{4}}\end{array},\end{array}$$

as a function of latitude (Fig. 3). Γ is Euler's Gamma function with $\mathrm{\Gamma }(x+\mathrm{1})=x\mathrm{\Gamma }\left(x\right)$. When we integrate this over the latitudes,

$$\begin{array}{}\text{(8)}& \begin{array}{rl}\stackrel{\mathrm{‾}}{T}& ={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}T\left(\mathit{\phi }\right)\cos \mathit{\phi }d\mathit{\phi }\\ & ={\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{\displaystyle \frac{\mathrm{\Gamma }(\mathrm{5}/\mathrm{8})}{\sqrt{\mathrm{2}\mathit{\pi }}\mathrm{\Gamma }(\mathrm{9}/\mathrm{8})}}\cdot \sqrt[\mathrm{4}]{{\displaystyle \frac{(\mathrm{1}-\mathit{\alpha })S}{\mathrm{4}\mathit{ϵ}\mathit{\sigma }}}}\underset{\sqrt{\mathit{\pi }}\mathrm{\Gamma }(\mathrm{9}/\mathrm{8})/\mathrm{\Gamma }(\mathrm{13}/\mathrm{8})}{\underbrace{\underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}(\cos \mathit{\phi }{)}^{\mathrm{5}/\mathrm{4}}d\mathit{\phi }}}\\ & =\underset{\frac{\sqrt{\mathrm{2}}}{\mathrm{4}}\frac{\mathrm{8}}{\mathrm{5}}=\mathrm{0.4}\sqrt{\mathrm{2}}}{\underbrace{{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}{\displaystyle \frac{\mathrm{1}}{\sqrt{\mathrm{2}}}}{\displaystyle \frac{\mathrm{\Gamma }(\mathrm{5}/\mathrm{8})}{\mathrm{\Gamma }(\mathrm{13}/\mathrm{8})}}}}\cdot \sqrt[\mathrm{4}]{{\displaystyle \frac{(\mathrm{1}-\mathit{\alpha })S}{\mathrm{4}\mathit{ϵ}\mathit{\sigma }}}}.\end{array}\end{array}$$

Therefore, the mean temperature is a factor $F=\mathrm{0.4}\sqrt{\mathrm{2}}\approx \mathrm{0.566}$ lower than 288 K, as stated by (3) and would be $\stackrel{\mathrm{‾}}{T}=\mathrm{163}$ K. The standard EBM in Fig. 1 has been imprinted into our thoughts and lectures. We should therefore be careful and pinpoint the reasons for the failure. What happens here is that the heat capacity of the Earth is neglected and there is a strong nonlinearity in the outgoing radiation. The local radiative equilibrium assumption cannot be used under these conditions.

2.2 The heat capacity and fast rotating body

If the heat capacity is taken into account, the energy balance reads as follows:

$$\begin{array}{}\text{(9)}& \begin{array}{lll}{C}_{\mathrm{p}}{\partial }_{t}T& =(\mathrm{1}-\mathit{\alpha })S& \cdot \cos \mathit{\phi }\cdot \cos \mathrm{\Theta }-\mathit{ϵ}\mathit{\sigma }{T}^{\mathrm{4}}\\ & & \text{for }-\mathit{\pi }/\mathrm{2}<\mathrm{\Theta }<\mathit{\pi }/\mathrm{2}\\ & =& -\mathit{ϵ}\mathit{\sigma }{T}^{\mathrm{4}}\\ & & \text{elsewhere }\end{array},\end{array}$$

with C_p representing the heat capacity multiplied with the depth of the atmosphere–ocean layer (C_p is on the order of 10⁷−10⁸ J K⁻¹ m⁻²). To simplify Eq. (9), the energy balance is integrated over the longitude and day, $\stackrel{\mathrm{̃}}{T}=\frac{\mathrm{1}}{\mathrm{2}\mathit{\pi }}\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}Td\mathrm{\Theta }.$ We assume that the exchange of operations in zonal and time mean and exponentials to the power of 4 is valid:

$$\begin{array}{}\text{(10)}& {\displaystyle \frac{\mathrm{1}}{\mathrm{2}\mathit{\pi }}}\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}{T}^{\mathrm{4}}\approx {\stackrel{\mathrm{̃}}{T}}^{\mathrm{4}},\end{array}$$

which is a good approximation in the climate system (Fig. A1) when analyzing hourly data temperature (Hersbach et al., 2020). Figure A1b indicates that the difference is lowest at low latitudes or marine-dominated regions where the diurnal variation of sea surface temperatures is small (e.g., Kawai and Kawamura, 2002; Stommel et al., 1969; Stuart-Menteth, et al. 2003; Ward, 2006). Therefore,

$$\begin{array}{}\text{(11)}& \begin{array}{rl}{C}_{\mathrm{p}}{\partial }_t\stackrel{\mathrm{̃}}{T}& =(\mathrm{1}-\mathit{\alpha })S\cos \mathit{\phi }\cdot {\displaystyle \frac{\mathrm{1}}{\mathrm{2}\mathit{\pi }}}\underset{\mathrm{2}}{\underbrace{\underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}\cos \mathrm{\Theta }d\mathrm{\Theta }}}-\mathit{ϵ}\mathit{\sigma }\underset{\approx {\stackrel{\mathrm{̃}}{T}}^{\mathrm{4}}}{\underbrace{{\displaystyle \frac{\mathrm{1}}{\mathrm{2}\mathit{\pi }}}\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}{T}^{\mathrm{4}}d\mathrm{\Theta }}}\\ & =(\mathrm{1}-\mathit{\alpha }){\displaystyle \frac{S}{\mathit{\pi }}}\cos \mathit{\phi }-\mathit{ϵ}\mathit{\sigma }{\stackrel{\mathrm{̃}}{T}}^{\mathrm{4}}.\end{array}\end{array}$$

It is emphasized that the assumption in Eq. (10) is made for the averaging of the diurnal cycle and longitude. It is completely different from the approach in Eqs. (3) and (5). Note that $\mathit{ϵ}\mathit{\sigma }{\stackrel{\mathrm{̃}}{T}}^{\mathrm{4}}$ has a pronounced latitudinal dependence, varying over 2 orders of magnitude. Treating this term as a constant number as in Eq. (3) is questionable.

The right-hand side of Eq. (11) is now time independent and the equilibrium solution is simply

$$\begin{array}{}\text{(12)}& \stackrel{\mathrm{̃}}{T}\left(\mathit{\phi }\right)=\sqrt[\mathrm{4}]{{\displaystyle \frac{\mathrm{4}}{\mathit{\pi }}}}\cdot \sqrt[\mathrm{4}]{{\displaystyle \frac{(\mathrm{1}-\mathit{\alpha })S}{\mathrm{4}\mathit{ϵ}\mathit{\sigma }}}}(\cos \mathit{\phi }{)}^{\mathrm{1}/\mathrm{4}},\end{array}$$

shown as the red line in Fig. 3. The global mean temperature is

$$\begin{array}{}\text{(13)}& \begin{array}{rl}\stackrel{\mathrm{‾}}{\stackrel{\mathrm{̃}}{T}}& =\sqrt[\mathrm{4}]{{\displaystyle \frac{\mathrm{4}}{\mathit{\pi }}}}\cdot \sqrt[\mathrm{4}]{{\displaystyle \frac{(\mathrm{1}-\mathit{\alpha })S}{\mathrm{4}\mathit{ϵ}\mathit{\sigma }}}}{\displaystyle \frac{\mathrm{1}}{\mathrm{2}}}\underset{\sqrt{\mathit{\pi }}\mathrm{\Gamma }(\mathrm{9}/\mathrm{8})/\mathrm{\Gamma }(\mathrm{13}/\mathrm{8})}{\underbrace{\underset{-\mathit{\pi }/\mathrm{2}}{\overset{\mathit{\pi }/\mathrm{2}}{\int }}(\cos \mathit{\phi }{)}^{\mathrm{5}/\mathrm{4}}d\mathit{\phi }}}\\ & =\sqrt{{\displaystyle \frac{\mathit{\pi }}{\mathrm{2}}}}{\displaystyle \frac{\mathrm{\Gamma }(\mathrm{9}/\mathrm{8})}{\mathrm{\Gamma }(\mathrm{13}/\mathrm{8})}}\cdot \sqrt[\mathrm{4}]{{\displaystyle \frac{(\mathrm{1}-\mathit{\alpha })S}{\mathrm{4}\mathit{ϵ}\mathit{\sigma }}}},\end{array}\end{array}$$

with the factor $G=\sqrt{\frac{\mathit{\pi }}{\mathrm{2}}}\frac{\mathrm{\Gamma }(\mathrm{9}/\mathrm{8})}{\mathrm{\Gamma }(\mathrm{13}/\mathrm{8})}\approx \mathrm{0.989}$ and $\stackrel{\mathrm{‾}}{\stackrel{\mathrm{̃}}{T}}=\mathrm{285}\text{ K}$.

Alternatively, one can obtain the numerical solution of Eq. (9), shown as the dashed brownish line in Fig. 3, where the diurnal cycle has been explicitly taken into account. Here, $C_{\mathrm{p}}=C_{\mathrm{p}}^{\mathrm{a}}$ has been chosen as the atmospheric heat capacity $C_{\mathrm{p}}^{\mathrm{a}}=c_{\mathrm{p}}{p}_{\mathrm{s}}/g=\mathrm{1004}\mathrm{J}{\mathrm{K}}^{-\mathrm{1}}{\mathrm{kg}}^{-\mathrm{1}}\times {\mathrm{10}}^{\mathrm{5}}\mathrm{Pa}/\left(\mathrm{9.81}\mathrm{m}{\mathrm{s}}^{-\mathrm{2}}\right)=\mathrm{1.02}\times {\mathrm{10}}^{\mathrm{7}}$ J K⁻¹ m⁻² which is the specific heat at constant pressure c_p times the total mass p_s∕g, p_s is the surface pressure, and g is the gravity. The global mean temperature $\stackrel{\mathrm{‾}}{T}$ is 286 K, again close to 288 K.

Figure 3 Latitudinal temperatures of the EBM with zero heat capacity (Eq. 7) shown in cyan (its mean shown as a dashed line), the global approach (Eq. 3) shown as a solid black line, and the zonal and time averaging (Eq. 12) shown in red. Please note that both curves are exactly $\mathrm{0}K=-\mathrm{273.15}$ ^∘C at the poles, but those values were not shown for graphical reasons. The dashed brownish curve shows the numerical solution by taking the zonal mean of Eq. (9).

The effect of heat capacity is systematically analyzed in Fig. 5. The temperatures are relatively insensitive for a wide range of C_p. We find a severe drop in temperatures for heat capacities below 0.01 of the atmospheric heat capacity $C_{\mathrm{p}}^{\mathrm{a}}$. Figure 4 shows the temperature dependence for different values of C_p and the length of the day, indicating a pronounced temperature drop during the night for low values of heat capacities and for (hypothetical) long days of 240 h instead of 24 h. We have chosen this feature for a particular latitude (here 45^∘ N). The analysis shows that the effective heat capacity is of great importance for the temperature, this depends on the atmospheric planetary boundary layer (how well mixed it is with small gradients in the vertical) and the depth of the mixed layer in the ocean, which will be analyzed later.

Figure 4 Temperature dependence on heat capacity (and rotation rate) when analyzing the daily mean temperature at 45^∘ N using Eq. (9).

Quite often, longwave radiation ϵσT⁴ is linearized in energy balance models. Indeed the linearization is performed around 0 ^∘C (North, 1975a, b; Chen et al., 1995; Lohmann and Gerdes, 1998; North and Kim, 2017) and is formulated as

$$\begin{array}{}\text{(14)}& A+B\cdot T.\end{array}$$

As the temperatures based on the local energy balance without a heat capacity would vary between T_min=0 K and $T_{\max }=\sqrt{\mathrm{2}}\cdot \mathrm{288}$ K =407 K, a linearization would be not permitted. Fig. A2 shows that a linear fit for the full range of these temperatures would provide high residuals (blue line) compared to the fit of observational ranges in temperature between −50 to 30 ^∘C (thin red line). The thick red line takes Budyko's (1969) empirical values for present climate conditions: A=203.3 W m⁻² and B=2.09 W m^−2 ∘C⁻¹. Therefore, the linearization implicitly assumes the above heat capacity and fast rotation arguments. If we assume a linearization, we can repeat the calculation Eq. (5) to get

$$\begin{array}{}\text{(15)}& \stackrel{\mathrm{‾}}{T}={\displaystyle \frac{(\mathrm{1}-\mathit{\alpha })S/\mathrm{4}-A}{B}},\end{array}$$

with $\stackrel{\mathrm{‾}}{T}=\mathrm{288}$ K taking Budyko's (1969) values for A and B.

2.3 Meridional heat transport

Equation (11) is the starting point for further investigation. One can easily include the meridional heat transport by diffusion and has been previously used in one-dimensional EBMs (e.g., Adem, 1965; Sellers, 1973; Budyko, 1969; North, 1975a, b). In the following, we will drop the tilde sign. Using a diffusion coefficient k, the meridional heat transport across a latitude is HT$=-k\mathrm{\nabla }T$. One can solve the EBM numerically.

$$\begin{array}{}\text{(16)}& C_{\mathrm{p}}{\partial }_{t}T=-\mathrm{\nabla }\cdot HT+(\mathrm{1}-\mathit{\alpha }){\displaystyle \frac{S}{\mathit{\pi }}}\cos \mathit{\phi }-\mathit{ϵ}\mathit{\sigma }{T}^{\mathrm{4}}\end{array}$$

The boundary condition is that the HT at the poles vanish. The values of k are in the range of earlier studies (North, 1975a, b; Stocker et al., 1992; Chen et al., 1995; Lohmann et al., 1996). Figure 6 shows the equilibrium solutions of Eq. (16) using different values of k (solid lines). The global mean temperature is not affected by the transport term because it depends only of global net radiative fluxes, not internal redistribution. Formally, the integration with boundary condition with zero heat transport at the poles provides no effect (note that ${\partial }_{y}T=\mathrm{0}$ at the North and South Pole). The same is true if we introduce zonal transports because of the cyclic boundary condition in θ direction.

Figure 5 Temperature depending on C_p when solving (Eq. 9) numerically. The reference heat capacity is the atmospheric heat capacity $C_{\mathrm{p}}^{\mathrm{a}}=\mathrm{1.02}\times {\mathrm{10}}^{\mathrm{7}}$ J K⁻¹ m⁻². The climate is insensitive to changes in heat capacity $C_{\mathrm{p}}\in [\mathrm{0.05}\cdot C_{\mathrm{p}}^{\mathrm{a}}$, $\mathrm{2}\cdot C_{\mathrm{p}}^{\mathrm{a}}]$.

Figure 6 Equilibrium temperature of Eq. (16) using different diffusion coefficients ($C_{\mathrm{p}}=C_{\mathrm{p}}^{\mathrm{a}}$) . The blue lines use 1.5×10⁶ m² s⁻¹ with no tilt (solid line), a tilt of 23.5^∘ (dotted line), and a tilt of 23.5^∘ with ice–albedo feedback using the representation of Sellers (1969) (dashed line). Except for the dashed line, the global mean values are identical to the value calculated in Eq. (13). Values are given in ^∘C.

Until now, we assumed that the Earth's axis of rotation was vertical with respect to the path of its orbit around the Sun. Instead Earth's axis is tilted off the vertical by about 23.5^∘. As the Earth orbits the Sun, the tilt causes one hemisphere to receive more direct sunlight and to have longer days. This is a redistribution of heat with more solar insolation at the poles and less at the Equator (formally it could be associated to an enhanced meridional heat transport HT). The resulting temperature is shown as the dotted blue line in Fig. 6. A spatially constant temperature in Eq. (1) can be formally seen as a system with infinite diffusion coefficient k→∞ (black line in Fig. 6).

The global mean temperatures are not affected by the tilt and the values are identical to the one calculated in Eq. (13). This is true even if we calculate the seasonal cycle (Berger and Loutre, 1991, 1997; Laepple and Lohmann, 2009). However, if we include nonlinearities such as the ice–albedo feedback (α as a function of T), the global mean value is changing (Budyko, 1969; Sellers, 1969; North, 1975a, b), e.g., the dashed blue line in Fig. 6. Such models can be improved by including an explicit spatial pattern with a seasonal cycle to study the long-term effects of climate on external forcing (Adem, 1981; North et al., 1983) or by adding noise mimicking the effect of short-term features on the long-term climate (Hasselmann, 1976; Lemke, 1977; Lohmann, 2018). A spatial explicit model would include the land–sea mask, heat capacities over land and the ocean, etc. This is, however, beyond the scope of the present paper. In the following we will analyze the effect of the heat capacity onto the seasonal cycle, and its effect on the global temperature.

Figure 7 Annual mean temperature depending on C_p when solving the seasonal resolved EBM (Eq. 17) numerically. For all solutions, we use $k=\mathrm{1.5}\times {\mathrm{10}}^{\mathrm{6}}$ m² s⁻¹¹, present-day orbital parameters and the ice–albedo feedback using the representation of Sellers (1969).

Figure 8 Seasonal amplitude of temperature for the two extreme scenarios in Fig. 7, indicating that the lower seasonality shown by the dashed black line relative to the blue line is linked to a warmer annual mean climate.

2.4 Effect of the seasonal cycle

As a logical next step, let us now include an explicit seasonal cycle in the EBM:

$$\begin{array}{}\text{(17)}& C_{\mathrm{p}}{\partial }_{t}T=-\mathrm{\nabla }\cdot \mathrm{HT}+(\mathrm{1}-\mathit{\alpha })S(\mathit{\phi },t)-\mathit{ϵ}\mathit{\sigma }{T}^{\mathrm{4}},\end{array}$$

with S(φ,t) being calculated daily (Berger and Loutre, 1991; 1997). Eq. (17) is calculated numerically for fixed diffusion coefficient $k=\mathrm{1.5}\times {\mathrm{10}}^{\mathrm{6}}$ m² s⁻¹ under present orbital conditions. Figure 7 indicates that the temperature gradient becomes flatter for large heat capacities. Furthermore, the mean temperature is affected by the choice of C_p. In the case of large heat capacity at high latitudes (for latitudes polewards of φ=50^∘, mimicking large mixed layer depths) and moderate elsewhere, we observe strong warming at high latitudes with moderate warming at low latitudes (dashed curve). This again indicates that we cannot neglect the time-dependent left-hand side in the energy balance equations.

Figure 9 Longwave radiation change with seasonal amplitude of temperature ΔT for different temperatures T.

Figure 10 Anomalous near-surface temperature for the vertical mixing experiment relative to the control climate. (a) Mean over boreal winter and austral summer (DJF). (b) Mean over austral winter and boreal summer (JJA). Shown is the 100-year mean after 900 years of integration using the Earth system model COSMOS. Values are given in ^∘C.

The question remains why the mean temperature in the dashed curve is much higher than the blue curve. Figure 8 shows the seasonal amplitude for the C_p scenarios as indicated by the blue and dashed black line. We see the large variation in the seasonal cycle $\mathrm{\Delta }T=T_{\mathrm{summer}}-T_{\mathrm{winter}}$ for the blue line in Fig. 8 as compared to the dashed line. A reduced seasonal cycle is responsible for significant warming due to ϵσ⋅T⁴. To estimate the order of magnitude, we can assume that the annual mean change in the net longwave radiation can be approximated by the mean of summer and winter values

$$\begin{array}{}\text{(18)}& \mathit{ϵ}\mathit{\sigma }\cdot \stackrel{\mathrm{‾}}{T^{\mathrm{4}}}\approx \mathit{ϵ}\mathit{\sigma }\cdot \mathrm{0.5}\cdot (T_{\mathrm{summer}}^{\mathrm{4}}+T_{\mathrm{winter}}^{\mathrm{4}}).\end{array}$$

If the seasonal cycle is damped or is zero in an extreme case, the net longwave radiation can be approximated by

$$\begin{array}{}\text{(19)}& \mathit{ϵ}\mathit{\sigma }\cdot {\stackrel{\mathrm{‾}}{T}}^{\mathrm{4}}\approx \mathit{ϵ}\mathit{\sigma }\cdot (\mathrm{0.5}\cdot [T_{\mathrm{summer}}+T_{\mathrm{winter}}]{)}^{\mathrm{4}}.\end{array}$$

A lower seasonal cycle provides for a lower net outgoing radiation, as shown by Eqs. (19) < (18).

Figure 9 shows the change in net outgoing longwave radiation with seasonal amplitude of temperature ΔT for different annual mean temperatures T. The temperature dependence on longwave radiation is relatively minor compared to ΔT. In the simulation with the larger seasonal contrast, the outgoing longwave radiation is up to 10 W m⁻² higher, yielding a colder climate (compare the blue and dashed black lines in Fig. 7). It is noted that this feature is missing in the linearized version (Eq. 14) of the outgoing radiation. Therefore, the change in seasonal amplitude does not directly influence the mean temperature in the linear version. In the following, we will analyze the change in seasonality and mean climate in a complex model.

2.5 Seasonal temperature changes in a complex model

A complex circulation model is used where the seasonal cycle is reduced by enhanced vertical mixing in the ocean. To make a rough estimate of the involved mixed layer, one can see that the effective heat capacity of the ocean is timescale dependent. A diffusive heat flux goes down the gradient of temperature and the convergence of this heat flux drives a ocean temperature tendency:

$$\begin{array}{}\text{(20)}& C_{\mathrm{p}}^{\mathrm{o}}{\partial }_{t}T=-{\partial }_z\left(k^{\mathrm{o}}{\partial }_{z}T\right),\end{array}$$

where $k_{\mathrm{v}}=k^{\mathrm{o}}/C_{\mathrm{p}}^{\mathrm{o}}$ is the oceanic vertical eddy diffusivity (in m² s⁻¹), and $C_{\mathrm{p}}^{\mathrm{o}}$ is the oceanic heat capacity relevant on the specific timescale. The vertical eddy diffusivity k_v can be estimated from climatological hydrographic data (Olbers et al., 1985; Munk and Wunsch, 1998) and varies roughly between 10⁻⁵ and 10⁻⁴ m² s⁻¹ depending on depth and region.

A scale analysis of Eq. (20) yields a characteristic depth scale h_T through

$$\begin{array}{}\text{(21)}& {\displaystyle \frac{\mathrm{\Delta }T}{\mathrm{\Delta }t}}=k_{\mathrm{v}}{\displaystyle \frac{\mathrm{\Delta }T}{h_T^{\mathrm{2}}}}⟶h_T=\sqrt{k_{\mathrm{v}}\mathrm{\Delta }t}.\end{array}$$

For the diurnal cycle, h_T is less that half a meter and the heat capacity generally less than that of the atmosphere. The seasonal mixed layer depth can be several hundred meters (e.g., de Boyer Montégut et al., 2004). As pointed out by Schwartz (2007), the effective heat capacity only reflects that portion of the global heat capacity that is coupled to the perturbation on the timescale of the perturbation. As an example, the effective heat capacity is subject to change in heat content in the context of global climate change induced by changes in gaseous components of the atmosphere on decadal to centennial timescales.

In order to test the effective heat capacity and mixing hypothesis, we employ the coupled climate model COSMOS, which was mainly developed at the Max Planck Institute for Meteorology in Hamburg (Jungclaus et al., 2010). The model contains explicit diurnal and seasonal cycles, it has no flux correction and has been successfully applied to test a variety of paleoclimate hypotheses, ranging from the Miocene climate (Knorr et al., 2011; Knorr and Lohmann, 2014; Stein et al., 2016), the Pliocene climate (Stepanek and Lohmann, 2012), and glacial (Zhang et al., 2013, 2014) and interglacial climates (Wei and Lohmann, 2012; Lohmann et al., 2013; Pfeiffer and Lohmann, 2016).

In order to mimic the effect of a higher effective heat capacity and deepened mixed layer depth, the vertical mixing coefficient is increased in the ocean, changing the values for the background vertical diffusivity (arbitrarily) by a factor of 25, providing a deeper thermocline. The mixing has a background value plus a mixing process strongly influenced by the shears of the mean currents. Although observations give a range of values of k_v for the ocean interior, models use simplified physics and prescribe a constant background value. The model uses a classical vertical eddy viscosity and diffusion scheme (Pacanowski and Philander, 1981). Orbital parameters are fixed to the present condition.

Figure 10 shows the anomalous near-surface temperature for the new vertical mixing experiment relative to the control climate (Wei and Lohmann, 2012). Both simulations were run over 1000 years of integration in order to receive a quasi-equilibrium at the surface. The differences are related to the last 100 years of the simulation. In the vertical mixing experiment, k_v was enhanced leading to more heat is taken up by the ocean producing equable climates with pronounced warming at polar latitudes (Fig. 10). Heat gained at the surface is diffused down the water column, and compared to the control simulation, the wind-induced Ekman cells in the upper part of the oceans intensified and deepened. The model indicates that the respective winter signal of high-latitude warming is most pronounced (Fig. 10), decreasing the seasonality, suggesting a common signal of pronounced warming and weaker seasonality, a feature already seen in the EBM (Figs. 8, 9).

3 Discussion

We can discuss the assumptions for obtaining global mean surface temperature of the Earth. We get realistic temperatures for the following cases.

• Equation (3), where we treat the temperature on the Earth as one number. This approach is written in most textbooks on climate. The averaging is problematic since T⁴ has a pronounced latitudinal dependence. Furthermore, the solution does not take the heat capacity and the Earth's rotation rate into account.

• Equation (13), where we assume a significant heat capacity and fast rotation of the Earth. The diurnal cycle is averaged out. The global mean temperature is similar to a Eq. (3), but multiplied with a factor $G=\sqrt{\frac{\mathit{\pi }}{\mathrm{2}}}\frac{\mathrm{\Gamma }(\mathrm{9}/\mathrm{8})}{\mathrm{\Gamma }(\mathrm{13}/\mathrm{8})}\approx \mathrm{0.989}$

• One can obtain the numerical solution of Eq. (9), where the diurnal cycle has been explicitly taken into account. The global mean temperature is similar to the previous solution and close to observations.

• If we linearize the outgoing radiation Eq. (14), which implicitly assumes the above heat capacity and fast rotation arguments, we get a realistic $\stackrel{\mathrm{‾}}{T}=\mathrm{288}$ K for present climate conditions Eq. (15) .

• The global mean temperatures in slightly more complex EBMs Eqs. (16) and (17) are in the range of observations for significant effective heat capacities.

• The global mean temperature is realistically simulated in a complex GCM with diurnal and seasonal cycle, as used in, e.g., Sect. 2.5 (Stepanek and Lohmann, 2012; Wei and Lohmann, 2012).

In the EBM it is found that a low effective heat capacity enhances the daily and seasonal cycle providing higher outgoing longwave radiation and cooling the Earth. In cases of zero heat capacity, the mean temperature is a factor $F=\mathrm{0.4}\sqrt{\mathrm{2}}\approx \mathrm{0.566}$ lower than 288 K, as stated in Eq. (3), and would be $\stackrel{\mathrm{‾}}{T}=\mathrm{163}$ K.

We can be put our findings into a general statement. Let us define $\stackrel{\mathrm{‾}}{\cdot }$ here as the averaging over an arbitrary time period (in our context the seasonal or diurnal cycle); thus, $\stackrel{\mathrm{‾}}{T^{\mathrm{4}}}>{\stackrel{\mathrm{‾}}{T}}^{\mathrm{4}}$, which is consistent with Hölder's inequality (Rodgers, 1888; Hölder 1889; Hardy et al., 1934, Kuptsov, 2001). For the diurnal cycle, $\stackrel{\mathrm{‾}}{T^{\mathrm{4}}}$ in Eq. (5) is greater than ${\stackrel{\mathrm{‾}}{T}}^{\mathrm{4}}$ as obtained from Eq. (8). We find that $\stackrel{\mathrm{‾}}{T}$ is a factor $F=\mathrm{0.4}\sqrt{\mathrm{2}}$ lower than the fourth root of $\stackrel{\mathrm{‾}}{T^{\mathrm{4}}}$. For the seasonal cycle, the climate with the lower seasonal cycle (dashed black line in Figs. 7 and 8) is warmer than state with a more pronounced seasonal cycle (blue line in Figs. 7 and 8). This is due to the strong dependence of the outgoing radiation on the seasonal range of temperatures (Fig. 9). The higher the daily/seasonal contrast, the more outgoing longwave radiation is increased. We examine whether this is related to the effective heat capacity of the system.

The effective heat capacity is connected to the perturbation on the timescale of the perturbation (Eq. 21). Previous studies have noted that changing the ocean mixed layer depth impacts the climatological annual mean temperature (Schneider and Zhu, 1998; Qiao et al., 2004; Donohoe et al., 2014; Wang et al., 2019). The increased heat capacity of the mixed layer reduced the magnitude of the annual cycle affecting the surface winds and upwelling, which may provide nonlinear effects (Wang et al., 2019). For the past, a strong warming at high latitudes is reconstructed for the Pliocene, Miocene, and Eocene periods (Markwick, 1994; Wolfe, 1994; Sloan and Rea, 1996; Huber et al., 2000; Shellito et al., 2003; Tripati et al., 2003; Mosbrugger et al., 2005; Utescher and Mosbrugger, 2007). It is a conundrum that the modeled high latitudes are not as warm as the reconstructions (e.g., Sloan and Rea, 1996; Huber et al., 2000; Mosbrugger et al., 2005; Knorr et al., 2011; Dowset et al., 2013). Inspired by the EBM and GCM results, we may think of a climate system as having a higher effective heat capacity producing a reduced seasonal cycle and flat temperature gradients. The changed vertical mixing coefficients could mimic possible effects like weak tidal dissipation or abyssal stratification (e.g., Lambeck 1977; Green and Huber, 2013), but its explicit physics are not evaluated here.

4 Conclusions

This paper revisits the relationship between the (global mean) surface temperature of the Earth and its radiation budget, as is frequently used in energy balance models (EBMs). The main point is that the effective heat capacity and its temporal variation over the daily/seasonal cycle needs to be taken into account when estimating surface temperature from the energy budget. EBMs provide a crucial tool in climate research, especially because they – confirmed by the results of the elaborate realistic climate models – describe the processes essential for the genesis of the global climate. EBMs are thus an admissible conceptual tool, due to their reduced complexity in relation to the essentials of “scientific understanding” represents (von Storch et al., 1999). This understanding states that the radiation balance on the ground and the absorption in the atmosphere are the essential factors for determining the temperature. The solution of the basic EBM says that the temperature is independent of the size of the Earth and the thermal characteristics but depends on the albedo, emissivity and solar constant.

The argument follows the conservation of energy: at steady state the Earth has to emit as much energy as it receives from the Sun. However, I argue that we shall explicitly emphasize the Earth as a rapidly rotating object with a significant heat capacity. Without these effects, the global mean temperature would be much lower and can be better used for objects like the Moon or Mercury (Vasavada et al., 1999; Lorenz, 2005), as they are slowly rotating bodies without significant heat capacity. The Earth system understanding states that the effect of the effective heat capacity is important for the radiation balance, other processes – like horizontal transport processes or the ice–albedo feedback – are only of secondary importance for the globally averaged temperature. Interestingly, the global mean temperature in the revised EBM is very close to the original proposed value.

As a basic feature, the net outgoing longwave radiation is reduced if the diurnal or seasonal cycle is reduced, providing for a significant warming. On long timescales, the effective heat capacity is linked to the mixed-layer depth of the ocean. A change in the mixed layer depth, which likely happened through glacial–interglacial cycles (e.g., Zhang et al., 2014), can therefore be an important driver constraining climate sensitivity (Köhler, et al., 2010). This could also be relevant for past and future climate change when the ocean stratification and mixing can change. The temperature dependence is indeed emphasized in a sensitivity study of climatological sea surface temperature to slab ocean model thickness (Wang et al., 2019). It might be that the more effective mixing provides an explanation for why high latitudes were much warmer than present and more equitable, in that the summer-to-winter range of temperature was much reduced (Sloan and Barron, 1990, Valdes et al., 1996; Sloan et al., 2001; Spicer et al. 2004). Interestingly, it has been suggested that the tight link between ocean temperature and CO₂ formed only during the Pliocene when the thermocline shoals and surface water became more sensitive to CO₂ (La Riviere et al., 2012), which is therefore of major importance for the understanding of the climate–carbon cycle (Wiebe and Weaver, 1999; Zachos et al., 2008; de Boer and Hogg, 2014).

It is concluded that climate studies should use improved representations of vertical mixing processes, including turbulence, tidal mixing, hurricanes, and wave breaking (e.g., Qiao et al., 2004; Huber et al., 2004; Simmons et al., 2004; Korty et al., 2008; Griffiths and Peltier, 2009; Green and Huber, 2013; Reichl and Hallberg, 2018). Global climate models treat ocean vertical mixing as static, although there is little reason to suspect this is correct (e.g., see Munk and Wunsch, 1998). In numerical modeling, the values are also constrained by the required numerical stability and fill gaps left by other parameterizations (e.g., Griffies, 2005). As a natural next step, one should analyze the ocean mixing and heat uptake (Luyten et al., 1983; Large et al., 1994; Nilsson, 1995) to understand past, present and future temperatures.

Appendix A

Figure A1 Longwave radiation ϵσT⁴ based on hourly near-surface temperatures, generated using Copernicus Climate Change Service information for 1979 CE following Hersbach et al. (2020). (a) It is seen that $\mathit{ϵ}\mathit{\sigma }\frac{\mathrm{1}}{\mathrm{2}\mathit{\pi }}\underset{\mathrm{0}}{\overset{\mathrm{2}\mathit{\pi }}{\int }}{T}^{\mathrm{4}}d\mathrm{\Theta }\approx \mathit{ϵ}\mathit{\sigma }{\stackrel{\mathrm{̃}}{T}}^{\mathrm{4}}$. (b) The differences are more than 2 orders of magnitude smaller than the mean values. The most pronounced differences are over the subtropics (deserts) and mid-to-high latitudes, where we have a pronounced daily cycle.

Figure A2 Black dots show ϵσT⁴. The blue line shows the linear fit for the range of temperatures 0 to 407 K. The thin red line shows the linear fit for temperatures between −50 to 30 ^∘C (range is shown as vertical dotted lines). The thick red line shows Budyko's (1969) linearization with A=203.3 W m⁻² and B=2.09 W m^−2 ∘C⁻¹. The regression coefficients for the blue line are A=321.4 W m⁻² and B=1.88 W m^−2 ∘C⁻¹, whereas for the thin red line they are A=199.5 W m⁻² and B=2.56 W m^−2 ∘C⁻¹.