# Modeling the global energy budget

## Introducing the zero-dimensional Energy Balance Model

____________
<a id='section1'></a>

## 1. Recap of the global energy budget
____________

Let's look again at the observations:



____________
<a id='section1'></a>

## 2. Tuning radiative fluxes to the observations
____________

### Recap of our simple greenhouse model

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:


```python
OLR_obs = 238.5  # in W/m2
sigma = 5.67E-8  # S-B constant
Ts_obs = 288.  # global average surface temperature
tau = OLR_obs / sigma / Ts_obs**4  # solve for tuned value of transmissivity
print(tau)
```

    0.6114139923687016


Let's now deal with the shortwave (solar) side of the energy budget.

### Absorbed Shortwave Radiation (ASR) and Planetary Albedo

Let's define a few terms.

#### Global mean insolation

From the observations, the area-averaged incoming solar radiation, or **insolation**, is 341.3 W m$^{-2}$.

Let's denote this quantity by $Q$.


```python
Q = 341.3  # the insolation
```

#### Planetary albedo

Some of the incoming radiation is not absorbed at all but simply reflected back to space. Let's call this quantity $F_{reflected}$

From observations we have:


```python
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:


```python
alpha = Freflected / Q
print(alpha)
```

    0.29856431292118374


That is, about 30% of the incoming radiation is reflected back to space.

#### Absorbed Shortwave Radiation

The **Absorbed Shortwave Radiation** or ASR is the part of the incoming sunlight that is *not* reflected back to space, i.e. that part that is absorbed somewhere within the Earth system.

Mathematically we write

$$ \text{ASR} = Q - F_{reflected}\\
= Q-\alpha Q \\
= (1-\alpha) Q $$

From the observations:


```python
ASRobserved = Q - Freflected
print(ASRobserved)
```

    239.4


As we noted last time, this number is *just slightly greater* than the observed OLR of 238.5 W m$^{-2}$.

____________
<a id='section3'></a>

## 3. Equilibrium temperature
____________


*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:

- What surface temperature do we need to have this balance?
- By how much would the temperature change in response to other changes in Earth system?
    - Changes in greenhouse gases
    - Changes in cloudiness
    - etc.




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:


```python
#  define a reusable function
def equilibrium_temperature(alpha, Q, tau):
    return ((1-alpha) *Q / (tau * sigma))**(1/4)

# call the function, passing arguments, and assign the return value to a new variable
Teq_obs = equilibrium_temperature(alpha, Q, tau)
print(Teq_obs)
```

    288.27131447889224


And this equilibrium temperature is *just slightly warmer* than 288 K.

____________
## 4. A climate change scenario
____________

Suppose that, due to global warming (changes in atmospheric composition and subsequent changes in cloudiness):

- The longwave transmissitivity decreases to $\tau = 0.57$ 
- The planetary albedo increases to $\alpha = 0.32$

What is the ***new equilibrium temperature***?

For this very simple model, we can work out the answer exactly:


```python
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_obs))
```

    The new equilibrium temperature is 291.10 K.
    The equilibrium temperature increased by about 2.8 K.


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).

### Python exercise

- Write **two** Python functions to calculate ASR and OLR for *arbitrary parameter values*.
- Verify the following:
    - With the new parameter values but the old temperature $T = 288$ K, is ASR greater or lesser than OLR? (Use formatted strings to print your results.)
    - Is the Earth gaining or losing energy?
    - How does your answer change if $T = 295$ K (or any other temperature greater than 291 K)?


```python

```


```python

