Modeling 1: your first model

In about fifteen minutes you write a process model as a TOML file, let Labinetix check it, generate its embedded kernel, and simulate it on the host. The model is a batch culture: biomass grows on one substrate with Monod kinetics. Nothing here talks to hardware.

This is the first of four tutorials:

1. Your first model (this page): author, check, generate, simulate.
2. Composition and balances: build a model from species and reactions, and let the elemental balances check it.
3. Estimation and identifiability: ask what data can tell you about the parameters, and estimate them.
4. From a model to a released package: the model’s generated package under taktwerk, bit-identical to the prediction, as an FMU, released, and back from its record.

Before you start

Do Getting started step 1 so that labinetix is on your PATH, and work from the root of your clone:

cd "$REPO_ROOT"
export PATH="$PWD/target/debug:$PATH"

Two directories matter here. models/ holds every model, one directory per model: models/<id>/authoring/<id>.toml is the file you write, models/<id>/generated/ is what Labinetix writes from it. models/ is the store the define commands work on (--store, default models, so commands run from the root of the clone need no flag).

1. Create the model

labinetix define new model_scaffold monod_batch

This writes models/monod_batch/authoring/monod_batch.toml with a placeholder model (a first-order lag), and generates it once. Open the file and replace its content with the batch culture:

# model_scaffold/monod_batch — a batch culture: biomass grows on one substrate with Monod
# kinetics and a constant yield, in a vessel that is neither fed nor drained.

[model]
name = "monod_batch"
contract_version = "0.1.0"
kind = "ode"
time_unit = "h"
doc = "Batch Monod growth: dx/dt = mu x, ds/dt = -mu x / y_xs, mu = mu_max s / (k_s + s)."

[parameters.mu_max]
unit = "1/h"
min = 0.0
default = 0.5
doc = "Maximum specific growth rate."

[parameters.k_s]
unit = "g/L"
min = 0.0
default = 0.1
doc = "Monod constant of the substrate."

[parameters.y_xs]
unit = "g/g"
min = 0.0
default = 0.5
doc = "Biomass yield on substrate."

[states.x]
unit = "g/L"
min = 0.0
doc = "Biomass concentration."

[states.s]
unit = "g/L"
min = 0.0
doc = "Substrate concentration."

[[locals]]
name = "mu"
unit = "1/h"
expr = "mu_max * monod(s, k_s)"
doc = "Specific growth rate."

[equations]
x = "mu * x"
s = "-mu * x / y_xs"

[outputs.biomass]
unit = "g/L"
expr = "x"
doc = "Biomass concentration, as measured."

[outputs.growth_rate]
unit = "1/h"
expr = "mu"
doc = "Specific growth rate."

What the tables mean:

Declaration order matters: it is the order of the state vector and of the generated interface.

Units are checked: define parses every unit (g/L, 1/h, mmol/(1e8 cells d); g_X is g, the tag names the species) and refuses a sum, a min/max, a comparison or a right-hand side whose units disagree in dimension or scale, naming the entry and both units. Each equation is in the state’s unit per time_unit. A bare number takes the unit its context needs (s + 1e-9); a symbol without a unit is unchecked.

Expressions use a small operator set that every embedded target can run: + - * /, pow, min, max, exp, log, sqrt, abs and ifelse(condition, a, b). Kinetic laws such as monod(c, k), haldane(c, k, ki), contois(s, x, k), moser(c, k, n), teissier(c, k), blackman(c, k) and inhibition(c, ki) are written as calls and expanded into that set: monod(s, k_s) is exactly s / (k_s + s).

2. Generate

labinetix define sync model_scaffold/monod_batch

sync parses and validates the model, and writes:

• models/monod_batch/generated/ir/monod_batch.json, the model as a canonical expression tree. Its hash is the model’s content hash: everything that pins this model pins that hash.
• models/monod_batch/generated/julia/monod_batch.jl, a ModelingToolkit view, as data.
• crates/kernels/core-models/generated/monod_batch.rs, the kernel: no_std Rust that computes the right-hand side, the outputs, and their exact Jacobians in the states, inputs and parameters, derived symbolically. This is the code that runs on a target.

It also adds the model to the kernels’ index, so the kernel compiles with every other one. sync commits what it wrote (sync(model_scaffold/monod_batch): regenerate views). The kernel’s right-hand side reads like the model:

#![allow(unused)]
fn main() {
fn rhs(&self, x: &[f64; NX], _u: &[f64; NU], p: &[f64; NP], _t: f64, dx: &mut [f64; NX]) {
    let l_mu = p[0] * (x[1] / (p[1] + x[1]));
    dx[0] = l_mu * x[0];
    dx[1] = -l_mu * x[0] / p[2];
}
}

3. When the model is wrong

Labinetix refuses a model it cannot run, before anything is generated, and names the file, the table and the symbol. Misspell the yield in the substrate equation (yxs for y_xs) and sync:

