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:

DataArray

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

class sympl.timedelta[source]

Difference between two datetime values.

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

sympl.set_constant(name, value, units)[source]

Sets the value of a constant.

Parameters:
  • name (str) – The name of the constant.
  • value (float) – The value to which the constant should be set.
  • units (str) – The units of the value given.
sympl.reset_constants()[source]

Reverts constants to their state when Sympl was originally imported. This includes removing any new constants, setting the original constants to their original values, and setting the condensible quantity to water.

Debugging and Logging Constants

You can get a string describing current constants by calling get_constants_string():

import sympl
print(sympl.get_constants_string())
sympl.get_constants_string()[source]
Returns:constant_string – A string listing all constants under each category of constants, with their current values and units.
Return type:str

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

sympl.set_condensible_name(name)[source]

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.
__repr__()[source]

Return repr(self).

__str__()[source]

Return str(self).

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.
__repr__()[source]

Return repr(self).

__str__()[source]

Return str(self).

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.
__init__()[source]

Initializes the Stepper object.

__repr__()[source]

Return repr(self).

__str__()[source]

Return str(self).

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.
__repr__()[source]

Return repr(self).

__str__()[source]

Return str(self).

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.
__repr__()[source]

Return repr(self).

__str__()[source]

Return str(self).

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.Monitor[source]
__repr__()[source]

Return repr(self).

__str__()[source]

Return str(self).

store(state)[source]

Stores the given state in the Monitor and performs class-specific actions.

Parameters:state (dict) – A model state dictionary.
Raises:InvalidStateError – If state is not a valid input for the DiagnosticComponent instance.
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.
write()[source]

Write all cached states to the NetCDF file, and clear the cache. This will append to any existing NetCDF file.

Raises:InvalidStateError – If cached states do not all have the same quantities as every other cached and written state.
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.
store(state)[source]

Updates the plot using the given state.

Parameters:state (dict) – A model state dictionary.

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

class sympl.MonitorComposite(*args)[source]
store(state)[source]

Stores the given state in the Monitor and performs class-specific actions.

Parameters:

state (dict) – A model state dictionary.

Raises:
  • KeyError – If a required quantity is missing from the state.
  • InvalidStateError – If state is not a valid input for a Monitor instance.

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

sympl.get_tracer_names()[source]
Returns:tracer_names – Tracer names in the order that they will appear in tracer arrays.
Return type:tuple of str
sympl.get_tracer_unit_dict()[source]
Returns:unit_dict – A dictionary whose keys are tracer quantity names as str and values are units of those quantities as str.
Return type: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.

sympl.units_are_same(unit1, unit2)[source]

Compare two unit strings for equality.

Parameters:
  • unit1 (str) –
  • unit2 (str) –
Returns:

units_are_same – True if the two input unit strings represent the same unit.

Return type:

bool

sympl.units_are_compatible(unit1, unit2)[source]

Determine whether a unit can be converted to another unit.

Parameters:
  • unit1 (str) –
  • unit2 (str) –
Returns:

units_are_compatible – True if the first unit can be converted to the second unit.

Return type:

bool

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.get_numpy_arrays_with_properties(state, property_dictionary)[source]
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.

  1. Fork the sympl repo on GitHub.

  2. Clone your fork locally:

    $ git clone git@github.com:your_name_here/sympl.git
    
  3. Install your local copy in development mode:

    $ cd sympl/
    $ python setup.py develop
    
  4. Create a branch for local development:

    $ git checkout -b name-of-your-bugfix-or-feature
    

    Now you can make your changes locally.

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

  6. 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
    
  7. Submit a pull request through the GitHub website.

Pull Request Guidelines

Before you submit a pull request, check that it meets these guidelines:

  1. The pull request should include tests.
  2. 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.
  3. 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.

Tips

To run a subset of tests:

$ py.test tests.test_timestepping

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 called stefan_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

Contributors

License

sympl is available under the open source BSD License.