```

____________

## 5. A time-dependent Energy Balance Model
____________

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.

### An energy balance **equation**

We will write the **total energy budget** of the Earth system as

\begin{align} 
\frac{dE}{dt} &= \text{net energy flux in to system} \\
 &= \text{flux in – flux out} \\
 &= \text{ASR} - \text{OLR}
\end{align}

where $E$ is the **enthalpy** or **heat content** of the total system.

We will express the budget per unit surface area, so each term above has units W m$^{-2}$

Note: any **internal exchanges** of energy between different reservoirs (e.g. between ocean, land, ice, atmosphere) do not appear in this budget – because $E$ is the **sum of all reservoirs**.

Also 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.



### An energy balance **model**

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:

$$ C \frac{dT_s}{dt} = \text{ASR} - \text{OLR} $$

where

- $C$ is the **heat capacity** of Earth system, in units of J m$^{-2}$ K$^{-1}$.
- $\frac{dT_s}{dt}$ is the rate of change of global average surface temperature.

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:

- Why is this a sensible assumption?
- What determines the heat capacity $C$?
- What are some limitations of this assumption?


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.


```python
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))
```

    The effective heat capacity is 4.0e+08 J/m2/K


### Solving the energy balance model

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:

1. Solve it analytically
2. Solve it numerically

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:

1. Analysis will often yield a deeper understanding of the behavior of the system
2. It gives us a benchmark against which to test the results of our numerical solutions.

____________

## 6. Representing time derivatives on a computer
____________

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}$$

- **On the computer there is no such thing as an instantaneous change.** 
- We are always dealing with *discrete quantities*.
- So we approximate the derivative with $\Delta T/ \Delta t$. 
- So long as we take the time interval $\Delta t$ "small enough", the approximation is valid and useful.
- (The meaning of "small enough" varies widely in practice. Let's not talk about it now)

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.

____________

## 7. Numerical solution of the Energy Balance Model
____________

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:


```python
dt = 60. * 60. * 24. * 365.   # one year expressed in seconds
```


```python
# 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) ) 
# above I am passing arguments using keywords so that I don't have to remember the order
print(T2)
```

    288.7678026614462


What happened? Why?

Try another timestep


```python
T1 = T2
T2 = T1 + dt / C * ( ASR(alpha=0.32) - OLR(T1, tau=0.57) )
print(T2)
```

    289.3479210238739


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:


```python
def step_forward(T):
    return T + dt / C * ( ASR(alpha=0.32) - OLR(T, tau=0.57) )
```

Try it out with an arbitrary temperature:


```python
step_forward(300.)
```




    297.658459884



Notice that our function calls other functions and variables we have already defined.

***
#### Python tip

Functions can access variables and other functions defined outside of the function. 

This is both very useful and occasionally confusing.
***

Now let's really harness the power of the computer by **making a loop** (and storing values in arrays).

### Python "For" Loops

Definite iteration loops are frequently referred to as **for** loops because for is the keyword that is used to introduce them in nearly all programming languages, including Python.

Definite iteration means that the number of repetitions is specified in advance. Later in the course we will introduce indefinite iteration, in which the code block executes until some condition is met.


```python
# Iterate through a list of numbers and execute statements
for n in [0, 1, 2, 3, 4]:
    print(n)

# Note the loop variable n takes on the value of the next element 
# in the collection each time through the loop.    
```

    0
    1
    2
    3
    4



```python
# Same thing, but use the built-in range() function
# range(<end>) returns an iterable that yields integers 
# starting with 0, up to but not including <end>:
for n in range(5):
    print(n)
```

    0
    1
    2
    3
    4


### Numpy arrays

[NumPy](https://numpy.org/) is the fundamental package for scientific computing with Python.

The fundamental data structure of NumPy is an **n-dimensional array**. N-dimensional means it can be 1-dimensional, 2-dimensional, 3-dimensional, and so on. 

As your programming skills grow, you will want to avoid for loops and intead use array programming to speed up operation runtime. Let's not worry about that now.

To access numpy, we must first import it.


```python
import numpy as np
```

The `linspace` function creates an array object that is an array of numbers evenly spaced between the start and end points.


```python
np.linspace(230,300,10)
```




    array([230.        , 237.77777778, 245.55555556, 253.33333333,
           261.11111111, 268.88888889, 276.66666667, 284.44444444,
           292.22222222, 300.        ])



An array can be a function argument.


```python
OLR(np.linspace(230,300,10), tau=0.57)
```




    array([ 90.44181279, 103.31014479, 117.50516903, 133.11508224,
           150.23091968, 168.9465551 , 189.35870079, 211.56690754,
           235.67356467, 261.7839    ])



The `zeros` function creates an array object that is an array of zeros of the specified size.


```python
np.zeros(10)
```




    array([0., 0., 0., 0., 0., 0., 0., 0., 0., 0.])



### Implementing For-loop and Array assignment in the Energy Balance Model


```python
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
    #Here we are calling a function inside a loop
    Tsteps[n+1] = step_forward( Tsteps[n] )
