Brian E. J. Rose, University at Albany
You really should be looking at The Climate Laboratory book by Brian Rose, where all the same content (and more!) is kept up to date.
Here you are likely to find broken links and broken code.
This document uses the interactive Jupyter notebook format. The notes can be accessed in several different ways:
github at https://github.com/brian-rose/ClimateModeling_coursewareAlso here is a legacy version from 2015.
Many of these notes make use of the climlab package, available at https://github.com/brian-rose/climlab
In [1]:
# Ensure compatibility with Python 2 and 3
from __future__ import print_function, division
Let's look again at the observations:
Last class we introduced a very simple model for the OLR or Outgoing Longwave Radiation to space:
$$ \text{OLR} = \tau \sigma T_s^4 $$where $\tau$ is the transmissivity of the atmosphere, a number less than 1 that represents the greenhouse effect of Earth's atmosphere.
We also tuned this model to the observations by choosing $ \tau \approx 0.61$.
More precisely:
In [ ]:
OLRobserved = 238.5 # in W/m2
sigma = 5.67E-8 # S-B constant
Tsobserved = 288. # global average surface temperature
tau = OLRobserved / sigma / Tsobserved**4 # solve for tuned value of transmissivity
print(tau)
Let's now deal with the shortwave (solar) side of the energy budget.
In [ ]:
Q = 341.3 # the insolation
In [ ]:
Freflected = 101.9 # reflected shortwave flux in W/m2
The planetary albedo is the fraction of $Q$ that is reflected.
We will denote the planetary albedo by $\alpha$.
From the observations:
In [ ]:
alpha = Freflected / Q
print(alpha)
That is, about 30% of the incoming radiation is reflected back to space.
From the observations:
In [ ]:
ASRobserved = Q - Freflected
print(ASRobserved)
As we noted last time, this number is just slightly greater than the observed OLR of 238.5 W m$^{-2}$.
This is one of the central concepts in climate modeling.
The Earth system is in energy balance when energy in = energy out, i.e. when
$$ \text{ASR} = \text{OLR} $$We want to know:
With our simple greenhouse model, we can get an exact solution for the equilibrium temperature.
First, write down our statement of energy balance:
$$ (1-\alpha) Q = \tau \sigma T_s^4 $$Rearrange to solve for $T_s$:
$$ T_s^4 = \frac{(1-\alpha) Q}{\tau \sigma} $$and take the fourth root, denoting our equilibrium temperature as $T_{eq}$:
$$ T_{eq} = \left( \frac{(1-\alpha) Q}{\tau \sigma} \right)^\frac{1}{4} $$Plugging the observed values back in, we compute:
In [ ]:
# define a reusable function!
def equilibrium_temperature(alpha,Q,tau):
return ((1-alpha)*Q/(tau*sigma))**(1/4)
Teq_observed = equilibrium_temperature(alpha,Q,tau)
print(Teq_observed)
And this equilibrium temperature is just slightly warmer than 288 K. Why?
For this very simple model, we can work out the answer exactly:
In [ ]:
Teq_new = equilibrium_temperature(0.32,Q,0.57)
# an example of formatted print output, limiting to two or one decimal places
print('The new equilibrium temperature is {:.2f} K.'.format(Teq_new))
print('The equilibrium temperature increased by about {:.1f} K.'.format(Teq_new-Teq_observed))
Most climate models are more complicated mathematically, and solving directly for the equilibrium temperature will not be possible!
Instead, we will be able to use the model to calculate the terms in the energy budget (ASR and OLR).
In [ ]:
The above exercise shows us that if some properties of the climate system change in such a way that the equilibrium temperature goes up, then the Earth system receives more energy from the sun than it is losing to space. The system is no longer in energy balance.
The temperature must then increase to get back into balance. The increase will not happen all at once! It will take time for energy to accumulate in the climate system. We want to model this time-dependent adjustment of the system.
In fact almost all climate models are time-dependent, meaning the model calculates time derivatives (rates of change) of climate variables.
We will write the total energy budget of the Earth system as
$$ \frac{dE}{dt} = (1-\alpha) Q - OLR $$Note: This is a generically true statement. We have just defined some terms, and made the (very good) assumption that the only significant energy sources are radiative exchanges with space.
This equation is the starting point for EVERY CLIMATE MODEL.
But so far, we don’t actually have a MODEL. We just have a statement of a budget. To use this budget to make a model, we need to relate terms in the budget to state variables of the atmosphere-ocean system.
For now, the state variable we are most interested in is temperature – because it is directly connected to the physics of each term above.
If we now suppose that
$$ E = C T_s $$where $T_s$ is the global mean surface temperature, and $C$ is a constant – the effective heat capacity of the atmosphere- ocean column.
then our budget equation becomes:
where
By adopting this equation, we are assuming that the energy content of the Earth system (atmosphere, ocean, ice, etc.) is proportional to surface temperature.
Important things to think about:
For our purposes here we are going to use a value of C equivalent to heating 100 meters of water:
$$C = c_w \rho_w H$$where
$c_w = 4 \times 10^3$ J kg$^{-1}$ $^\circ$C$^{-1}$ is the specific heat of water,
$\rho_w = 10^3$ kg m$^{-3}$ is the density of water, and
$H$ is an effective depth of water that is heated or cooled.
In [ ]:
c_w = 4E3 # Specific heat of water in J/kg/K
rho_w = 1E3 # Density of water in kg/m3
H = 100. # Depth of water in m
C = c_w * rho_w * H # Heat capacity of the model
print('The effective heat capacity is {:.1e} J/m2/K'.format(C))
This is a first-order Ordinary Differential Equation (ODE) for $T_s$ as a function of time. It is also our very first climate model.
To solve it (i.e. see how $T_s$ evolves from some specified initial condition) we have two choices:
Option 1 (analytical) will usually not be possible because the equations will typically be too complex and non-linear. This is why computers are our best friends in the world of climate modeling.
HOWEVER it is often useful and instructive to simplify a model down to something that is analytically solvable when possible. Why? Two reasons:
Recall that the derivative is the instantaneous rate of change. It is defined as
$$ \frac{dT}{dt} = \lim_{\Delta t\rightarrow 0} \frac{\Delta T}{\Delta t}$$So we write our model as
$$ C \frac{\Delta T}{\Delta t} \approx \text{ASR} - \text{OLR}$$where $\Delta T$ is the change in temperature predicted by our model over a short time interval $\Delta t$.
We can now use this to make a prediction:
Given a current temperature $T_1$ at time $t_1$, what is the temperature $T_2$ at a future time $t_2$?
We can write
$$ \Delta T = T_2-T_1 $$$$ \Delta t = t_2-t_1 $$and so our model says
$$ C \frac{T_2-T_1}{\Delta t} = \text{ASR} - \text{OLR} $$Which we can rearrange to solve for the future temperature:
$$ T_2 = T_1 + \frac{\Delta t}{C} \left( \text{ASR} - \text{OLR}(T_1) \right) $$We now have a formula with which to make our prediction!
Notice that we have written the OLR as a function of temperature. We will use the current temperature $T_1$ to compute the OLR, and use that OLR to determine the future temperature.
The quantity $\Delta t$ is called a timestep. It is the smallest time interval represented in our model.
Here we're going to use a timestep of 1 year:
In [ ]:
dt = 60. * 60. * 24. * 365. # one year expressed in seconds
In [ ]:
# Try a single timestep, assuming we have working functions for ASR and OLR
T1 = 288.
T2 = T1 + dt / C * ( ASR(alpha=0.32) - OLR(T1, tau=0.57) )
print(T2)
What happened? Why?
Try another timestep
In [ ]:
T1 = T2
T2 = T1 + dt / C * ( ASR(alpha=0.32) - OLR(T1, tau=0.57) )
print(T2)
Warmed up again, but by a smaller amount.
But this is tedious typing. Time to define a function to make things easier and more reliable:
In [ ]:
def step_forward(T):
return T + dt / C * ( ASR(alpha=0.32) - OLR(T, tau=0.57) )
Try it out with an arbitrary temperature:
In [ ]:
step_forward(300.)
Notice that our function calls other functions and variables we have already defined.
Now let's really harness the power of the computer by making a loop (and storing values in arrays):
In [ ]:
import numpy as np
numsteps = 20
Tsteps = np.zeros(numsteps+1)
Years = np.zeros(numsteps+1)
Tsteps[0] = 288.
for n in range(numsteps):
Years[n+1] = n+1
Tsteps[n+1] = step_forward( Tsteps[n] )
print(Tsteps)
What did we just do?
In [ ]:
# a special instruction for the Jupyter notebook
# Display all plots inline in the notebook
%matplotlib inline
# import the plotting package
import matplotlib.pyplot as plt
In [ ]:
plt.plot( Years, Tsteps )
plt.xlabel('Years')
plt.ylabel('Global mean temperature (K)');
Note how the temperature adjusts smoothly toward the equilibrium temperature, that is, the temperature at which ASR = OLR.
If the planetary energy budget is out of balance, the temperature must change so that the OLR gets closer to the ASR!
The adjustment is actually an exponential decay process: The rate of adjustment slows as the temperature approaches equilibrium.
The temperature gets very very close to equilibrium but never reaches it exactly.
plt.plot(x,y), where x and y are arrays of the same size. But we must import it first.This is actually not native Python, but uses a special graphics library called matplotlib.
Just about all of our notebooks will start with this:
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
In [ ]:
We are going to linearize the equation for small perturbations away from this equilibrium.
Let $T_s = T_{eq} + T_s^\prime$ and restrict our solution to $T_s^\prime << T_{eq}$.
Note this this is not a big restriction! For example, a 10 degree warming or cooling is just $\pm$3.4% of the absolute equilibrium temperature.
and the budget for the perturbation temperature thus becomes
$$C \frac{d T_s^\prime}{d t} = -\lambda_0 T_s^\prime$$where we define
$$\lambda_0 = 4 \tau \sigma T_{eq}^3 $$Putting in our observational values, we get
In [ ]:
lambda_0 = 4 * sigma * tau * Teq_observed**3
# This is an example of formatted text output in Python
print( 'lambda_0 = {:.2f} W m-2 K-1'.format(lambda_0) )
This is actually our first estimate of what is often called the Planck feedback. It is the tendency for a warm surface to cool by increased longwave radiation to space.
It may also be refered to as the "no-feedback" climate response parameter. As we will see, $\lambda_0$ quantifies the sensitivity of the climate system in the absence of any actual feedback processes.
Now define
$$ t^* = \frac{C}{\lambda_0} $$This is a positive constant with dimensions of time (seconds). With these definitions the temperature evolves according to
$$ \frac{d T_s^\prime}{d t} = - \frac{T_s^\prime}{t^*}$$This is one of the simplest ODEs. Hopefully it looks familiar to most of you. It is the equation for an exponential decay process.
We can easily solve for the temperature evolution by integrating from an initial condition $T_s^\prime(0)$:
$$ \int_{T_s^\prime(0)}^{T_s^\prime(t)} \frac{d T_s^\prime}{T_s^\prime} = -\int_0^t \frac{dt}{t^*}$$$$\ln \bigg( \frac{T_s^\prime(t)}{T_s^\prime(0)} \bigg) = -\frac{t}{t^*}$$$$T_s^\prime(t) = T_s^\prime(0) \exp \bigg(-\frac{t}{t^*} \bigg)$$I hope that the mathematics is straightforward for everyone in this class. If not, go through it carefully and make sure you understand each step.
Our model says that surface temperature will relax toward its equilibrium value over a characteristic time scale $t^*$. This is an e-folding time – the time it takes for the perturbation to decay by a factor $1/e = 0.37$
What should this timescale be for the climate system?
To estimate $t^*$ we need a value for the effective heat capacity $C$.
Our "quick and dirty" estimate above used 100 meters of water to set this heat capacity.
We will revisit this question later in the course. For now, let’s just continue assuming $H = 100$ m (a bit deeper than the typical depth of the surface mixed layer in the oceans).
Now calculate the e-folding time for the surface temperature:
In [ ]:
tstar = C / lambda_0 # Calculated value of relaxation time constant
seconds_per_year = 60.*60.*24.*365.
print( 'The e-folding time is {:1.2e} seconds or about {:1.0f} years.'.format(tstar, tstar / seconds_per_year))
This is a rather fast timescale relative to other processes that can affect the planetary energy budget.
But notice that the climate feedback parameter $\lambda$ is smaller, the timescale gets longer. We will come back to this later.
In [2]:
%load_ext version_information
%version_information numpy, matplotlib
Out[2]:
The author of this notebook is Brian E. J. Rose, University at Albany.
It was developed in support of ATM 623: Climate Modeling, a graduate-level course in the Department of Atmospheric and Envionmental Sciences
Development of these notes and the climlab software is partially supported by the National Science Foundation under award AGS-1455071 to Brian Rose. Any opinions, findings, conclusions or recommendations expressed here are mine and do not necessarily reflect the views of the National Science Foundation.
In [ ]: