Brian E. J. Rose, University at Albany
climlab
¶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
# Ensure compatibility with Python 2 and 3
from __future__ import print_function, division
climlab
¶climlab
is a python package for process-oriented climate modeling.
It is based on a very general concept of a model as a collection of individual,
interacting processes. climlab
defines a base class called Process
, which
can contain an arbitrarily complex tree of sub-processes (each also some
sub-class of Process
). Every climate process (radiative, dynamical,
physical, turbulent, convective, chemical, etc.) can be simulated as a stand-alone
process model given appropriate input, or as a sub-process of a more complex model.
New classes of model can easily be defined and run interactively by putting together an
appropriate collection of sub-processes.
climlab
is an open-source community project. The latest code can always be found on github
:
https://github.com/brian-rose/climlab
You can install climlab
by doing
conda install -c conda-forge climlab
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
import climlab
Recall that we have worked with a zero-dimensional Energy Balance Model
$$ C \frac{dT_s}{dt} = (1-\alpha) Q - \tau \sigma T_s^4 $$Here we are going to implement this exact model using climlab
.
Yes, we have already written code to implement this model, but we are going to repeat this effort here as a way of learning how to use climlab
.
There are tools within climlab
to implement much more complicated models, but the basic interface will be the same.
# create a zero-dimensional domain with a single surface temperature
state = climlab.surface_state(num_lat=1, # a single point
water_depth = 100., # 100 meters slab of water (sets the heat capacity)
)
state
AttrDict({'Ts': Field([[32.]])})
Here we have created a dictionary called state
with a single item called Ts
:
state['Ts']
Field([[32.]])
This dictionary holds the state variables for our model -- which is this case is a single number! It is a temperature in degrees Celsius.
For convenience, we can access the same data as an attribute (which lets us use tab-autocomplete when doing interactive work):
state.Ts
Field([[32.]])
It is also possible to see this state
dictionary as an xarray.Dataset
object:
climlab.to_xarray(state)
<xarray.Dataset> Dimensions: (depth: 1, depth_bounds: 2, lat: 1, lat_bounds: 2) Coordinates: * depth (depth) float64 50.0 * lat (lat) float64 0.0 * depth_bounds (depth_bounds) float64 0.0 100.0 * lat_bounds (lat_bounds) float64 -90.0 90.0 Data variables: Ts (depth, lat) float64 32.0
# create the longwave radiation process
olr = climlab.radiation.Boltzmann(name='OutgoingLongwave',
state=state,
tau = 0.612,
eps = 1.,
timestep = 60*60*24*30.)
# Look at what we just created
print(olr)
climlab Process of type <class 'climlab.radiation.boltzmann.Boltzmann'>. State variables and domain shapes: Ts: (1, 1) The subprocess tree: OutgoingLongwave: <class 'climlab.radiation.boltzmann.Boltzmann'>
# create the shortwave radiation process
asr = climlab.radiation.SimpleAbsorbedShortwave(name='AbsorbedShortwave',
state=state,
insolation=341.3,
albedo=0.299,
timestep = 60*60*24*30.)
# Look at what we just created
print(asr)
climlab Process of type <class 'climlab.radiation.absorbed_shorwave.SimpleAbsorbedShortwave'>. State variables and domain shapes: Ts: (1, 1) The subprocess tree: AbsorbedShortwave: <class 'climlab.radiation.absorbed_shorwave.SimpleAbsorbedShortwave'>
# couple them together into a single model
ebm = olr + asr
# Give the parent process name
ebm.name = 'EnergyBalanceModel'
# Examine the model object
print(ebm)
climlab Process of type <class 'climlab.process.time_dependent_process.TimeDependentProcess'>. State variables and domain shapes: Ts: (1, 1) The subprocess tree: EnergyBalanceModel: <class 'climlab.process.time_dependent_process.TimeDependentProcess'> OutgoingLongwave: <class 'climlab.radiation.boltzmann.Boltzmann'> AbsorbedShortwave: <class 'climlab.radiation.absorbed_shorwave.SimpleAbsorbedShortwave'>
The object called ebm
here is the entire model -- including its current state (the temperature Ts
) as well as all the methods needed to integrated forward in time!
The current model state, accessed two ways:
ebm.state
AttrDict({'Ts': Field([[32.]])})
ebm.Ts
Field([[32.]])
Here is some internal information about the timestep of the model:
print(ebm.time['timestep'])
print(ebm.time['steps'])
2592000.0 0
This says the timestep is 2592000 seconds (30 days!), and the model has taken 0 steps forward so far.
To take a single step forward:
ebm.step_forward()
ebm.Ts
Field([[31.61786227]])
The model got colder!
To see why, let's look at some useful diagnostics computed by this model:
ebm.diagnostics
{'OLR': Field([[300.896072]]), 'ASR': 239.25130000000004}
This is another dictionary, now with two items. They should make sense to you.
Just like the state
variables, we can access these diagnostics
variables as attributes:
ebm.OLR
Field([[300.896072]])
ebm.ASR
239.25130000000004
So why did the model get colder in the first timestep?
What do you think will happen next?
Let's look at how the model adjusts toward its equilibrium temperature.
Exercise:
for
loop, take 500 steps forward with this model
Suppose we want to investigate the effects of a small decrease in the transmissitivity of the atmosphere tau
.
Previously we used the zero-dimensional model to investigate a hypothetical climate change scenario in which:
tau
decreases to 0.57How would we do that using climlab
?
Recall that the model is comprised of two sub-components:
for name, process in ebm.subprocess.items():
print(name)
print(process)
OutgoingLongwave climlab Process of type <class 'climlab.radiation.boltzmann.Boltzmann'>. State variables and domain shapes: Ts: (1, 1) The subprocess tree: OutgoingLongwave: <class 'climlab.radiation.boltzmann.Boltzmann'> AbsorbedShortwave climlab Process of type <class 'climlab.radiation.absorbed_shorwave.SimpleAbsorbedShortwave'>. State variables and domain shapes: Ts: (1, 1) The subprocess tree: AbsorbedShortwave: <class 'climlab.radiation.absorbed_shorwave.SimpleAbsorbedShortwave'>
The parameter tau
is a property of the OutgoingLongwave
subprocess:
ebm.subprocess['OutgoingLongwave'].tau
0.612
and the parameter albedo
is a property of the AbsorbedShortwave
subprocess:
ebm.subprocess['AbsorbedShortwave'].albedo
0.299
Let's make an exact clone of our model and then change these two parameters:
ebm2 = climlab.process_like(ebm)
print(ebm2)
climlab Process of type <class 'climlab.process.time_dependent_process.TimeDependentProcess'>. State variables and domain shapes: Ts: (1, 1) The subprocess tree: EnergyBalanceModel: <class 'climlab.process.time_dependent_process.TimeDependentProcess'> OutgoingLongwave: <class 'climlab.radiation.boltzmann.Boltzmann'> AbsorbedShortwave: <class 'climlab.radiation.absorbed_shorwave.SimpleAbsorbedShortwave'>
ebm2.subprocess['OutgoingLongwave'].tau = 0.57
ebm2.subprocess['AbsorbedShortwave'].albedo = 0.32
Now our model is out of equilibrium and the climate will change!
To see this without actually taking a step forward:
# Computes diagnostics based on current state but does not change the state
ebm2.compute_diagnostics()
ebm2.ASR - ebm2.OLR
Field([[-46.76117229]])
Shoud the model warm up or cool down?
Well, we can find out:
ebm2.Ts
Field([[31.61786227]])
ebm2.step_forward()
ebm2.Ts
Field([[31.32798841]])
Often we want to integrate a model forward in time to equilibrium without needing to store information about the transient state.
climlab
offers convenience methods to do this easily:
ebm3 = climlab.process_like(ebm2)
ebm3.integrate_years(50)
Integrating for 608 steps, 18262.11 days, or 50 years. Total elapsed time is 50.10373938170343 years.
# What is the current temperature?
ebm3.Ts
Field([[17.94837835]])
# How close are we to energy balance?
ebm3.ASR - ebm3.OLR
Field([[-0.00021699]])
# We should be able to accomplish the exact same thing with explicit timestepping
for n in range(608):
ebm2.step_forward()
ebm2.Ts
Field([[17.94837835]])
ebm2.ASR - ebm2.OLR
Field([[-0.00021699]])
We will be using climlab
extensively throughout this course. Lots of examples of more advanced usage are found here in the course notes. Here are some links to other resources:
climlab
at the 2018 AMS Python symposium (January 2018)%load_ext version_information
%version_information numpy, matplotlib, climlab
Loading extensions from ~/.ipython/extensions is deprecated. We recommend managing extensions like any other Python packages, in site-packages.
Software | Version |
---|---|
Python | 3.7.3 64bit [Clang 4.0.1 (tags/RELEASE_401/final)] |
IPython | 7.6.0 |
OS | Darwin 17.7.0 x86_64 i386 64bit |
numpy | 1.16.4 |
matplotlib | 3.1.1 |
climlab | 0.7.5 |
Wed Jul 03 14:49:48 2019 EDT |
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. ____________