Extending mass-action semantics, Part 2

Building a categorical model in CatColab

CatColab
modeling
Author
Published

2026-04-18

Abstract

In the prequel to this post, we saw some benefits of taking a categorical approach to modelling, but also how the process can be fiddly and error-prone, even for a small model. In this short post we describe the implementation of the extended mass-action semantics in CatColab and revisit the model from last time.

Note

This blog post is a direct continuation of Part 1.

We left off last time saying how nice it would be if CatColab supported a generalised version of mass-action semantics for Petri nets. To cut straight to the punchline, as of the recent CatColab v0.5 release, Petri nets (as well as stock-flow diagrams) now support this “unbalanced” mass-action dynamics. In the analyses for mass-action dynamics, there is a setting menu with an option to toggle conserve mass which, if unticked, will allow you to pick a rate granularity of either per place or per transition. The help pages describe this in much more detail.

If you just want to see the results, I’ve re-implemented the model from Part 1. Below you can view both the model and the derived equations, or you can play around with them yourself at this link. The key thing to note here is that the equations are derived from the model, which means that you can add or remove transitions and immediately see how the equations change. This solves the main problem that I pointed out last time: how easy it is to make a mistake when translating a Petri net into the mass-action equations by hand.

Before proceeding, I want to point out one thing that I find particularly exciting. Although CatColab is still under very active development, it has now reached the stage where I could to go from idea, to research, to implementation, to having a live public-facing feature in a matter of weeks. That’s because there is now a solid foundation of tooling and examples in the codebase for people to build on, in terms of both foundational categorical logic and UI components and general software infrastructure. And of course this is only the case due to a gargantuan effort from all the people working on CatColab.

GitHub tells me that my implementation of unbalanced mass-action semantics consists of 1,059 lines of code added, and 927 lines of code deleted. If we pretend that positive code and negative code cancel out, then that’s under 150 lines of code — not too bad! The Petri net (and stock-flow) help pages on CatColab contain documentation for this new semantics.

For those interested in a small walkthrough of some of the key parts of the implementation, found in mass_action.rs, you can read more below!

1 Defining some types

