Welcome to `musica`!#

This tutorial will walk you through getting an installation working and running a box model using the musica python package.

NOTE: We’ll use musica to refer to the python package developed as part of the larger MUlti-Scale Infrastructure for Chemistry and Aerosols (MUSICA) Project. To learn more about the MUSICA Project, see here.

Installation (optional)#

You probably arrived at this tutorial through binder. If so, please skip this section as we’ve taken care of installing everything for you.

If you want to install musica locally, it is available from pypi. Below are instructions that will install musica’s dependencies alongside all of the dependencies requires for this tutorial

pip install musica[tutorial]

We do recommend you install into some virtual environment. We’ve included instructions for both `uv <https://docs.astral.sh/uv/>`__, and `conda <https://docs.conda.io/projects/conda/en/latest/user-guide/getting-started.html>`__

uv (recommended)#

uv python install 3.13
uv venv --python 3.13

source .venv/bin/activate   # macOS/Linux
.venv\Scripts\activate      # Windows

uv pip install musica[tutorial] # bash
uv pip install 'musica[tutorial]' # zsh

Conda#

conda create --name musica python=3.13
conda activate musica
pip install musica[tutorial] # bash
pip install 'musica[tutorial]' # zsh

CLI#

The musica pacakge comes with a cli tool, musica-cli. This can be useful to copy out some of our examples, or printing the version. You can read more about its usage and options here. To verify your installation worked, you can ask it to print the version number

musica-cli --version

You should see something like

musica 0.15.0 (MICM 3.12.0, TUV-x 0.15.0, CARMA 4.0.0)

Using `musica`#

Let’s start off by importing musica and getting the versions in python.

Version checks#

[ ]:

import musica

[ ]:

print(musica.__version__)

Each of musica’s componenents can be imported as a separate package. There are currently three:

micm: ODE solvers and gas-phase chemistry package.
tuv-x**: Photolysis rate constant calculator.
carma**: Aerosol model.

NOTE: The unreliability of fortran compilers prevents inclusion of tuv-x and carma in musica on some systems (including Windows). Some of the cells below will fail if you are on one these systems. As legacy code is ported to modern languages, we expect to be able to include all musica components in all deployments.

[ ]:

print(musica.micm.__version__)

[ ]:

print(musica.tuvx.__version__)

[ ]:

print(musica.carma.__version__)

To run chemistry you only need micm, which is our ODE solver. There are separate tutorials for tuv-x and carma so we will forget about them for now.

A simple gas-phase system from a configuration#

There are two ways to configure musica:

using a valid mechanism configuration in json/yaml
using the musica API.

We’ll use the configuration file approach in this first tutorial chapter. In the next chapter we’ll present the in-code approach, which we’ll use in the remaining examples.

We’ve provided several packaged example configurations with musica to make their usage easy. We will run a simple 3-species system below. We’ll start by importing the required parts of musica.

musica, this will give us access to all parts of musica, including MICM, our gas-phase solver
Parser, which can read our mechanism files
find_config_path, a utility function that will make it easy to load our example

[ ]:

import musica
from musica.mechanism_configuration import Parser
from musica.utils import find_config_path

Before we build our model, let’s take a look at the configuration file to see how it’s organized

[ ]:

import json

with open(find_config_path("v1", "analytical", "config.json")) as f:
    config = json.load(f)

print(json.dumps(config, indent=2))

This simple mechanism comprises:

Three species (A, B, and C) in the gas-phase
Two Arrhenius reactions
- A -> B
- B -> C

The format of the musica configuration is intended to describe the science without introducing details of how musica software is implemented in code.

Now, let’s build a real model

[ ]:

parser = Parser()
mechanism = parser.parse(find_config_path("v1", "analytical", "config.json"))

[ ]:

for reaction in mechanism.reactions:
    print(f"[{reaction.type.name}] {reaction.to_equation()}")

Above, we had our parser read the mechanism file and return a valid mechanism. We can see the same reactions from the configuration file, now part of a mechanism instance.

We can now use this to create a MICM solver and a state.

[ ]:

solver = musica.MICM(mechanism=mechanism)
state = solver.create_state()

Why do we need a solver and a state?

Conceptually, they fill two distinct roles:

The state is a container for time-varying properties for your system of interest (temperature, pressure, species concentrations, etc.)
The solver is the logic needed to advance a state in time. In our case we get the default MICM solver, which is a Rosenbrock 3-stage solver with a predetermined set of solver parameters. We could get a different ODE solver, or set different solver parameters by including additional arguments to the MICM() constructor.

Practically, having a separate solver and state allows you to create as many state containers as you want. Each state can all be advanced in time using the same solver instance. Each state instance can also be configured to represent multiple grid cells, which you’ll read more about in the next tutorial.

For now, we’ll keep things simple with a single one-grid-cell state. Let’s set the initial conditions and move the chemistry forward!

[ ]:

state.set_conditions(temperatures=[298.15], pressures=[101325])
state.set_concentrations(
    {
        'A' : [1.0],
        'B' : [0.0]
    }
)

Now all that’s left is to solve and collect the concentrations as they change over time. For now, we’ll print them out. Later tutorials have more code to make some pretty graphs.

solver.solve takes in a state and timestep to solve for. Let’s solve for 10 seconds and see how things change

[ ]:

def print_concentrations(t, state):
    cs = state.get_concentrations()
    print(f"{t}\t| {cs['A'][0]:.4f}\t| {cs['B'][0]:.4f}")

print_concentrations(0, state)
for t in range(1, 11):
    solver.solve(state, 1)
    print_concentrations(t, state)

And that’s your first MUSICA box model! Read on to learn how to define mechanisms in-code for solving and to do more some more complex tasks.

Welcome to musica!#