labinetix define sync model_scaffold/monod_batch

labinetix: define sync: invalid entity state: …/models/monod_batch/authoring/monod_batch.toml: [equations] `s`: references `yxs`, which no table declares

The same holds for an operator outside the set, a local used before it is defined, a state without an equation, or bounds that contradict each other. Fix the name and sync again.

define check verifies the whole store: every generated file is what its source generates today, and every pin still matches. It prints nothing when all is well. Run it before you commit by hand:

labinetix define check

4. A scenario

A scenario says how to run the model: the solver, the initial state, input values and changes, and parameter values that differ from the defaults. It lives beside the model and pins the model’s content hash, so it can never silently run a changed model. Find the hash in the kernel (CONTENT_HASH) or the model’s lock file, and write models/monod_batch/scenarios/batch.toml:

[scenario]
model = "monod_batch"
content_hash = "sha256:…"      # the model's content hash
doc = "10 g/L substrate, 0.1 g/L biomass, 12 h."

[solver]
method = "rk4"                 # rk4 or euler: what a target runs; ros2 for stiff models
step = 0.1                     # model time per step (h)
substeps = 4
steps = 120
record_every = 10

[initial]
x = 0.1
s = 10.0

Inputs go in [inputs] (their value at the start), changes in [[schedule]] (step, input, value), overrides in [parameters]. A discrete model’s solver is method = "discrete" with step = 1 and substeps = 1: one step of the model per step, so a row’s t is the step index.

5. Simulate

The simulator steps the compiled kernel, so rebuild labinetix after every define sync that changed a model:

cargo build -p labinetix
labinetix analyze simulate models/monod_batch/scenarios/batch.toml

{
  "kind": "labinetix_model_prediction",
  "model": "monod_batch",
  "content_hash": "sha256:8cd6…",
  "scenario": "batch",
  "scenario_hash": "sha256:b8cb…",
  "columns": {"x": ["x","s"], "u": [], "y": ["biomass","growth_rate"]},
  "rows": [
    {"t": 0.0, "x": [0.1,10.0], "u": [], "y": [0.1,0.49504950495049505]},
    {"t": 1.0, "x": [0.16405317914198603,9.871893641716031], "u": [], "y": [0.16405317914198603,0.4949859072111608]},
    …
    {"t": 12.0, "x": [5.100000000000007,3.1297485988535e-84], "u": [], "y": [5.100000000000007,1.56487429942675e-83]}
  ]
}

The culture grows at nearly mu_max until the substrate runs out, around 9 h, and ends at 0.1 + 0.5 · 10 = 5.1 g/L biomass. The simulation steps the generated kernel through the same integrator a target uses, so this is what the target computes.

Save it beside the scenario as batch.prediction.json (--out) and commit both: the test suite re-runs every committed scenario and requires the same bytes, so a change to the model or the toolchain that moves a prediction cannot go unnoticed.

labinetix analyze simulate models/monod_batch/scenarios/batch.toml \
  --out models/monod_batch/scenarios/batch.prediction.json

6. Change the model

Edit the authoring, raise contract_version as the change deserves (a new parameter appended last is a minor change, a removed output a major one), sync and rebuild. Even a new parameter default moves the content hash, so the scenario refuses to run until you re-pin it, after you have looked at what changed:

labinetix analyze simulate models/monod_batch/scenarios/batch.toml

labinetix: scenario `batch` pins sha256:8cd6…, and model `monod_batch` is sha256:ad8f…: re-pin it after reviewing the change

7. Models ported from a paper

Some models in models/ are written without composition, as this one is, and reproduce numbers their source prints. Each header cites the table and page of every default and names any deviation from the printed equations; a test in crates/studio/analyze-model/tests/ checks the committed prediction against a second implementation and against the printed numbers:

• models/hybridoma_fed_batch: hybridoma cells on glucose and glutamine with lactate, ammonia, death and antibody (de Tremblay et al. 1992, re-simulated by Roubos et al. 1997). Its constant-feed scenario ends at Roubos’s 280.4 mg of antibody, within 0.2 %.
• models/penicillin_fed_batch: the penicillin benchmark of Lim et al. with Contois growth and Haldane-inhibited production (Ochoa 2016). Its luus-20-stages scenario replays Luus’s optimal feed (1993) and meets every state of his Table II to 2e-4 and his 87.95 g of penicillin.
• models/zymomonas_fed_batch: Z. mobilis ethanol fed-batch with threshold and power-law inhibition, in which growth consumes no glucose (Wang and Shyu 1997). Its table-1-idp-sqp scenario replays the paper’s optimal feed to 78.635 kg of ethanol against the printed 78.6302 kg; its initial biomass is fitted to the paper’s Fig. 1 and labelled so.

Next

Composition and balances builds a model from species and reactions instead of writing every balance by hand, and lets Labinetix check that the reactions conserve atoms and electrons.