Let’s begin by writing down the different types of mass-action semantics that we want to differentiate between. We are interested in three types, each more expressive than the last:

  • Balanced. Transitions conserve mass, and are described by a single transition rate parameter r_T. So a transition [A,B]\to\boxed{T}\to[X,Y] requires one parameter r_T, and gives the equations \left\{ \begin{aligned} \dot{A} &= -r_T AB \\\dot{B} &= -r_T AB \\\dot{X} &= \phantom{-}r_T AB \\\dot{Y} &= \phantom{-}r_T AB \end{aligned} \right.

  • Unbalanced, per transition. Transitions do not necessarily conserve mass, and are described a consumption rate parameter \kappa_T and a production rate parameter \rho_T. So a transition [A,B]\to\boxed{T}\to[X,Y] requires two parameters (\kappa_T,\rho_T), and gives the equations \left\{ \begin{aligned} \dot{A} &= -\kappa_T AB \\\dot{B} &= -\kappa_T AB \\\dot{X} &= \phantom{-}\rho_T AB \\\dot{Y} &= \phantom{-}\rho_T AB \end{aligned} \right.

  • Unbalanced, per place. Transitions do not necessarily conserve mass, and are described a consumption rate parameter \kappa_T^A for each input place A and a production rate parameter \rho_T^X for each output place X. So a transition [A,B]\to\boxed{T}\to[X,Y] requires four parameters (\kappa_T^A,\kappa_T^B,\rho_T^X,\rho_T^Y), and gives the equations \left\{ \begin{aligned} \dot{A} &= -\kappa_T^A AB \\\dot{B} &= -\kappa_T^B AB \\\dot{X} &= \phantom{-}\rho_T^X AB \\\dot{Y} &= \phantom{-}\rho_T^Y AB \end{aligned} \right.

Note that these semantics make sense for Petri nets, but the first two also make sense for stock-flow diagrams. Because of this, we sometimes switch between the words “transition” and “flow”.

We can write down a simple description of this hierarchy in Rust types. First off, we specify the three types of mass-action, using RateGranularity to describe the distinction between the per-place and per-transition cases.

pub enum MassConservationType {
    Balanced,
    Unbalanced(RateGranularity),
}

pub enum RateGranularity {
    PerTransition,
    PerPlace,
}

Next we describe what the parameters will look like for each flow. In CatColab, QualifiedName is the type of a name referring to an object or morphism declared in a model. (It is called “qualifed” because a name may consist of multiple segments when a model is constructed by hierarchical composition.) In particular, every place and transition in a Petri net is identified by a QualifiedName.

Recall that, in the balanced case, we don’t need to worry about the direction of flows (e.g. whether a place is an input or an output to the transition); in the unbalanced per-transition case, we don’t need to worry about distinguishing between any of the input places (resp. output places); in the unbalanced per-place case, we need to worry about everything! Note that the Direction type has the opposite terminology from what we might expect: a transition A\to\boxed{T}\to B (or a flow A\Rightarrow B) gives rise to an incoming flow to the output B and an outgoing flow from the input A.

pub enum FlowParameter {
    Balanced {
        transition: QualifiedName,
    },
    Unbalanced {
        direction: Direction,
        parameter: RateParameter,
    },
}

pub enum RateParameter {
    PerTransition {
        transition: QualifiedName,
    },
    PerPlace {
        transition: QualifiedName,
        place: QualifiedName,
    },
}

pub enum Direction {
    IncomingFlow,
    OutgoingFlow,
}

2 Building the system of equations

With the types all set up, the process of actually building the equations from a given model is rather routine, and consists of matching against the mass-action type and then using all of the built-in functionality of catlog (the core package of CatColab).

The first step is to simply create an equation for each object (i.e. place or stock) of the form \dot{A}=0, so that we can add all the contributions to it when we later iterate over the morphisms (i.e. transitions or flows).

let mut sys = PolynomialSystem::new();

for ob in model.ob_generators_with_type(&self.stock_ob_type) {
    sys.add_term(ob, Polynomial::zero());
}

Next we’ll iterate over all the pairs (flow, term), where flow is a morphism in the model and term is the monomial built by multiplying together all of the variables corresponding to the input places (or, in the case of stock-flow diagrams, the variable corresponding to the input stock and those of any input links). For this we use a helper function flow_monomials which we won’t look at in any detail. As we said above, all we need to do here is match against MassConservationType, so let’s start by just sketching out the shape of the iteration.

let terms: Vec<_> = self.flow_monomials(model).into_iter().collect();

for (flow, term) in terms {
    let dom = model.mor_generator_dom(&flow).unwrap_basic();
    let cod = model.mor_generator_cod(&flow).unwrap_basic();
    match mass_conservation_type {
        MassConservationType::Balanced => {
            ...
        }
        MassConservationType::Unbalanced(granularity) => {
            ...
        }
    }
}

For the case of MassConservationType::Balanced, we simply need to build the monomial given by multiplying term by a single rate parameter, i.e. of of type FlowParameter::Balanced.

MassConservationType::Balanced => {
    let term: Polynomial<_, _, _> = [(
        Parameter::generator(FlowParameter::Balanced { transition: mor }),
        term.clone(),
    )]
    .into_iter()
    .collect();

    for input in inputs {
        sys.add_term(input.unwrap_generator(), -term.clone());
    }

    for output in outputs {
        sys.add_term(output.unwrap_generator(), term.clone());
    }
}

For MassConservationType::Unbalanced(_), we need to treat the per-transition and per-place cases differently. In the former, our parameter will be of type RateParameter::PerTransition, and in the latter of type RateParameter::PerPlace. To avoid repetition, we’ll just look at iterating over the input places, since the code for output places is essentially identical: we just change Direction from OutgoingFlow to IncomingFlow and add instead of adding -input_term to our system we add output_term.

MassConservationType::Unbalanced(granularity) => {
    for input in inputs {
        let input_term: Polynomial<_, _, _> = match granularity {
            RateGranularity::PerTransition => [(
                Parameter::generator(FlowParameter::Unbalanced {
                    direction: Direction::OutgoingFlow,
                    parameter: RateParameter::PerTransition {
                        transition: mor.clone(),
                    },
                }),
                term.clone(),
            )],
            RateGranularity::PerPlace => [(
                Parameter::generator(FlowParameter::Unbalanced {
                    direction: Direction::OutgoingFlow,
                    parameter: RateParameter::PerPlace {
                        transition: mor.clone(),
                        place: input.clone().unwrap_generator(),
                    },
                }),
                term.clone(),
            )],
        }
        .into_iter()
        .collect();

        sys.add_term(input.unwrap_generator(), -input_term.clone());
    }
    for output in outputs {
      ...
    }

    sys.normalize()
}

3 Writing a test

There is a growing collection of test models in catlog that we use for unit and regression testing. One of them is a small Petri net representing a catalytic reaction:

Thanks to all the work on the DoubleTT text elaborator, we can construct this model without too much boilerplate:

let th = Rc::new(th_sym_monoidal_category());
let model = tt::modelgen::Model::from_text(
    &th.into(),
    "[
        x : Object,
        y : Object,
        c : Object,
        f : (Hom Object)[@tensor [x, c], @tensor [y, c]],
    ]",
);
let model = model.unwrap().as_modal().unwrap();

We can now write a test that applies unbalanced, per-place, mass-action semantics to this model and verifies that the output is what we expect, namely \left\{ \begin{aligned} \dot{x} &= -\kappa_f^x cx \\\dot{y} &= \phantom{-}\rho_r^y cx \\\dot{c} &= (\rho_f^c - \kappa_f^c) cx \end{aligned} \right. Note that the “catalyst” c is not actually left unchanged unless f is balanced with respect to c, i.e. unless \rho_f^c = \kappa_f^c.

#[test]
fn catalysis_dynamics() {
    let th = Rc::new(th_sym_monoidal_category());
    let model = catalyzed_reaction(th);
    let sys = PetriNetMassActionAnalysis::default().build_system(
        &model,
        MassConservationType::Unbalanced(RateGranularity::PerPlace),
    );
    let expected = expect!([r#"
        dx = -(x->[f]) c x
        dy = ([f]->y) c x
        dc = (([f]->c) - (c->[f])) c x
    "#]);
    expected.assert_eq(&sys.to_string());
}

Here, instead of LaTeX, we’re just using a simple ASCII printer for parameters, which outputs e.g. (x->[f]) instead of \kappa_f^x and ([f]->y) instead of \rho_f^y.

Running cargo test tells us that everything works just as we want. Nice!

Leaving a comment will set a cookie in your browser. For more information, see our cookies policy.