Sympl: A System for Modelling Planets¶
sympl is an open source project aims to enable researchers and other users to write understandable, modular, accessible Earth system and planetary models in Python. It is meant to be used in combination with other packages that provide model components in order to write model scripts. Its source code can be found on GitHub.
New users should read the Quickstart. This framework is meant to be used along with model toolkits like CliMT to write models. See the paper on Sympl and CliMT for a good overview and some examples! You may also want to check out the CliMT documentation.
Documentation¶
Overview: Why Sympl?¶
Traditional atmospheric and Earth system models can be hard to read and change for many reasons. Sympl tries to learn from the past experience of these models to speed up research and improve accessibility.
At its core, Sympl defines a model in terms of a state that gets changed in sequence by components of a model (like the radiation scheme, or dynamical core). Each of those components as well-defined and documented inputs and outputs, and code in Sympl will automatically handle unit and dimensionality conversions (such as dimension orderings) to give components the inputs they need.
Sympl defines a framework of Python object interfaces (APIs) that can be combined to create a model. This has many benefits:
- Objects can use code written in any language that can be called from Python, including Fortran, C, C++, Julia, Matlab, and others.
- Each object, such as a radiation parameterization, has a clearly documented interface and can be understood without looking at any other part of a model’s code. Certain interfaces have been designed to force model code to self-document, such as having inputs and outputs as properties of a scheme.
- Objects can be swapped out with other compatible objects. For example, Sympl makes it trivial to change the type of time stepping used on a prognostic scheme.
- Code can be re-used between different types of models. For instance, an atmospheric general circulation model, numerical weather prediction model, and large-eddy simulation could all use the same RRTM radiation object.
- Already-existing documentation for Sympl can tell your users how to configure and run your model. You will likely spend less time writing documentation, but end up with a better documented model. As long as you write docstrings, you’re good to go!
Sympl also contains a number of commonly used objects, such as time steppers and NetCDF output objects.
So is Sympl a model?¶
Sympl is not a model itself. In particular, physical parameterizations and dynamical cores are not present in Sympl. This code instead can be found in other projects that make use of Sympl.
Sympl is meant to be a community ecosystem that allows researchers and other users to use and combine components from a number of different sources. By keeping model physics/dynamics code outside of Sympl itself, researchers can own and maintain their own models. The framework API ensures that models using Sympl are clear and accessible, and allows components from different models and packages to be used alongside one another.
Then where’s the model?¶
Check out CliMT as an example. We highly recommend reading the paper on Sympl and CliMT.
In Sympl, the “model” in the traditional sense of Fortran models is essentially your run script. You use a Python run script instead of a Fortran module like main.f90. This may sound scary, but the idea is that the python run script is as easy (or easier) to understand than a configuration file. Sympl makes choices that force those run scripts to be easier to understand. You can see examples in the above paper.
A “model developer” in the traditional sense would write a toolkit package containing model components that are used by those run scripts. See CliMT for an example of this. That toolkit should also come with example run scripts using its components.
In a way, when you configure the model you are writing the model itself. This is reasonable in Sympl because the model run script should be accessible and readable by users with basic knowledge of programming (even users who don’t know Python). By being readable, the model run script tells others clearly and precisely how you configured and ran your model.
If someone wants to write a model in the traditional way, where their Python run script is never changed and instead you configure the model by changing a configuration file, they can do that, too! Read about the particular model you’re using for details.
The API¶
In a Sympl model, the model state is contained within a “state dictionary”. This is a Python dictionary whose keys are strings indicating a quantity, and values are DataArrays with the values of those quantities. The one exception is “time”, which is stored as a timedelta or datetime-like object, not as a DataArray. The DataArrays also contain information about the units of the quantity, and the grid it is located on. At the start of a model script, the state dictionary should be set to initial values. Code to do this may be present in other packages, or you can write this code yourself. The state and its initialization is discussed further in Model State.
The state dictionary is evolved by TendencyStepper
and
Stepper
objects. These types of objects take in the state
and a timedelta object that indicates the time step, and return the next
model state. TendencyStepper
objects do this by wrapping
TendencyComponent
objects, which calculate tendencies using the
state dictionary. We should note that the meaning of “Stepper” in Sympl is
slightly different than its traditional definition. Here an “Stepper” object is
one that calculates the new state directly from the current state, or any
object that requires the timestep to calculate the new state, while
“TendencyComponent” objects are ones that calculate tendencies without using the
timestep. If a TendencyStepper
or Stepper
object needs to use multiple time steps in its calculation, it does so by
storing states it was previously given until they are no longer needed.
The state is also calculated using DiagnosticComponent
objects which
determine diagnostic quantities at the current time from the current state,
returning them in a new dictionary. This type of object is particularly useful
if you want to write your own online diagnostics.
The state can be stored or viewed using Monitor
objects.
These take in the model state and do something with it, such as storing it in
a NetCDF file, or updating an interactive plot that is being shown to the user.
Quickstart¶
Here we have an example of how Sympl might be used to construct a model run script, with explanations of what’s going on. Here is the full model script we will be looking at (we break it into smaller pieces below):
from model_package import (
get_initial_state, Radiation, BoundaryLayer, DeepConvection,
ImplicitDynamics)
from sympl import (
AdamsBashforth, PlotFunctionMonitor, UpdateFrequencyWrapper,
datetime, timedelta)
def my_plot_function(fig, state):
ax = fig.add_subplot(1, 1, 1)
ax.set_xlabel('longitude')
ax.set_ylabel('latitude')
ax.set_title('Lowest model level air temperature (K)')
im = ax.pcolormesh(
state['air_temperature'].to_units('degK').values[0, :, :],
vmin=260.,
vmax=310.)
cbar = fig.colorbar(im)
plot_monitor = PlotFunctionMonitor(my_plot_function)
state = get_initial_state(nx=256, ny=128, nz=64)
state['time'] = datetime(2000, 1, 1)
physics_stepper = AdamsBashforth([
UpdateFrequencyWrapper(Radiation(), timedelta(hours=2)),
BoundaryLayer(),
DeepConvection(),
])
implicit_dynamics = ImplicitDynamics()
timestep = timedelta(minutes=30)
while state['time'] < datetime(2010, 1, 1):
physics_diagnostics, state_after_physics = physics_stepper(state, timestep)
dynamics_diagnostics, next_state = implicit_dynamics(state_after_physics, timestep)
state.update(physics_diagnostics)
state.update(dynamics_diagnostics)
plot_monitor.store(state)
next_state['time'] = state['time'] + timestep
state = next_state
Importing Packages¶
At the beginning of the script we have import statements:
from model_package import (
get_initial_state, Radiation, BoundaryLayer, DeepConvection,
ImplicitDynamics)
from sympl import (
AdamsBashforth, PlotFunctionMonitor, UpdateFrequencyWrapper,
datetime, timedelta)
These grant access to the objects that will be used to construct the model,
and are dependent on the model package you are using. Here, the names
model_package
, get_initial_state
, Radiation
,
BoundaryLayer
, DeepConvection
, and
ImplicitDynamics
are placeholders, and do not refer to
an actual existing package.
Defining a PlotFunctionMonitor¶
Here we define a plotting function, and use it to create a
Monitor
using PlotFunctionMonitor
:
def my_plot_function(fig, state):
ax = fig.add_subplot(1, 1, 1)
ax.set_xlabel('longitude')
ax.set_ylabel('latitude')
ax.set_title('Lowest model level air temperature (K)')
im = ax.pcolormesh(
state['air_temperature'].to_units('degK').values[0, :, :],
vmin=260.,
vmax=310.)
cbar = fig.colorbar(im)
plot_monitor = PlotFunctionMonitor(my_plot_function)
That Monitor
will be used to produce an animated plot of the lowest model
level air temperature as the model runs. Here we assume that the first axis
is the vertical axis, and that the lowest level is at the lowest index, but
this depends entirely on your model. The [0, :, :]
part might be different
for your model.
Initialize the Model State¶
To initialize the model, we need to create a dictionary which contains the model state. The way this is done is model-dependent. Here we assume there is a function that was defined by the model_package package which handles this for us:
state = get_initial_state(nx=256, ny=128, nz=64)
state['time'] = datetime(2000, 1, 1)
An initialized state
is a dictionary whose keys are strings (like
‘air_temperature’) and values are DataArray
objects, which
store not only the data but also metadata like units. The one exception
is the “time” quantity which is either a datetime
-like or timedelta
-like
object. Here we are calling sympl.datetime()
to initialize time,
rather than directly creating a Python datetime. This is because
sympl.datetime()
can support a number of calendars using the
netcdftime package, if installed, unlike the built-in datetime
which only
supports the Proleptic Gregorian calendar.
You can read more about the state
, including sympl.datetime()
in
Model State.
Initialize Components¶
Now we need the objects that will process the state to move it forward in time. Those are the “components”:
physics_stepper = AdamsBashforth([
UpdateFrequencyWrapper(Radiation(), timedelta(hours=2)),
BoundaryLayer(),
DeepConvection(),
])
implicit_dynamics = ImplicitDynamics()
AdamsBashforth
is a TendencyStepper
, which is
created with a set of TendencyComponent
components.
The TendencyComponent
components we have here are Radiation
,
BoundaryLayer
, and DeepConvection
. Each of these carries information about
what it takes as inputs and provides as outputs, and can be called with a model
state to return tendencies for a set of quantities. The
TendencyStepper
uses this information to step the model state
forward in time.
The UpdateFrequencyWrapper
applied to the Radiation
object
is an object that acts like a TendencyComponent
but only computes
its output if at least a certain amount of model time has passed since the last
time the output was computed. Otherwise, it returns the last computed output.
This is commonly used in atmospheric models to avoid doing radiation
calculations (which are very expensive) every timestep, but it can be applied
to any TendencyComponent.
The ImplicitDynamics
class is a Stepper
object, which
steps the model state forward in time in the same way that a TendencyStepper
would, but doesn’t use TendencyComponent
objects in doing so.
The Main Loop¶
With everything initialized, we have the part of the code where the real computation is done – the main loop:
timestep = timedelta(minutes=30)
while state['time'] < datetime(2010, 1, 1):
physics_diagnostics, state_after_physics = physics_stepper(state, timestep)
dynamics_diagnostics, next_state = implicit_dynamics(state_after_physics, timestep)
state.update(physics_diagnostics)
state.update(dynamics_diagnostics)
plot_monitor.store(state)
next_state['time'] = state['time'] + timestep
state = next_state
In the main loop, a series of component calls update the state, and the figure
presented by plot_monitor
is updated. The code is meant to be as
self-explanatory as possible. It is necessary to manually set the time of the
next state at the end of the loop. This is not done automatically by
TendencyStepper
and Stepper
objects, because
in many models you may want to update the state with multiple such objects
in a sequence over the course of a single time step.
Frequently Asked Questions¶
Isn’t Python too slow for Earth System Models?¶
Not in general. Most model run time is spent within code such as the dynamical core, radiation parameterization, and other physics parameterizations. These components can be written in your favorite compiled language like Fortran or C, and then run from within Python. For new projects where you’re writing a component from scratch, we recommend Cython, as it allows you to write typed Python code which gets converted into C code and then compiled. Sympl is designed so that only overhead tasks need to be written in Python.
If 90% of a model’s run time is spent within this computationally intensive, compiled code, and the other 10% is spent in overhead code, then that overhead code taking 3x as long to run would only increase the model’s run time by 1/5th.
But the run time of a model isn’t the only important aspect, you also have to consider time spent programming a model. Poorly designed and documented code can cost weeks of researcher time. It can also take a long time to perform tasks that Sympl makes easy, like porting a component from one model to another. Time is also saved when others have to read and understand your model code.
In short, the more your work involves configuring and developing models, the more time you will save, at the cost of slightly slower model runs. But in the end, what is the cost of your sanity?
What calendar is my model using?¶
Hopefully the section on Choice of Datetime can clear this up.
Installation¶
Latest release¶
To install Sympl, run this command in your terminal:
$ pip install sympl
This is the preferred method to install Sympl, as it will always install the most recent release.
If you don’t have pip installed, this Python installation guide can guide you through the process.
From sources¶
The sources for Sympl can be downloaded from the Github repo.
You can either clone the public repository:
$ git clone git://github.com/mcgibbon/sympl
Or download the tarball:
$ curl -OL https://github.com/mcgibbon/sympl/tarball/master
Once you have a copy of the source, you can install it with:
$ python setup.py install
If you are looking to modify the code, you can install it with:
$ python setup.py develop
This configures the package so that Python points to the current directory instead of copying files. Then when you make modifications to the source code in that directory, they are automatically used by any new Python sessions.
Model State¶
In a Sympl model, physical quantities are stored in a state dictionary. This is
a Python dict with keys that are strings, indicating the quantity name, and
values are DataArray
objects. The DataArray
is a slight modification of the DataArray
object from xarray. It
maintains attributes when it is on the left hand side of addition or
subtraction, and contains a helpful method for converting units. Any
information about the grid the data is using that components need should be
put as attributes in the attrs
of the DataArray
objects. Deciding on
these attributes (if any) is mostly up to the component developers. However,
in order to use the TendencyStepper objects and several helper functions from Sympl,
it is required that a “units” attribute is present.
-
class
sympl.
DataArray
(data, coords=None, dims=None, name=None, attrs=None, encoding=None, fastpath=False)[source]¶ -
__add__
(other)[source]¶ If this DataArray is on the left side of the addition, keep its attributes when adding to the other object.
-
__sub__
(other)[source]¶ If this DataArray is on the left side of the subtraction, keep its attributes when subtracting the other object.
-
to_units
(units)[source]¶ Convert the units of this DataArray, if necessary. No conversion is performed if the units are the same as the units of this DataArray. The units of this DataArray are determined from the “units” attribute in attrs.
Parameters: units (str) – The desired units.
Raises: ValueError
– If the units are invalid for this object.KeyError
– If this object does not have units information in its attrs.
Returns: converted_data – A DataArray containing the data from this object in the desired units, if possible.
Return type:
-
There is one quantity which is not stored as a DataArray
, and
that is “time”. Time must be stored as a datetime or timedelta-like object.
Code to initialize the state is intentionally not present in Sympl, since this depends heavily on the details of the model you are running. You may find helper functions to create an initial state in model packages, or you can write your own. For example, below you can see code to initialize a state with random temperature and pressure on a lat-lon grid (random values are used for demonstration purposes only, and are not recommended in a real model).
from datetime import datetime
import numpy as np
from sympl import DataArray, add_direction_names
n_lat = 64
n_lon = 128
n_height = 32
add_direction_names(x='lat', y='lon', z=('mid_levels', 'interface_levels'))
state = {
"time": datetime(2000, 1, 1),
"air_temperature": DataArray(
np.random.rand(n_lat, n_lon, n_height),
dims=('lat', 'lon', 'mid_levels'),
attrs={'units': 'degK'}),
"air_pressure": DataArray(
np.random.rand(n_lat, n_lon, n_height),
dims=('lat', 'lon', 'mid_levels'),
attrs={'units': 'Pa'}),
"air_pressure_on_interface_levels": DataArray(
np.random.rand(n_lat, n_lon, n_height + 1),
dims=('lat', 'lon', 'interface_levels'),
attrs=('units': 'Pa')),
}
The call to add_direction_names()
tells Sympl
what dimension names correspond to what directions. This information is used
by components to make sure the axes are in the right order.
Choice of Datetime¶
The built-in datetime
object in Python (as used above) assumes the
proleptic Gregorian calendar, which extends the Gregorian calendar back
infinitely. Sympl provides a datetime()
function which returns
a datetime-like object, and allows a variety of different calendars. If a
calendar other than ‘proleptic_gregorian’ is specified, one of the classes from
the netcdftime package will be used. Of course, this requires that it is
installed! If it’s not, you will get an error, and should pip install netcdftime
.
Sympl also includes timedelta
for convenience. This is just
the default Python timedelta
.
To repeat, the calendar your model is using depends entirely on what
object you’re using to store time in the state dictionary, and the default one
uses the proleptic Gregorian calendar used by the default Python datetime
.
-
sympl.
datetime
(year, month, day, hour=0, minute=0, second=0, microsecond=0, tzinfo=None, calendar='proleptic_gregorian')[source]¶ Retrieves a datetime-like object with the requested calendar. Calendar types other than proleptic_gregorian require the netcdftime module to be installed.
Parameters: - year (int,) –
- month (int,) –
- day (int,) –
- hour (int, optional) –
- minute (int, optional) –
- second (int, optional) –
- microsecond (int, optional) –
- tzinfo (datetime.tzinfo, optional) – A timezone informaton class, such as from pytz. Can only be used with ‘proleptic_gregorian’ calendar, as netcdftime does not support timezones.
- calendar (string, optional) – Should be one of ‘proleptic_gregorian’, ‘no_leap’, ‘365_day’, ‘all_leap’, ‘366_day’, ‘360_day’, ‘julian’, or ‘gregorian’. Default is ‘proleptic_gregorian’, which returns a normal Python datetime. Other options require the netcdftime module to be installed.
Returns: datetime – The requested datetime. May be a Python datetime, or one of the datetime-like types in netcdftime.
Return type: datetime-like
Naming Quantities¶
If you are a model user, the names of your quantities should coincide with the names used by the components you are using in your model. Basically, the components you are using dictate what quantity names you must use. If you are a model developer, we have a set of guidelines for naming quantities.
Note
The following is intended for model developers.
All quantity names should be verbose, and fully descriptive. Within a component you can set a quantity to an abbreviated variable, such as
theta = state['air_potential_temperature']
This ensures that your code is self-documenting. It is immediately apparent to anyone reading your code that theta refers to potential temperature of air, even if they are not familiar with theta as a common abbreviation.
We strongly recommend using the standard names according to CF conventions. In addition to making sure your code is self-documenting, this helps make sure that different components are compatible with one another, since they all need to use the same name for a given quantity in the model state.
If your quantity is on vertical interface levels, you should name it using the form “<name>_on_interface_levels”. If this is not specified, it is assumed that the quantity is on vertical mid levels. This is necessary because the same quantity may be specified on both mid and interface levels in the same model state.
When in doubt about names, look at what other components have been written that use the same quantity. If it looks like their name is verbose and follows the CF conventions then you should probably use the same name.
Constants¶
Configuration is an important part of any modelling framework. In Sympl, component-specific configuration is given to components directly. However, configuration values that may be shared by more than one component are stored as constants. Good examples of these are physical constants, such as gravitational_acceleration, or constants specifying processor counts.
Getting and Setting Constants¶
You can retrieve and set constants using get_constant()
and
set_constant()
. set_constant()
will
allow you to set constants regardless of whether a value is already defined
for that constant, allowing you to add new constants we haven’t thought of.
The constant library can be reverted to its original state when Sympl is
imported by calling reset_constants()
.
-
sympl.
get_constant
(name, units)[source]¶ Retrieves the value of a constant.
Parameters: - name (str) – The name of the constant.
- units (str) – The units requested for the returned value.
Returns: value – The value of the constant in the requested units.
Return type: float
Debugging and Logging Constants¶
You can get a string describing current constants by calling get_constants_string()
:
import sympl
print(sympl.get_constants_string())
Condensible Quantities¶
For Earth system modeling, water is used as a condensible compound. By
default, condensible quantities such as ‘density_of_ice’ and
‘heat_capacity_of_liquid_phase’ are aliases for the corresponding value for
water. If you would like to use a different condensible compound, you can
use the set_condensible_name()
function. For example:
import sympl
sympl.set_condensible_name('carbon_dioxide')
sympl.get_constant('heat_capacity_of_solid_phase', 'J kg^-1 K^-1')
will set the condensible compound to carbon dioxide, and then get the heat capacity of solid carbon dioxide (if it has been set). For example, the constant name ‘heat_capacity_of_solid_phase’ would then be an alias for ‘heat_capacity_of_solid_carbon_dioxide’.
When setting the value of an alias, the value of the aliased quantity is the one which will be altered. For example, if you run
import sympl
sympl.set_constant('heat_capacity_of_liquid_phase', 1.0, 'J kg^-1 K^-1')
you would change the heat capacity of liquid water (since water is the default condensible compound).
Default Constants¶
The following constants are available in Sympl by default:
-
class
sympl._core.constants.
ConstantList
[source]¶ - Planetary
gravitational_acceleration: 9.80665 m s^-2
planetary_radius: 6371000.0 m
planetary_rotation_rate: 7.292e-05 s^-1
seconds_per_day: 86400.0
- Chemical
heat_capacity_of_water_vapor_at_constant_pressure: 1846.0 J kg^-1 K^-1
density_of_liquid_water: 1000.0 kg m^-3
gas_constant_of_water_vapor: 461.5 J kg^-1 K^-1
latent_heat_of_vaporization_of_water: 2500000.0 J kg^-1
heat_capacity_of_liquid_water: 4185.5 J kg^-1 K^-1
latent_heat_of_fusion_of_water: 333550.0 J kg^-1
heat_capacity_of_solid_water_as_ice: 2108.0 J kg^-1 K^-1
heat_capacity_of_solid_water_as_snow: 2108.0 J kg^-1 K^-1
thermal_conductivity_of_solid_water_as_ice: 2.22 W m^-1 K^-1
thermal_conductivity_of_solid_water_as_snow: 0.2 W m^-1 K^-1
thermal_conductivity_of_liquid_water: 0.57 W m^-1 K^-1
density_of_solid_water_as_ice: 916.7 kg m^-3
density_of_solid_water_as_snow: 100.0 kg m^-3
freezing_temperature_of_liquid_water: 273.0 K
specific_enthalpy_of_water_vapor: 2500.0 J kg^-1
density_of_snow: 100.0 kg m^-3
heat_capacity_of_snow: 2108.0 J kg^-1 K^-1
heat_capacity_of_ice: 2108.0 J kg^-1 K^-1
density_of_ice: 916.7 kg m^-3
thermal_conductivity_of_ice: 2.22 W m^-1 K^-1
thermal_conductivity_of_snow: 0.2 W m^-1 K^-1
- Condensible
density_of_liquid_phase: 1000.0 kg m^-3
heat_capacity_of_liquid_phase: 4185.5 J kg^-1 K^-1
heat_capacity_of_vapor_phase: 1846.0 J kg^-1 K^-1
specific_enthalpy_of_vapor_phase: 2500.0 J kg^-1
gas_constant_of_vapor_phase: 461.5 J kg^-1 K^-1
latent_heat_of_condensation: 2500000.0 J kg^-1
latent_heat_of_fusion: 333550.0 J kg^-1
density_of_solid_phase_as_ice: 916.7 kg m^-3
density_of_solid_phase_as_snow: 100.0 kg m^-3
heat_capacity_of_solid_phase_as_ice: 2108.0 J kg^-1 K^-1
heat_capacity_of_solid_phase_as_snow: 2108.0 J kg^-1 K^-1
thermal_conductivity_of_solid_phase_as_ice: 2.22 W m^-1 K^-1
thermal_conductivity_of_solid_phase_as_snow: 0.2 W m^-1 K^-1
thermal_conductivity_of_liquid_phase: 0.57 W m^-1 K^-1
freezing_temperature_of_liquid_phase: 273.0 K
enthalpy_of_fusion: 333550.0 J kg^-1
latent_heat_of_vaporization: 2500000.0 J kg^-1
- Physical
stefan_boltzmann_constant: 5.670367e-08 W m^-2 K^-4
avogadro_constant: 6.022140857e+23 mole^-1
speed_of_light: 299792458.0 m s^-1
boltzmann_constant: 1.38064852e-23 J K^-1
loschmidt_constant: 2.6516467e+25 m^-3
universal_gas_constant: 8.3144598 J mole^-1 K^-1
planck_constant: 6.62607004e-34 J s
- Atmospheric
heat_capacity_of_dry_air_at_constant_pressure: 1004.64 J kg^-1 K^-1
gas_constant_of_dry_air: 287.0 J kg^-1 K^-1
thermal_conductivity_of_dry_air: 0.026 W m^-1 K^-1
reference_air_pressure: 101320.0 Pa
reference_air_temperature: 300.0 degK
- Stellar
stellar_irradiance: 1367.0 W m^-2
solar_constant: 1367.0 W m^-2
Timestepping¶
TendencyStepper
objects use time derivatives from
TendencyComponent
objects to step a model state forward in time.
They are initialized using any number of TendencyComponent
objects.
from sympl import AdamsBashforth
time_stepper = AdamsBashforth(MyPrognostic(), MyOtherPrognostic())
Once initialized, a TendencyStepper
object has a very similar
interface to the Stepper
object.
from datetime import timedelta
time_stepper = AdamsBashforth(MyPrognostic())
timestep = timedelta(minutes=10)
diagnostics, next_state = time_stepper(state, timestep)
state.update(diagnostics)
The returned diagnostics
dictionary contains diagnostic quantities from
the timestep of the input state
, while next_state
is the state
dictionary for the next timestep. It is possible that some of the arrays in
diagnostics
may be the same arrays as were given in the input state
,
and that they have been modified. In other words, state
may be modified by
this call. For instance, the time filtering necessary when using Leapfrog
time stepping means the current model state has to be modified by the filter.
It is only after calling the TendencyStepper
and getting the
diagnostics that you will have a complete state with all diagnostic quantities.
This means you will sometimes want to pass state
to your
Monitor
objects after calling
the TendencyStepper
and getting next_state
.
Warning
TendencyStepper
objects do not, and should not,
update ‘time’ in the model state.
Keep in mind that for split-time models, multiple TendencyStepper
objects might be called in in a single pass of the main loop. If each one
updated state['time']
, the time would be moved forward more than it should.
For that reason, TendencyStepper
objects do not update
state['time']
.
There are also
Stepper
objects which evolve the state forward in time
without the use of TendencyComponent objects. These function exactly the same as a
TendencyStepper
once they are created, but do not accept
TendencyComponent
objects when you create them. One example might
be a component that condenses all supersaturated moisture over some time period.
Stepper
objects are generally used for parameterizations
that work by determining the target model state in some way, or involve
limiters, and cannot be represented as a TendencyComponent
.
-
class
sympl.
TendencyStepper
(*args, **kwargs)[source]¶ An object which integrates model state forward in time.
It uses TendencyComponent and DiagnosticComponent objects to update the current model state with diagnostics, and to return the model state at the next timestep.
-
diagnostic_properties
¶ A dictionary whose keys are quantities for which values for the old state are returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
output_properties
¶ A dictionary whose keys are quantities for which values for the new state are returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
prognostic
¶ A composite of the TendencyComponent and ImplicitPrognostic objects used by the TendencyStepper.
Type: ImplicitTendencyComponentComposite
-
prognostic_list
¶ A list of TendencyComponent objects called by the TendencyStepper. These should be referenced when determining what inputs are necessary for the TendencyStepper.
Type: list of TendencyComponent and ImplicitPrognosticComponent
-
tendencies_in_diagnostics
¶ A boolean indicating whether this object will put tendencies of quantities in its diagnostic output.
Type: bool
-
time_unit_name
¶ The unit to use for time differencing when putting tendencies in diagnostics.
Type: str
-
time_unit_timedelta
¶ A timedelta corresponding to a single time unit as used for time differencing when putting tendencies in diagnostics.
Type: timedelta
-
name
¶ A label to be used for this object, for example as would be used for Y in the name “X_tendency_from_Y”.
Type: string
-
__call__
(state, timestep)[source]¶ Retrieves any diagnostics and returns a new state corresponding to the next timestep.
Parameters: - state (dict) – The current model state.
- timestep (timedelta) – The amount of time to step forward.
Returns: - diagnostics (dict) – Diagnostics from the timestep of the input state.
- new_state (dict) – The model state at the next timestep.
-
__init__
(*args, **kwargs)[source]¶ Initialize the TendencyStepper.
Parameters: - *args (TendencyComponent or ImplicitTendencyComponent) – Objects to call for tendencies when doing time stepping.
- tendencies_in_diagnostics (bool, optional) – A boolean indicating whether this object will put tendencies of quantities in its diagnostic output. Default is False. If set to True, you probably want to give a name also.
- name (str) – A label to be used for this object, for example as would be used for Y in the name “X_tendency_from_Y”. By default the class name is used.
-
array_call
(state, timestep)[source]¶ Gets diagnostics from the current model state and steps the state forward in time according to the timestep.
Parameters: - state (dict) – A numpy array state dictionary. Instead of data arrays, should include numpy arrays that satisfy the input_properties of this object.
- timestep (timedelta) – The amount of time to step forward.
Returns: - diagnostics (dict) – Diagnostics from the timestep of the input state, as numpy arrays.
- new_state (dict) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the timestep after input state, as numpy arrays.
-
-
class
sympl.
AdamsBashforth
(*args, **kwargs)[source]¶ A TendencyStepper using the Adams-Bashforth scheme.
-
__init__
(*args, **kwargs)[source]¶ Initialize an Adams-Bashforth time stepper.
Parameters: - *args (TendencyComponent or ImplicitTendencyComponent) – Objects to call for tendencies when doing time stepping.
- order (int, optional) – The order of accuracy to use. Must be between 1 and 4. 1 is the same as the Euler method. Default is 3.
-
-
class
sympl.
Leapfrog
(*args, **kwargs)[source]¶ A TendencyStepper using the Leapfrog scheme.
This scheme calculates the values at time $t_{n+1}$ using the derivatives at $t_{n}$ and values at $t_{n-1}$. Following the step, an Asselin filter is applied to damp the computational mode that results from the scheme and maintain stability. The Asselin filter brings the values at $t_{n}$ (and optionally the values at $t_{n+1}$, according to Williams (2009)) closer to the mean of the values at $t_{n-1}$ and $t_{n+1}$.
-
__init__
(*args, **kwargs)[source]¶ Initialize a Leapfrog time stepper.
Parameters: - *args (TendencyComponent or ImplicitTendencyComponent) – Objects to call for tendencies when doing time stepping.
- asselin_strength (float, optional) – The filter parameter used to determine the strength of the Asselin filter. Default is 0.05.
- alpha (float, optional) – Constant from Williams (2009), where the midpoint is shifted by alpha*influence, and the right point is shifted by (1-alpha)*influence. If alpha is 1 then the behavior is that of the classic Robert-Asselin time filter, while if it is 0.5 the filter will conserve the three-point mean. Default is 0.5.
References
Williams, P., 2009: A Proposed Modification to the Robert-Asselin Time Filter. Mon. Wea. Rev., 137, 2538–2546, doi: 10.1175/2009MWR2724.1.
-
Component Types¶
In Sympl, computation is mainly performed using TendencyComponent
,
DiagnosticComponent
, and Stepper
objects.
Each of these types, once initialized, can be passed in a current model state.
TendencyComponent
objects use the state to return tendencies and
diagnostics at the current time. DiagnosticComponent
objects
return only diagnostics from the current time. Stepper
objects will take in a timestep along with the state, and then return the
next state as well as modifying the current state to include more diagnostics
(it is similar to a TendencyStepper
in how it is called).
In specific cases, it may be necessary to use a ImplicitTendencyComponent
object, which is discussed at the end of this section.
These classes themselves (listed in the previous paragraph) are not ones you can initialize (e.g. there is no one ‘prognostic’ scheme), but instead should be subclassed to contain computational code relevant to the model you’re running.
In addition to the computational functionality below, all components have “properties” for their inputs and outputs, which are described in the section Input/Output Properties.
Details on the internals of components and how to write them are in the section on Writing Components.
TendencyComponent¶
As stated above, TendencyComponent
objects use the state to return
tendencies and diagnostics at the current time. In a full model, the tendencies
are used by a time stepping scheme (in Sympl, a TendencyStepper
)
to determine the values of quantities at the next time.
You can call a TendencyComponent
directly to get diagnostics and
tendencies like so:
radiation = RRTMRadiation()
diagnostics, tendencies = radiation(state)
diagnostics
and tendencies
in this case will both be dictionaries,
similar to state
. Even if the TendencyComponent
being called
does not compute any diagnostics, it will still return an empty
diagnostics dictionary.
Usually, you will call a TendencyComponent object through a
TendencyStepper
that uses it to determine values at the next
timestep.
-
class
sympl.
TendencyComponent
(tendencies_in_diagnostics=False, name=None)[source]¶ -
input_properties
¶ A dictionary whose keys are quantities required in the state when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
tendency_properties
¶ A dictionary whose keys are quantities for which tendencies are returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
diagnostic_properties
¶ A dictionary whose keys are diagnostic quantities returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
tendencies_in_diagnostics
¶ A boolean indicating whether this object will put tendencies of quantities in its diagnostic output based on first order time differencing of output values.
Type: bool
-
name
¶ A label to be used for this object, for example as would be used for Y in the name “X_tendency_from_Y”.
Type: string
-
__call__
(state)[source]¶ Gets tendencies and diagnostics from the passed model state.
Parameters: state (dict) – A model state dictionary satisfying the input_properties of this object.
Returns: - tendencies (dict) – A dictionary whose keys are strings indicating state quantities and values are the time derivative of those quantities in units/second at the time of the input state.
- diagnostics (dict) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state.
Raises: KeyError
– If a required quantity is missing from the state.InvalidStateError
– If state is not a valid input for the TendencyComponent instance.
-
__init__
(tendencies_in_diagnostics=False, name=None)[source]¶ Initializes the Stepper object.
Parameters: - tendencies_in_diagnostics (bool, optional) – A boolean indicating whether this object will put tendencies of quantities in its diagnostic output.
- name (string, optional) – A label to be used for this object, for example as would be used for Y in the name “X_tendency_from_Y”. By default the class name in lowercase is used.
-
array_call
(state)[source]¶ Gets tendencies and diagnostics from the passed model state.
Parameters: state (dict) – A model state dictionary. Instead of data arrays, should include numpy arrays that satisfy the input_properties of this object. Returns: - tendencies (dict) – A dictionary whose keys are strings indicating state quantities and values are the time derivative of those quantities in units/second at the time of the input state, as numpy arrays.
- diagnostics (dict) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state, as numpy arrays.
-
-
class
sympl.
ConstantTendencyComponent
(tendencies, diagnostics=None, **kwargs)[source]¶ Prescribes constant tendencies provided at initialization.
-
input_properties
¶ A dictionary whose keys are quantities required in the state when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
tendency_properties
¶ A dictionary whose keys are quantities for which tendencies are returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
diagnostic_properties
¶ A dictionary whose keys are diagnostic quantities returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
input_scale_factors
¶ A (possibly empty) dictionary whose keys are quantity names and values are floats by which input values are scaled before being used by this object.
Type: dict
-
tendency_scale_factors
¶ A (possibly empty) dictionary whose keys are quantity names and values are floats by which tendency values are scaled before being returned by this object.
Type: dict
-
diagnostic_scale_factors
¶ A (possibly empty) dictionary whose keys are quantity names and values are floats by which diagnostic values are scaled before being returned by this object.
Type: dict
-
update_interval
¶ If not None, the component will only give new output if at least a period of update_interval has passed since the last time new output was given. Otherwise, it would return that cached output.
Type: timedelta
-
tendencies_in_diagnostics
¶ A boolean indicating whether this object will put tendencies of quantities in its diagnostic output based on first order time differencing of output values.
Type: boo
-
name
¶ A label to be used for this object, for example as would be used for Y in the name “X_tendency_from_Y”.
Type: string
Note
Any arrays in the passed dictionaries are not copied, so that if you were to modify them after passing them into this object, it would also modify the values inside this object.
-
__init__
(tendencies, diagnostics=None, **kwargs)[source]¶ Parameters: - tendencies (dict) – A dictionary whose keys are strings indicating state quantities and values are the time derivative of those quantities in units/second to be returned by this TendencyComponent.
- diagnostics (dict, optional) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities to be returned by this TendencyComponent. By default an empty dictionary is used.
- input_scale_factors (dict, optional) – A (possibly empty) dictionary whose keys are quantity names and values are floats by which input values are scaled before being used by this object.
- tendency_scale_factors (dict, optional) – A (possibly empty) dictionary whose keys are quantity names and values are floats by which tendency values are scaled before being returned by this object.
- diagnostic_scale_factors (dict, optional) – A (possibly empty) dictionary whose keys are quantity names and values are floats by which diagnostic values are scaled before being returned by this object.
- update_interval (timedelta, optional) – If given, the component will only give new output if at least a period of update_interval has passed since the last time new output was given. Otherwise, it would return that cached output.
- tendencies_in_diagnostics (bool, optional) – A boolean indicating whether this object will put tendencies of quantities in its diagnostic output.
- name (string, optional) – A label to be used for this object, for example as would be used for Y in the name “X_tendency_from_Y”. By default the class name in lowercase is used.
-
array_call
(state)[source]¶ Gets tendencies and diagnostics from the passed model state.
Parameters: state (dict) – A model state dictionary. Instead of data arrays, should include numpy arrays that satisfy the input_properties of this object. Returns: - tendencies (dict) – A dictionary whose keys are strings indicating state quantities and values are the time derivative of those quantities in units/second at the time of the input state, as numpy arrays.
- diagnostics (dict) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state, as numpy arrays.
-
-
class
sympl.
RelaxationTendencyComponent
(quantity_name, units, **kwargs)[source]¶ Applies Newtonian relaxation to a single quantity.
The relaxation takes the form \(\frac{dx}{dt} = - \frac{x - x_{eq}}{\tau}\) where \(x\) is the quantity being relaxed, \(x_{eq}\) is the equilibrium value, and \(\tau\) is the timescale of the relaxation.
-
input_properties
¶ A dictionary whose keys are quantities required in the state when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
tendency_properties
¶ A dictionary whose keys are quantities for which tendencies are returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
diagnostic_properties
¶ A dictionary whose keys are diagnostic quantities returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
input_scale_factors
¶ A (possibly empty) dictionary whose keys are quantity names and values are floats by which input values are scaled before being used by this object.
Type: dict
-
tendency_scale_factors
¶ A (possibly empty) dictionary whose keys are quantity names and values are floats by which tendency values are scaled before being returned by this object.
Type: dict
-
diagnostic_scale_factors
¶ A (possibly empty) dictionary whose keys are quantity names and values are floats by which diagnostic values are scaled before being returned by this object.
Type: dict
-
update_interval
¶ If not None, the component will only give new output if at least a period of update_interval has passed since the last time new output was given. Otherwise, it would return that cached output.
Type: timedelta
-
tendencies_in_diagnostics
¶ A boolean indicating whether this object will put tendencies of quantities in its diagnostic output based on first order time differencing of output values.
Type: boo
-
name
¶ A label to be used for this object, for example as would be used for Y in the name “X_tendency_from_Y”.
Type: string
-
__init__
(quantity_name, units, **kwargs)[source]¶ Parameters: - quantity_name (str) – The name of the quantity to which Newtonian relaxation should be applied.
- units (str) – The units of the relaxed quantity as to be used internally when computing tendency. Can be any units convertible from the actual input you plan to use.
- input_scale_factors (dict, optional) – A (possibly empty) dictionary whose keys are quantity names and values are floats by which input values are scaled before being used by this object.
- tendency_scale_factors (dict, optional) – A (possibly empty) dictionary whose keys are quantity names and values are floats by which tendency values are scaled before being returned by this object.
- diagnostic_scale_factors (dict, optional) – A (possibly empty) dictionary whose keys are quantity names and values are floats by which diagnostic values are scaled before being returned by this object.
- update_interval (timedelta, optional) – If given, the component will only give new output if at least a period of update_interval has passed since the last time new output was given. Otherwise, it would return that cached output.
-
array_call
(state)[source]¶ Gets tendencies and diagnostics from the passed model state.
Parameters: state (dict) – A model state dictionary as numpy arrays. Below, (quantity_name) refers to the quantity_name passed at initialization. The state must contain:
- (quantity_name)
- equilibrium_(quantity_name), unless this was passed at initialisation time in which case that value is used
- (quantity_name)_relaxation_timescale, unless this was passed at initialisation time in which case that value is used
Returns: - tendencies (dict) – A dictionary whose keys are strings indicating state quantities and values are the time derivative of those quantities in units/second at the time of the input state, as numpy arrays.
- diagnostics (dict) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state, as numpy arrays.
-
DiagnosticComponent¶
DiagnosticComponent
objects use the state to return quantities
(‘diagnostics’) from the same timestep as the input state. You can call a
DiagnosticComponent
directly to get diagnostic quantities like so:
diagnostic_component = MyDiagnostic()
diagnostics = diagnostic_component(state)
You should be careful to check in the documentation of the particular
DiagnosticComponent
you are using to see whether it modifies the
state
given to it as input. DiagnosticComponent
objects in charge
of updating ghost cells in particular may modify the arrays in the input
dictionary, so that the arrays in the returned diagnostics
dictionary are
the same ones as were sent as input in the state
. To make it clear that
the state is being modified when using such objects, we recommend using a
syntax like:
state.update(diagnostic_component(state))
-
class
sympl.
DiagnosticComponent
[source]¶ -
input_properties
¶ A dictionary whose keys are quantities required in the state when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
diagnostic_properties
¶ A dictionary whose keys are diagnostic quantities returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
__call__
(state)[source]¶ Gets diagnostics from the passed model state.
Parameters: state (dict) – A model state dictionary satisfying the input_properties of this object.
Returns: diagnostics – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state.
Return type: dict
Raises: KeyError
– If a required quantity is missing from the state.InvalidStateError
– If state is not a valid input for the TendencyComponent instance.
-
array_call
(state)[source]¶ Gets diagnostics from the passed model state.
Parameters: state (dict) – A model state dictionary. Instead of data arrays, should include numpy arrays that satisfy the input_properties of this object. Returns: diagnostics – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state, as numpy arrays. Return type: dict
-
-
class
sympl.
ConstantDiagnosticComponent
(diagnostics, **kwargs)[source]¶ Yields constant diagnostics provided at initialization.
-
input_properties
¶ A dictionary whose keys are quantities required in the state when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
diagnostic_properties
¶ A dictionary whose keys are diagnostic quantities returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
input_scale_factors
¶ A (possibly empty) dictionary whose keys are quantity names and values are floats by which input values are scaled before being used by this object.
Type: dict
-
diagnostic_scale_factors
¶ A (possibly empty) dictionary whose keys are quantity names and values are floats by which diagnostic values are scaled before being returned by this object.
Type: dict
-
update_interval
¶ If not None, the component will only give new output if at least a period of update_interval has passed since the last time new output was given. Otherwise, it would return that cached output.
Type: timedelta
Note
Any arrays in the passed dictionaries are not copied, so that if you were to modify them after passing them into this object, it would also modify the values inside this object.
-
__init__
(diagnostics, **kwargs)[source]¶ Parameters: - diagnostics (dict) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities. The values in the dictionary will be returned when this DiagnosticComponent is called.
- input_scale_factors (dict, optional) – A (possibly empty) dictionary whose keys are quantity names and values are floats by which input values are scaled before being used by this object.
- diagnostic_scale_factors (dict, optional) – A (possibly empty) dictionary whose keys are quantity names and values are floats by which diagnostic values are scaled before being returned by this object.
- update_interval (timedelta, optional) – If given, the component will only give new output if at least a period of update_interval has passed since the last time new output was given. Otherwise, it would return that cached output.
-
array_call
(state)[source]¶ Gets diagnostics from the passed model state.
Parameters: state (dict) – A model state dictionary. Instead of data arrays, should include numpy arrays that satisfy the input_properties of this object. Returns: diagnostics – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state, as numpy arrays. Return type: dict
-
Stepper¶
Stepper
objects use a state and a timestep to return the next
state, and update the input state with any relevant diagnostic quantities. You
can call an Stepper object like so:
from datetime import timedelta
implicit = MyImplicit()
timestep = timedelta(minutes=10)
diagnostics, next_state = implicit(state, timestep)
state.update(diagnostics)
The returned diagnostics
dictionary contains diagnostic quantities from
the timestep of the input state
, while next_state
is the state
dictionary for the next timestep. It is possible that some of the arrays in
diagnostics
may be the same arrays as were given in the input state
,
and that they have been modified. In other words, state
may be modified by
this call. For instance, the object may need to update ghost cells in the
current state. Or if an object provides ‘cloud_fraction’ as a diagnostic, it
may modify an existing ‘cloud_fraction’ array in the input state if one is
present, instead of allocating a new array.
-
class
sympl.
Stepper
(tendencies_in_diagnostics=False, name=None)[source]¶ -
input_properties
¶ A dictionary whose keys are quantities required in the state when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
diagnostic_properties
¶ A dictionary whose keys are quantities for which values for the old state are returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
output_properties
¶ A dictionary whose keys are quantities for which values for the new state are returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
tendencies_in_diagnostics
¶ A boolean indicating whether this object will put tendencies of quantities in its diagnostic output based on first order time differencing of output values.
Type: bool
-
time_unit_name
¶ The unit to use for time differencing when putting tendencies in diagnostics.
Type: str
-
time_unit_timedelta
¶ A timedelta corresponding to a single time unit as used for time differencing when putting tendencies in diagnostics.
Type: timedelta
-
name
¶ A label to be used for this object, for example as would be used for Y in the name “X_tendency_from_Y”.
Type: string
-
__call__
(state, timestep)[source]¶ Gets diagnostics from the current model state and steps the state forward in time according to the timestep.
Parameters: - state (dict) – A model state dictionary satisfying the input_properties of this object.
- timestep (timedelta) – The amount of time to step forward.
Returns: - diagnostics (dict) – Diagnostics from the timestep of the input state.
- new_state (dict) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the timestep after input state.
Raises: KeyError
– If a required quantity is missing from the state.InvalidStateError
– If state is not a valid input for the Stepper instance for other reasons.
-
__init__
(tendencies_in_diagnostics=False, name=None)[source]¶ Initializes the Stepper object.
Parameters: - tendencies_in_diagnostics (bool, optional) – A boolean indicating whether this object will put tendencies of quantities in its diagnostic output based on first order time differencing of output values.
- name (string, optional) – A label to be used for this object, for example as would be used for Y in the name “X_tendency_from_Y”. By default the class name in lowercase is used.
-
array_call
(state, timestep)[source]¶ Gets diagnostics from the current model state and steps the state forward in time according to the timestep.
Parameters: - state (dict) – A numpy array state dictionary. Instead of data arrays, should include numpy arrays that satisfy the input_properties of this object.
- timestep (timedelta) – The amount of time to step forward.
Returns: - diagnostics (dict) – Diagnostics from the timestep of the input state, as numpy arrays.
- new_state (dict) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the timestep after input state, as numpy arrays.
-
Input/Output Properties¶
You may have noticed when reading the documentation for the classes above that
there are a number of attributes with names like input_properties
for
components. These attributes give a fairly complete description of the inputs
and outputs of the component.
You can access them like this (for an example TendencyComponent
class RRTMRadiation
):
radiation = RRTMRadiation()
radiation.input_properties
radiation.diagnostic_properties
radiation.tendency_properties
Input¶
All components have input_properties, because they all take inputs. This
attribute (like all the other properties attributes) is a python dict
,
or “dictionary” (if you are unfamiliar with these, please read the Python
documentation for dicts).
An example input_properties would be
{
'air_temperature': {
'dims': ['*', 'z'],
'units': 'degK',
},
'vertical_wind': {
'dims': ['*', 'z'],
'units': 'm/s',
'match_dims_like': ['air_temperature']
}
}
Each entry in the input_properties dictionary is a quantity that the object
requires as an input, and its value is another dictionary that tells you how
the object uses that quantity. The units
property is the units used
internally in the object. You don’t need to pass in the quantity with those
those units, as long as the units can be converted, but if you do use the same
units in the input state it will avoid the computational cost of
converting units.
The dims
property can be more confusing, but is very useful. It says what
dimensions the component uses internally for those quantities. The component
requires that you give it quantities that can be transformed into those
internal dimensions, but it can take care of that transformation itself. In
this example, it will transform the arrays for both quantities to put the
vertical dimension last, and collect all the other dimensions into a single
first dimension. If you pass this object arrays that have their vertical
dimension last, it may speed up the computation, depending on the component
(but not for all components!).
So what are ‘*’ and ‘z’ anyways? These are wildcard dimensions. ‘z’ will
match any dimension that is vertical, while ‘*’ will match any dimension that
is not specified somewhere else in the dims
list. There are also ‘x’ and
‘y’ for horizontal dimensions. The directional matches are given to Sympl
using the functions set_direction_names()
or
add_direction_names()
. If you’re using someone else’s package
for a component, it is likely that they call these functions for you, so you
don’t have to (and if you’re writing such a package, you should use
add_direction_names()
).
If a component is using a wildcard it means it doesn’t care very much about those directions. For example, a column component like radiation will simply call itself on each column of the domain, so it doesn’t care about the specifics of what the non-vertical dimensions are, as long as the desired quantities are co-located.
That’s where match_dims_like
comes in. This property says the object
requires all shared wildcard dimensions between the two quantity match the
same dimensions as the other
specified quantity. In this case, it will ensure that vertical_wind
is on
the same grid as air_temperature
.
Let’s consider a slight variation on the earlier example:
{
'air_temperature': {
'dims': ['*', 'mid_levels'],
'units': 'degK',
},
'vertical_wind': {
'dims': ['*', 'interface_levels'],
'units': 'm/s',
'match_dims_like': ['air_temperature']
}
}
This version requires that air_temperature
be on the mid_levels
vertical
grid, while vertical_wind
is on the interface_levels
. It still requires
that all other dimensions are the same between the two quantities, so that they
are on the same horizontal grid (if they have a horizontal grid).
Outputs¶
There are a few output property dictionaries in Sympl: tendency_properties
,
diagnostic_properties
, and output_properties
. They are all formatted
the same way with the same properties, but tell you about the tendencies,
diagnostics, or next state values that are output by the component,
respectively.
Here’s an example output dictionary:
tendency_properties = {
'air_temperature': {
'dims_like': 'air_temperature',
'units': 'degK/s',
}
}
In tendency_properties, the quantity names specify the quantities for which
tendencies are given. The units
are the units of the output value, which
is also put in the output DataArray
as the units
attribute.
dims_like
is telling you that the output array will have the same dimensions
as the array you gave it for air_temperature
as an input. If you pass it
an air_temperature
array with (‘latitude’, ‘longitude’, ‘mid_levels’) as
its axes, it will return an array with (‘latitude’, ‘longitude’, ‘mid_levels’)
for the temperature tendency. If dims_like
is not specified in the
tendency_properties
dictionary, it is assumed to be the matching quantity
in the input, but for the other quantities dims_like
must always be
explicitly defined. For instance, if the object as a diagnostic_properties
equal to:
diagnostic_properties = {
'cloud_fraction': {
'dims_like': 'air_temperature',
'units': '',
}
}
that the object will output cloud_fraction
in its diagnostics on the
same grid as air_temperature
, in dimensionless units.
ImplicitTendencyComponent¶
Warning
This component type should be avoided unless you know you need it, for reasons discussed in this section.
In addition to the component types described above, computation may be performed by a
ImplicitTendencyComponent
. This class should be avoided unless you
know what you are doing, but it may be necessary in certain cases. An
ImplicitTendencyComponent
, like a TendencyComponent
,
calculates tendencies, but it does so using both the model state and a timestep.
Certain components, like ones handling advection using a spectral method, may
need to derive tendencies from an Stepper
object by
representing it using an ImplicitTendencyComponent
.
The reason to avoid using an ImplicitTendencyComponent
is that if
a component requires a timestep, it is making internal assumptions about how
you are timestepping. For example, it may use the timestep to ensure that all
supersaturated water is condensed by the end of the timestep using an assumption
about the timestepping. However, if you use a TendencyStepper
which does not obey those assumptions, you may get unintended behavior, such as
some supersaturated water remaining, or too much water being condensed.
For this reason, the TendencyStepper
objects included in Sympl
do not wrap ImplicitTendencyComponent
components. If you would like
to use this type of component, and know what you are doing, it is pretty easy
to write your own TendencyStepper
to do so (you can base the code
off of the code in Sympl), or the model you are using might already have
components to do this for you.
If you are wrapping a parameterization and notice that it needs a timestep to
compute its tendencies, that is likely not a good reason to write an
ImplicitTendencyComponent
. If at all possible you should modify the
code to compute the value at the next timestep, and write an
Stepper
component. You are welcome to reach out to the
developers of Sympl if you would like advice on your specific situation! We’re
always excited about new wrapped components.
-
class
sympl.
ImplicitTendencyComponent
(tendencies_in_diagnostics=False, name=None)[source]¶ -
input_properties
¶ A dictionary whose keys are quantities required in the state when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
tendency_properties
¶ A dictionary whose keys are quantities for which tendencies are returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
diagnostic_properties
¶ A dictionary whose keys are diagnostic quantities returned when the object is called, and values are dictionaries which indicate ‘dims’ and ‘units’.
Type: dict
-
tendencies_in_diagnostics
¶ A boolean indicating whether this object will put tendencies of quantities in its diagnostic output based on first order time differencing of output values.
Type: bool
-
name
¶ A label to be used for this object, for example as would be used for Y in the name “X_tendency_from_Y”.
Type: string
-
__call__
(state, timestep)[source]¶ Gets tendencies and diagnostics from the passed model state.
Parameters: - state (dict) – A model state dictionary satisfying the input_properties of this object.
- timestep (timedelta) – The time over which the model is being stepped.
Returns: - tendencies (dict) – A dictionary whose keys are strings indicating state quantities and values are the time derivative of those quantities in units/second at the time of the input state.
- diagnostics (dict) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state.
Raises: KeyError
– If a required quantity is missing from the state.InvalidStateError
– If state is not a valid input for the TendencyComponent instance.
-
__init__
(tendencies_in_diagnostics=False, name=None)[source]¶ Initializes the Stepper object.
Parameters: - tendencies_in_diagnostics (bool, optional) – A boolean indicating whether this object will put tendencies of quantities in its diagnostic output.
- name (string, optional) – A label to be used for this object, for example as would be used for Y in the name “X_tendency_from_Y”. By default the class name in lowercase is used.
-
array_call
(state, timestep)[source]¶ Gets tendencies and diagnostics from the passed model state.
Parameters: - state (dict) – A model state dictionary. Instead of data arrays, should include numpy arrays that satisfy the input_properties of this object.
- timestep (timedelta) – The time over which the model is being stepped.
Returns: - tendencies (dict) – A dictionary whose keys are strings indicating state quantities and values are the time derivative of those quantities in units/second at the time of the input state, as numpy arrays.
- diagnostics (dict) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state, as numpy arrays.
-
Tracer Properties¶
You may notice some components have properties that mention tracers. uses_tracers is a boolean that tells you whether the component makes use of tracers or not, while tracer_dims contains the dimensions of tracer arrays used internally by the object. For more on tracers as a user, see Tracers, and for more on tracers as a component author, see Writing Components.
Monitors¶
Monitor
objects store states in some way, whether it is by
displaying the new state on a plot that is shown to the user, updating
information on a web server, or saving the state to a file. They are called
like so:
monitor = MyMonitor()
monitor.store(state)
The Monitor
will take advantage of the ‘time’ key in the
state
dictionary in order to determine the model time of the state. This is
particularly important for a Monitor
which outputs a series
of states to disk.
-
class
sympl.
NetCDFMonitor
(filename, time_units='seconds', store_names=None, write_on_store=False, aliases=None)[source]¶ A Monitor which caches stored states and then writes them to a NetCDF file when requested.
-
__init__
(filename, time_units='seconds', store_names=None, write_on_store=False, aliases=None)[source]¶ Parameters: - filename (str) – The file to which the NetCDF file will be written.
- time_units (str, optional) – The units in which time will be stored in the NetCDF file. Time is stored as an integer number of these units. Default is seconds.
- store_names (iterable of str, optional) – Names of quantities to store. If not given, all quantities are stored.
- write_on_store (bool, optional) – If True, stored changes are immediately written to file. This can result in many file open/close operations. Default is to write only when the write() method is called directly.
- aliases (dict) – A dictionary of string replacements to apply to state variable names before saving them in netCDF files.
-
store
(state)[source]¶ Caches the given state. If write_on_store=True was passed on initialization, also writes to file. Normally a call to the write() method is required to write to file.
Parameters: state (dict) – A model state dictionary. Raises: InvalidStateError
– If state is not a valid input for the DiagnosticComponent instance.
-
-
class
sympl.
PlotFunctionMonitor
(plot_function, interactive=True)[source]¶ A Monitor which uses a user-defined function to draw figures using model state.
-
__init__
(plot_function, interactive=True)[source]¶ Initialize a PlotFunctionMonitor.
Parameters: - plot_function (func) – A function plot_function(fig, state) that draws the given state onto the given (initially clear) figure.
- interactive (bool, optional) – If true, matplotlib’s interactive mode will be enabled, allowing plot animation while other computation is running.
-
Composites¶
There are a set of objects in Sympl that wrap multiple components into a single
object so they can be called as if they were one component. There is one each
for TendencyComponent
, DiagnosticComponent
, and
Monitor
. These can be used to simplify code, so that
the way you call a list of components is the same as the way you would
call a single component. For example, instead of writing:
tendency_component_list = [
MyTendencyComponent(),
MyOtherTendencyComponent(),
YetAnotherTendencyComponent(),
]
all_diagnostics = {}
total_tendencies = {}
for tendency_component in tendency_component_list:
tendencies, diagnostics = tendency_component(state)
# this should actually check to make sure nothing is overwritten,
# but this code does not
total_tendencies.update(tendencies)
for name, value in tendencies.keys():
if name not in total_tendencies:
total_tendencies[name] = value
else:
total_tendencies[name] += value
for name, value in diagnostics.items():
all_diagnostics[name] = value
You could write:
tendency_component_composite = TendencyComponentComposite([
MyTendencyComponent(),
MyOtherTendencyComponent(),
YetAnotherTendencyComponent(),
])
tendencies, diagnostics = tendency_component_composite(state)
This second call is much cleaner. It will also automatically detect whether
multiple components are trying to write out the same diagnostic, and raise
an exception if that is the case (so no results are being silently
overwritten). You can get similar simplifications for
DiagnosticComponent
and Monitor
.
Note
TendencyComponentComposites are mainly useful inside of TimeSteppers, so if you’re only writing a model script it’s unlikely you’ll need them.
API Reference¶
-
class
sympl.
TendencyComponentComposite
(*args)[source]¶ -
__call__
(state)[source]¶ Gets tendencies and diagnostics from the passed model state.
Parameters: state (dict) – A model state dictionary.
Returns: - tendencies (dict) – A dictionary whose keys are strings indicating state quantities and values are the time derivative of those quantities in units/second at the time of the input state.
- diagnostics (dict) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state.
Raises: KeyError
– If a required quantity is missing from the state.InvalidStateError
– If state is not a valid input for a TendencyComponent instance.
-
__init__
(*args)[source]¶ Parameters: *args – The components that should be wrapped by this object.
Raises: SharedKeyError
– If two components compute the same diagnostic quantity.InvalidPropertyDictError
– If two components require the same input or compute the same output quantity, and their dimensions or units are incompatible with one another.
-
array_call
(state)[source]¶ Gets tendencies and diagnostics from the passed model state.
Parameters: state (dict) – A model state dictionary. Instead of data arrays, should include numpy arrays that satisfy the input_properties of this object. Returns: - tendencies (dict) – A dictionary whose keys are strings indicating state quantities and values are the time derivative of those quantities in units/second at the time of the input state, as numpy arrays.
- diagnostics (dict) – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state, as numpy arrays.
-
component_class
¶ alias of
sympl._core.base_components.TendencyComponent
-
-
class
sympl.
DiagnosticComponentComposite
(*args)[source]¶ -
__call__
(state)[source]¶ Gets diagnostics from the passed model state.
Parameters: state (dict) – A model state dictionary.
Returns: diagnostics – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state.
Return type: dict
Raises: KeyError
– If a required quantity is missing from the state.InvalidStateError
– If state is not a valid input for a DiagnosticComponent instance.
-
array_call
(state)[source]¶ Gets diagnostics from the passed model state.
Parameters: state (dict) – A model state dictionary. Instead of data arrays, should include numpy arrays that satisfy the input_properties of this object. Returns: diagnostics – A dictionary whose keys are strings indicating state quantities and values are the value of those quantities at the time of the input state, as numpy arrays. Return type: dict
-
component_class
¶ alias of
sympl._core.base_components.DiagnosticComponent
-
Tracers¶
In an Earth system model, “tracer” refers to quantities that are passively moved around in a model without actively interacting with a component. Generally these are moved around by a dynamical core or subgrid advection scheme. It is possible for components to do something else to tracers (let us know if you think of something!) but for now let’s assume that’s what’s going on.
If a component moves around tracers, it will have its uses_tracers
property
set to True
, and will also have a tracer_dims
property set.
You can tell Sympl that you want components to move a tracer around by
registering it with register_tracer()
.
-
sympl.
register_tracer
(name, units)[source]¶ Parameters: - name (str) – Quantity name to register as a tracer.
- units (str) – Unit string of that quantity.
To see the current list of registered tracers, you can call
get_tracer_names()
or
get_tracer_unit_dict()
.
Units¶
While nearly all unit functionality is handled internally within Sympl, it does expose a few helper functions which may be useful to you.
-
sympl.
is_valid_unit
(unit_string)[source]¶ Returns True if the unit string is recognized, and False otherwise.
Writing Components¶
Note
This section is intended for model developers. If you intend to use only components that are already written, you can probably ignore it.
Perhaps the best way to learn how to write components is to read components someone else has written. For example, you can look at the CliMT project. Here we will go over a couple examples of physically simple, made-up components to talk about the parts of their code.
Writing an Example¶
Let’s start with a TendencyComponent component which relaxes temperature towards some target temperature. We’ll go over the sections of this example step-by-step below.
from sympl import (
TendencyComponent, get_numpy_arrays_with_properties,
restore_data_arrays_with_properties)
class TemperatureRelaxation(TendencyComponent):
input_properties = {
'air_temperature': {
'dims': ['*'],
'units': 'degK',
},
'vertical_wind': {
'dims': ['*'],
'units': 'm/s',
'match_dims_like': ['air_temperature']
}
}
diagnostic_properties = {}
tendency_properties = {
'air_temperature': {
'dims_like': 'air_temperature',
'units': 'degK/s',
}
}
def __init__(self, damping_timescale_seconds=1., target_temperature_K=300.):
self._tau = damping_timescale_seconds
self._T0 = target_temperature_K
def array_call(self, state):
tendencies = {
'air_temperature': (state['air_temperature'] - self._T0)/self._tau,
}
diagnostics = {}
return tendencies, diagnostics
Imports¶
There are a lot of parts to that code, so let’s go through some of them step-by-step. First we have to import objects and functions from Sympl that we plan to use. The import statement should always go at the top of your file so that it can be found right away by anyone reading your code.
from sympl import (
TendencyComponent, get_numpy_arrays_with_properties,
restore_data_arrays_with_properties)
Define an Object¶
Once these are imported, there’s this line:
class TemperatureRelaxation(TendencyComponent):
This is the syntax for defining an object in Python. TemperatureRelaxation
will be the name of the new object. The TendencyComponent
in parentheses is telling Python that TemperatureRelaxation
is a subclass of
TendencyComponent
. This tells Sympl that it can expect your object
to behave like a TendencyComponent
.
Define Attributes¶
The next few lines define attributes of your object:
input_properties = {
'air_temperature': {
'dims': ['*'],
'units': 'degK',
},
'eastward_wind': {
'dims': ['*'],
'units': 'm/s',
'match_dims_like': ['air_temperature']
}
}
diagnostic_properties = {}
tendency_properties = {
'air_temperature': {
'dims_like': 'air_temperature',
'units': 'degK/s',
}
}
Note
‘eastward_wind’ wouldn’t normally make sense as an input for this object, it’s only included so we can talk about match_dims_like.
These attributes will be attributes both of the class object you’re defining and of any instances of that object. That means you can access them using:
TemperatureRelaxation.input_properties
or on an instance, as when you do:
prognostic = TemperatureRelaxation()
prognostic.input_properties
These properties are described in Component Types. They are very useful! They clearly document your code. Here we can see that air_temperature will be used as a 1-dimensional flattened array in units of degrees Kelvin. Sympl uses these properties to automatically acquire arrays in the dimensions and units that you need, and to automatically convert your output back into a form consistent with the dimensions of the model state. It will warn you if you create extra outputs which are not defined in the properties, or if there is an output defined in the properties that is missing.
It is possible that some of these attributes won’t be known until you
create the object (they may depend on things passed in on initialization).
If that’s the case, you can write the __init__
method (see below) so that
it sets any relevant properties like self.input_properties
to have the
correct values.
Initialization Method¶
Next we see a method being defined for this class, which may seem to have a weird name:
def __init__(self, damping_timescale_seconds=1., target_temperature_K=300.):
self._tau = damping_timescale_seconds
self._T0 = target_temperature_K
This is the function that is called when you create an instance of your object.
All methods on objects take in a first argument called self
. You don’t see
it when you call those methods, it gets added in automatically. self
is
a variable that refers to the object on which the method is being called -
it’s the object itself! When you store attributes on self, as we see in this
code, they stay there. You can access them when the object is called later.
Notice some things about the way variables have been named in this __init__
.
The parameters are fairly verbose names which almost fully describe what they
are (apart from the units, which are in the documentation string). This is
best because it is entirely clear what these values are when others are using
your object. You write code for people, not computers! Compilers write code for
computers.
Then we take these inputs and store them as attributes with shorter names. This is also optimal. What these attributes mean is clearly defined in the two lines:
self._tau = damping_timescale_seconds
self._T0 = target_temperature_K
Obviously self._tau
is the damping timescale, and self._T0
is the
target temperature for the relaxation. Now you can use these shorter variables
in the actual code to keep long lines for equations short, knowing that your
variables are well-documented.
The Computation¶
That brings us to the array_call
method. In Sympl components, this is the
method which takes in a state dictionary as numpy arrays (not DataArray
)
and returns dictionaries with numpy array outputs.
def array_call(self, state):
tendencies = {
'air_temperature': (state['air_temperature'] - self._T0)/self._tau,
}
diagnostics = {}
return tendencies, diagnostics
Sympl will automatically handle taking in the input state of DataArray
objects and converting it to the form defined by the input_properties
of your
component. It will convert units to ensure the numbers are in the specified
units, and it will reshape the data to give it the shape specified in dims
.
For example, if dims is ['*', 'mid_levels']
then it will give you a
2-dimensional array whose second axis is the vertical on mid levels, and first
axis is a flattening of any other dimensions. The match_dims_like
property on air_pressure
tells Sympl that any wildcard-matched dimensions
(‘*’) should be the same between the two
quantities, meaning they’re on the same grid for those wildcards. You can still,
however, have one be on say ‘mid_levels’ and another on ‘interface_levels’ if
those dimensions are explicitly listed.
After you return dictionaries of numpy arrays, Sympl will convert these outputs
back to DataArray
objects. In this example, it takes the tendencies
dictionary and converts the value for ‘air_temperature’ from a numpy array to a
DataArray
that has the same dimensions as air_temperature
had in the
input state. That means that you could pass this object a state with whatever
dimensions you want, whether it’s ('longitude', 'latitude', 'mid_levels')
, or
('interface_levels',)
or ('station_number', 'planet_number')
, etc.
and this component will be able to take in that state, and return a
tendency dictionary with the same dimensions (and order) that the model uses!
And internally you can work with a simple 1-dimensional array. This is
particularly useful for writing pointwise components using ['*']
or column
components with, for example, ['*', 'mid_levels']
or
['interface_levels', '*']
.
You can read more about properties in the section Input/Output Properties.
-
sympl.
restore_data_arrays_with_properties
(raw_arrays, output_properties, input_state, input_properties, ignore_names=None, ignore_missing=False)[source]¶ Parameters: - raw_arrays (dict) – A dictionary whose keys are quantity names and values are numpy arrays containing the data for those quantities.
- output_properties (dict) – A dictionary whose keys are quantity names and values are dictionaries with properties for those quantities. The property “dims” must be present for each quantity not also present in input_properties. All other properties are included as attributes on the output DataArray for that quantity, including “units” which is required.
- input_state (dict) – A state dictionary that was used as input to a component for which DataArrays are being restored.
- input_properties (dict) – A dictionary whose keys are quantity names and values are dictionaries with input properties for those quantities. The property “dims” must be present, indicating the dimensions that the quantity was transformed to when taken as input to a component.
- ignore_names (iterable of str, optional) – Names to ignore when encountered in output_properties, will not be included in the returned dictionary.
- ignore_missing (bool, optional) – If True, ignore any values in output_properties not present in raw_arrays rather than raising an exception. Default is False.
Returns: out_dict – A dictionary whose keys are quantities and values are DataArrays corresponding to those quantities, with data, shapes and attributes determined from the inputs to this function.
Return type: dict
Raises: InvalidPropertyDictError
– When an output property is specified to have dims_like an input property, but the arrays for the two properties have incompatible shapes.
Aliases¶
Note
Using aliases isn’t necessary, but it may make your code easier to read if you have long quantity names
Let’s say if instead of the properties we set before, we have
input_properties = {
'air_temperature': {
'dims': ['*'],
'units': 'degK',
'alias': 'T',
},
'eastward_wind': {
'dims': ['*'],
'units': 'm/s',
'match_dims_like': ['air_temperature']
'alias': 'u',
}
}
The difference here is we’ve set ‘T’ and ‘u’ to be aliases for ‘air_temperature’ and ‘eastward_wind’. What does that mean? Well, in the computational code, we can write:
def array_call(self, state):
tendencies = {
'T': (state['T'] - self._T0)/self._tau,
}
diagnostics = {}
return tendencies, diagnostics
Instead of using ‘air_temperature’ in the raw_arrays and raw_tendencies dictionaries, we can use ‘T’. This doesn’t matter much for a name as short as air_temperature, but it might matter for longer names like ‘correlation_of_eastward_wind_and_liquid_water_potential_temperature_on_interface_levels’.
Also notice that even though the alias is set in input_properties, it is also
used when restoring DataArrays. If there is an output that is not
also an input, the alias could instead be set in diagnostic_properties
,
tendency_properties
, or output_properties
, wherever is relevant.
Using Tracers¶
Note
This feature is mostly used in dynamical cores. If you don’t think you need this, you probably don’t.
Sympl’s base components have some features to automatically create tracer arrays
for use by dynamical components. If an Stepper
,
TendencyComponent
, or ImplicitTendencyComponent
component specifies uses_tracers = True
and sets tracer_dims
, this
feature is enabled.
class MyDynamicalCore(Stepper):
uses_tracers = True
tracer_dims = ['tracer', '*', 'mid_levels']
[...]
tracer_dims
is a list or tuple in the form of a dims
attribute on one of its
inputs, and must have a “tracer” dimension. This dimension refers to which
tracer (you could call it “tracer number”).
Once this feature is enabled, the state
passed to array_call
on the
component will include a quantity called “tracers” with the dimensions
specified by tracer_dims
. It will also be required that these tracers
are used in the output. For a Stepper
component, “tracers”
must be present in the output state, and for a TendencyComponent
or
ImplicitTendencyComponent
component “tracers” must be present in
the tendencies, with the same dimensions as the input “tracers”.
On these latter two components, you should also specify a
tracer_tendency_time_unit
property, which refers to the time part of the
tendency unit. For example, if the input tracer is in units of g m^-3
, and
tracer_tendency_time_unit
is “s”, then the output tendency will be in
units of g m^-3 s^-1
. This value is set as “s” (or seconds) by default.
class MyDynamicalCore(TendencyComponent):
uses_tracers = True
tracer_dims = ['tracer', '*', 'mid_levels']
tracer_tendency_time_unit = 's'
[...]
Memory Management¶
Warning
This section contains fairly advanced topics. If you find it confusing, that’s because the behavior is confusing.
Arrays¶
If possible, you should try to be aware of when there are two code references to the same in-memory array. This can help avoid some common bugs. Let’s start with an example. Say you create a ConstantTendencyComponent object like so:
>>> import numpy as np
>>> from sympl import ConstantTendencyComponent, DataArray
>>> array = DataArray(
np.ones((5, 5, 10)),
dims=('lon', 'lat', 'lev'), attrs={'units': 'K/s'})
>>> tendencies = {'air_temperature': array}
>>> tendency_component = ConstantTendencyComponent(tendencies)
This is all fine so far. But it’s important to know that now array
is the
same array stored inside tendency_component
:
>>> out_tendencies, out_diagnostics = tendency_component({})
>>> out_tendencies['air_temperature'] is array # same place in memory
True
So if you were to modify array
, it would change the output given by
tendency_component:
>>> array[:] = array * 5.
>>> out_tendencies, out_diagnostics = tendency_component({})
>>> out_tendencies['air_temperature'] is array
True
>>> np.all(out_tendencies['air_temperature'].values == array.values)
True
When in doubt, assume that any array you put into a component when it is initialized should not be modified any more, unless changing the values in the component is intentional.
However, this code would not modify the array in tendency_component
:
>>> array = array * 5.
>>> out_tendencies, out_diagnostics = tendency_component({})
>>> out_tendencies['air_temperature'] is array
False
>>> np.all(out_tendencies['air_temperature'].values == array.values)
False
What’s the difference? We took away the [:]
on the left hand side of the
assignment operator. when [:]
is included, python modifies the array on the
left hand side, but when it’s not included it tells the python variable name
“array” to refer to what is on the right hand side. These are subtly different
things - one involves modifying the memory that array
already refers to,
the other involves telling array
to refer to a different place in memory.
More precisely, having array =
tells python
that you want to change what the variable array
refers to, and set it to
be the thing on the right hand side, while array[:] =
tells python to
call the __setitem__(key, value)
method of array
with the contents
of the square parentheses as the key and the right hand side as the value.
Interestingly, array = array * 5.
has different behavior from
array *= 5.
. The first one will change what array
refers to, as before,
while the second one will modify array
in-place without changing the
reference. Writing array *= 5
is the same as writing array[:] = array * 5'
.
All similarly written operations (-=
, +=
, /=
, etc.) are
in-place operations.
Contributing¶
Contributions are welcome, and they are greatly appreciated! Every little bit helps. You can contribute in many ways:
Types of Contributions¶
Usage in Publications¶
If you use Sympl to perform research, your publication is a valuable resource for others looking to learn the ways they can leverage Sympl’s capabilities. If you have used Sympl in a publication, please let us know so we can add it to the list.
Working on projects that use Sympl¶
Sympl is only as useful as the components it has available. You can make Sympl more useful for others by contributing to model projects which use Sympl, or by writing/wrapping model components and deploying them in your own Python packages.
Presenting Sympl to Others¶
Sympl is meant to be an accessible, community-driven tool. You can help the community of users grow and be more effective in many ways, such as:
- Running a workshop
- Offering to be a resource for others to ask questions
- Presenting research that uses Sympl
If you or someone you know is contributing to the Sympl community by presenting it or assisting others with the model, please let us know so we can add that person to the contributors list.
Report Bugs¶
Report bugs at https://github.com/mcgibbon/sympl/issues.
If you are reporting a bug, please include:
- Your operating system name and version.
- Any details about your local setup that might be helpful in troubleshooting.
- Detailed steps to reproduce the bug.
Fix Bugs¶
Look through the GitHub issues for bugs. Anything tagged with “bug” and “help wanted” is open to whoever wants to implement it.
Implement Features¶
Look through the GitHub issues for features. Anything tagged with “enhancement” and “help wanted” is open to whoever wants to implement it.
Write Documentation¶
Sympl could always use more documentation. You could:
- Clean up or add to the official Sympl docs and docstrings.
- Write useful and clear examples that are missing from the examples folder.
- Create a Jupyter notebook that uses Sympl and share it with others.
- Prepare reproducible model scripts to distribute with a paper using Sympl.
- Anything else that communicates useful information about Sympl.
Submit Feedback¶
The best way to send feedback is to file an issue at https://github.com/mcgibbon/sympl/issues.
If you are proposing a feature:
- Explain in detail how it would work.
- Keep the scope as narrow as possible, to make it easier to implement.
- Remember that this is a volunteer-driven project, and that contributions are welcome :)
Get Started!¶
Ready to contribute? Here’s how to set up sympl for local development.
Fork the sympl repo on GitHub.
Clone your fork locally:
$ git clone git@github.com:your_name_here/sympl.git
Install your local copy in development mode:
$ cd sympl/ $ python setup.py develop
Create a branch for local development:
$ git checkout -b name-of-your-bugfix-or-feature
Now you can make your changes locally.
When you’re done making changes, check that your changes pass flake8 and the tests, including testing other Python versions with tox:
$ flake8 sympl tests $ python setup.py test or py.test $ tox
To get flake8 and tox, just pip install them.
Commit your changes and push your branch to GitHub:
$ git add . $ git commit -m "Your detailed description of your changes." $ git push origin name-of-your-bugfix-or-feature
Submit a pull request through the GitHub website.
Pull Request Guidelines¶
Before you submit a pull request, check that it meets these guidelines:
- The pull request should include tests.
- If the pull request adds functionality, the docs should be updated. Put your new functionality into a function with a docstring, and add the feature to the list in README.rst.
- The pull request should work for Python 2.7, 3.4 and 3.5. Check https://travis-ci.org/mcgibbon/sympl/pull_requests and make sure that the tests pass for all supported Python versions.
Style¶
In the Sympl code, we follow PEP 8 style guidelines (tested by flake8). You can test style by running “tox -e flake8” from the root directory of the repository. There are some exceptions to PEP 8:
- All lines should be shorter than 80 characters. However, lines longer than this are permissible if this increases readability (particularly for lines representing complicated equations).
- Space should be assigned around arithmetic operators in a way that maximizes readability. For some cases, this may mean not including whitespace around certain operations to make the separation of terms clearer, e.g. “Cp*T + g*z + Lv*q”.
- While state dictionary keys are full and verbose, within components they may be assigned to shorter names if it makes the code clearer.
- We can take advantage of known scientific abbreviations for quantities within components (e.g. “T” for “air_temperature”) even thought they do not follow pothole_case.
What’s New¶
v0.4.0¶
- Stepper, DiagnosticComponent, ImplicitTendencyComponent, and TendencyComponent base classes were modified to include functionality that was previously in ScalingWrapper, UpdateFrequencyWrapper, and TendencyInDiagnosticsWrapper. The functionality of TendencyInDiagnosticsWrapper is now to be used in Stepper and TendencyStepper objects.
- Composites now have a component_list attribute which contains the components being composited.
- TimeSteppers now have a prognostic_list attribute which contains the prognostics used to calculate tendencies.
- TimeSteppers from sympl can now handle ImplicitTendencyComponent components.
- Added a check for netcdftime having the required objects, to fall back on not using netcdftime when those are missing. This is because most objects are missing in older versions of netcdftime (that come packaged with netCDF4) (closes #23).
- TimeSteppers should now be called with individual Prognostics as args, rather than a list of components, and will emit a warning when lists are given.
- TimeSteppers now have input, output, and diagnostic properties as attributes. These are handled entirely by the base class.
- TimeSteppers now allow you to put tendencies in their diagnostic output. This is done using first-order time differencing.
- Composites now have properties dictionaries.
- Updated basic components to use new component API.
- Components enforce consistency of output from array_call with properties dictionaries, raising ComponentMissingOutputError or ComponentExtraOutputError respectively if outputs do not match.
- Added a priority order of property types for determining which aliases are returned by get_component_aliases.
- Fixed a bug where TendencyStepper objects would modify the arrays passed to them by TendencyComponent objects, leading to unexpected value changes.
- Fixed a bug where constants were missing from the string returned by get_constants_string, particularly any new constants (issue #27)
- Fixed a bug in NetCDFMonitor which led to some aliases being skipped.
- Modified class checking on components so that components which satisfy the component’s API will be recognized as instances using isinstance(obj, Class). Right now this only checks for the presence and lack of presence of component attributes, and correct signature of __call__. Later it may also check properties dictionaries for consistency, or perform other checks.
- Fixed a bug where ABCMeta was not being used in Python 3.
- Added initialize_numpy_arrays_with_properties which creates zero arrays for an output properties dictionary.
- Added reference_air_temperature constant.
- Fixed bug where degrees Celcius or Fahrenheit could not be used as units on inputs because it would lead to an error.
- Added combine_component_properties as a public function.
- Added some unit helper functions (units_are_same, units_are_compatible, is_valid_unit) to public API.
- Added tracer-handling funcitonality to component base classes.
Breaking changes¶
- Implicit, Timestepper, Prognostic, ImplicitPrognostic, and Diagnostic objects have been renamed to TendencyStepper, Stepper, TendencyComponent, ImplicitTendencyComponent, and DiagnosticComponent. These changes are also reflected in subclass names.
- inputs, outputs, diagnostics, and tendencies are no longer attributes of components. In order to get these, you should use e.g. input_properties.keys()
- properties dictionaries are now abstract methods, so subclasses must define them. Previously they defaulted to empty dictionaries.
- Base classes now raise InvalidPropertyDictError when output property units conflict with input property units (which probably indicates that they’re wrong).
- Components should now be written using a new array_call method rather than __call__. __call__ will automatically unwrap DataArrays to numpy arrays to be passed into array_call based on the component’s properties dictionaries, and re-wrap to DataArrays when done.
- TimeSteppers should now be written using a _call method rather than __call__. __call__ wraps _call to provide some base class functionality, like putting tendencies in diagnostics.
- ScalingWrapper, UpdateFrequencyWrapper, and TendencyInDiagnosticsWrapper have been removed. The functionality of these wrappers has been moved to the component base types as methods and initialization options.
- ‘time’ now must be present in the model state dictionary. This is strictly required for calls to DiagnosticComponent, TendencyComponent, ImplicitTendencyComponent, and Stepper components, and may be strictly required in other ways in the future
- Removed everything to do with directional wildcards. Currently ‘*’ is the only wildcard dimension. ‘x’, ‘y’, and ‘z’ refer to their own names only.
- Removed the combine_dimensions function, which wasn’t used anywhere and no longer has much purpose without directional wildcards
- RelaxationTendencyComponent no longer allows caching of equilibrium values or timescale. They must be provided through the input state. This is to ensure proper conversion of dimensions and units.
- Removed ComponentTestBase from package. All of its tests except for output caching are now performed on object initialization or call time.
- “*” matches are now enforced to be the same across all quantities of a component, such that the length of the “*” axis will be the same for all quantities. Any missing dimensions that are present on other quantities will be created and broadcast to achieve this.
- dims_like is obsolete as a result, and is no longer used. dims should be used instead. If present, dims from input properties will be used as default.
- Components will now raise an exception when __call__ of the component base class (e.g. Stepper, TendencyComponent, etc.) if the __init__ method of the base class has not been called, telling the user that the component __init__ method should make a call to the superclass init.
v0.3.2¶
- Exported get_constants_string to the public API
- Added “aliases” kwarg to NetCDFMonitor, allowing the monitor to shorten variable names when writing to netCDF
- Added get_component_aliases() to get a dictionary of quantity aliases from a list of Components (used by NetCDFMonitor to shorten variable names)
- Added tests for NetCDFMonitor aliases and get_component_aliases()
Breaking changes¶
- tendencies in diagnostics are now named as X_tendency_from_Y, instead of tendency_of_X_due_to_Y. The idea is that it’s shorter, and can easily be shortened more by aliasing “tendency” to “tend”
v0.3.1¶
- Fixed botched deployment, see v0.3.0 for the real changes
v0.3.0¶
- Modified component class checking to look at the presence of properties
- Added ScalingWrapper
- Fixed bug in TendencyInDiagnosticsWrapper where tendency_diagnostics_properties were being copied into input_properties
- Modified component class checking to look at the presence of properties attributes instead of checking type when verifying component class.
- Removed Python 3.4 from Travis CI testing
- added some more constants to default_constants related to conductivity of water in all phases and phase changes of water.
- increased the verbosity of the error output on shape mismatch in restore_data_arrays_with_properties
- corrected heat capacity of snow and ice to be floats instead of ints
- Added get_constant function as the way to retrieve constants
- Added ImplicitTendencyComponent as a new component type. It is like a TendencyComponent, but its call signature also requires that a timestep be given.
- Added TimeDifferencingWrapper, which turns an Stepper into an ImplicitTendencyComponent by applying first-order time differencing.
- Added set_condensible_name as a way of changing what condensible aliases (for example, density_of_solid_phase) refer to. Default is ‘water’.
- Moved wrappers to their own file (out from util.py).
- Corrected str representation of DiagnosticComponent to say DiagnosticComponent instead of Stepper.
- Added a function reset_constants to reset the constants library to its initial state.
- Added a function datetime which accepts calendar as a keyword argument, and returns datetimes from netcdftime when non-default calendars are used. The dependency on netcdftime is optional, the other calendars just won’t work if it isn’t installed
- Added a reference to the built-in timedelta for convenience.
Breaking changes¶
- Removed default_constants from the public API, use get_constant and set_constant instead.
- Removed replace_none_with_default. Use get_constant instead.
- set_dimension_names has been removed, use set_direction_names instead.
0.2.1¶
- Fixed value of planetary radius, added specific heat of water vapor.
- Added function set_constant which provides an easy interface for setting values in the default_constants dictionary. Users can already set them manually by creating DataArray objects. This automates the DataArray creation, which should make user code cleaner.
0.2.0¶
- Added some more physical constants.
- Added readthedocs support.
- Overhaul of documentation.
- Docstrings now use numpy style instead of Google style.
- Expanded tests.
- Added function to put prognostic tendencies in diagnostic output.
- NetCDFMonitor is actually working now, and has tests.
- There are now helper functions for automatically extracting required numpy arrays with correct dimensions and units from input state dictionaries. See the note about _properties attributes in Breaking changes below.
- Added base object for testing components
- Renamed set_dimension_names to set_direction_names, set_dimension_names is now deprecated and gives a warning. add_direction_names was added to append to the dimension list instead of replacing it.
Breaking changes¶
- The constant
stefan_boltzmann
is now calledstefan_boltzmann_constant
to maintain consistency with other names. - Removed add_dicts_inplace from public API
- combine_dimensions will raise exceptions in a few more cases where it should do so. Particularly, if there is an extra dimension in the arrays.
- Default out_dims is removed from combine_dimensions.
- input_properties, tendency_properties, etc. dictionaries have been added to components, which contain information about the units and dimensions required for those arrays, and can include more properties as required by individual projects. This makes it possible to extract appropriate numpy arrays from a model state in an automated fashion based on these properties, significantly reducing boilerplate code. These dictionaries need to be defined by subclasses, instead of the old “inputs”, “outputs” etc. lists which are auto-generated from these new dictionaries.
- Class wrapping now works by inheritance, instead of by monkey patching methods.
- All Exception classes (e.g. SharedKeyException) have been renamed to “Error” classes (e.g. SharedKeyError) to be consistent with normal Python naming conventions
0.1.1 (2017-01-05)¶
- First release on PyPI.
Credits¶
Development Lead¶
- Jeremy McGibbon <mcgibbon@uw.edu>
Contributors¶
- Joy Monteiro <joy.monteiro@misu.su.se>
License¶
sympl is available under the open source BSD License.