print(Tsteps)
```

    [288.         288.76780266 289.34792102 289.78523685 290.11433323
     290.36166675 290.54736768 290.68669049 290.79115953 290.86946109
     290.92813114 290.97208122 291.00499865 291.02964965 291.0481083
     291.06192909 291.07227674 291.08002371 291.08582346 291.09016532
     291.09341571]


What did we just do?

- Created an array of zeros
- Set the initial temperature to 288 K
- Repeated our time step 20 times
- Stored the results of each time step into the array

***
#### Python tip

Use square bracket [ ] to refer to elements of an array or list. Use round parentheses ( ) for function arguments. 
***

### Plotting the result

Now let's draw a picture of our result!


```python
# 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
```


```python
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.

***
#### Python tip

We can easily make simple graphs with the function `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 graphics library called [matplotlib](https://matplotlib.org). This is the workhorse of scientific plotting in Python, and we will be using it all the time!

Just about all of our notebooks will start with this:
```
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
```
***

____________
## 8. Summary
____________

- We looked at the flows of energy in and out of the Earth system. 
- These are determined by radiation at the top of the Earth's atmosphere.
- Any imbalance between shortwave absorption (ASR) and longwave emission (OLR) drives a change in temperature
- Using this idea, we built a climate model!
- This **Zero-Dimensional Energy Balance Model** solves for the global, annual mean surface temperature $T_s$
- Two key assumptions:
    - Energy content of the Earth system varies proportionally to $T_s$
    - The OLR increases as $\tau \sigma T_s^4$ (our simple greenhouse model)
- Earth (or any planet) has a well-defined **equilibrium temperature** $T_{eq}$ at which ASR = OLR, because of the *temperature dependence of the outgoing longwave radiation*.

- If $T_s < T_{eq}$, the model will warm up.
- We can represent the continous warming process on the computer using discrete timesteps.
- We can plot the result.

____________

## Credits

This notebook is part of [The Climate Laboratory](https://brian-rose.github.io/ClimateLaboratoryBook), an open-source textbook developed and maintained by [Brian E. J. Rose](http://www.atmos.albany.edu/facstaff/brose/index.html), University at Albany. It has been modified by [Nicole Feldl](http://nicolefeldl.com), UC Santa Cruz.

It is licensed for free and open consumption under the
[Creative Commons Attribution 4.0 International (CC BY 4.0)](https://creativecommons.org/licenses/by/4.0/) license.

Development of these notes and the [climlab software](https://github.com/brian-rose/climlab) 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.
____________


```python

```

1c08e320eb72f3ad3fa082c654609eb768c6edd1
ipynb
Jupyter Notebook
content/courseware/zero-dim-ebm.ipynb
nfeldl/ClimateLaboratoryBook
05eb0395c0e07d3724e6569e160fbefc9829a990
CC-BY-4.0
2021-08-25T13:02:15.000Z
2021-08-25T13:02:15.000Z
content/courseware/zero-dim-ebm.ipynb
nfeldl/ClimateLaboratoryBook
05eb0395c0e07d3724e6569e160fbefc9829a990
CC-BY-4.0
content/courseware/zero-dim-ebm.ipynb
nfeldl/ClimateLaboratoryBook
05eb0395c0e07d3724e6569e160fbefc9829a990
CC-BY-4.0
2021-07-21T20:43:20.000Z
2021-08-25T13:02:16.000Z
Qwen/Qwen-72B
1. YES
2. YES