# Process¶

**Prepared by:****Covered topics:**0. Recap

1. Understanding Process/Processes/processes

2. Making a simple biokinetic model

3. Advanced features

**Video demo:**

To run tutorials in your browser, go to this Binder page.

You can also watch a video demo on YouTube (part 1, part 2) (subscriptions & likes appreciated!).

## 0. Recap¶

Check out this video here about `qsdsan`

, make sure you understand the UML (unified modeling language) diagram.

```
[1]:
```

```
from IPython.display import Image
Image(url='https://lucid.app/publicSegments/view/c8de361f-7292-47e3-8870-d6f668691728/image.png', width=800)
```

```
[1]:
```

This tutorial will focus on the *Process* module in `qsdsan`

. This module was developed to provide modeling capacity for dynamic conversions of waste streams in a `SanUnit`

. With this module, simulations can describe not only a system’s mass and energy flows at steady states but also the changes of these flows over time.

```
[2]:
```

```
import qsdsan as qs
print(f'This tutorial was made with qsdsan v{qs.__version__}.')
```

Back to top

## 1. Understanding `Process`

/`Processes`

/`processes`

¶

Similar to the Component module, the *Process* module includes three main classes: `Process`

, `Processes`

, and `CompiledProcesses`

.

The `Process`

class is used to model individual dynamic processes with user-defined stoichiometry and reaction kinetics. The `Processes`

class is simply a collection of multiple `Process`

objects. Compilation of `Processes`

yields `CompiledProcesses`

, which provides additional attributes of the collection. It is often used to implement (pseudo-)mechanistic models with parellel biological and/or physiochemical processes in dynamic simulations of treatment systems.

```
[3]:
```

```
Image(url='https://lucid.app/publicSegments/view/2c231fa2-6065-46b9-83af-a790ce38b6c0/image.png', width=600)
```

In this section, we will focus on understanding the structures of these classes by going over some existing processes in `qsdsan`

. Next section will guide you through the creation of your own `Processes`

subclass.

### 1.1. `Processes`

vs. `processes`

¶

You might remember `SanUnit`

vs. `sanunits`

. `Processes`

and `processes`

have similar distinctions.

```
[4]:
```

```
# If you check
qs.Process
# or
# qs.Processes
# or
# qs.CompiledProcesses
```

This means `Process`

/`Processes`

/`CompiledProcesses`

is a class in the `_process`

module (i.e., the codes defining theses classes are in the `_process.py`

file within the `qsdsan`

folder).

```
[5]:
```

```
# If you check
qs.processes
```

This shows that `processes`

is actually a folder. (Why? See here.) It contains the process models (i.e., subclasses of `Process`

or `CompiledProcesses`

) that have been included in `qsdsan`

.

```
[6]:
```

```
# To see the list of objects that can be directly imported from this folder
qs.processes.__all__
```

The list above include common process models like `ASM1`

, `ASM2d`

, `ADM1`

, which are implemented as subclasses of `CompiledProcesses`

. You can also directly import utility functions that create `CompiledComponents`

objects corresponding to the state variables in those process models (i.e., `create_<ModelName>_cmps`

) among others.

`Decay`

and `KineticReaction`

are classes for simple processes whose state variables can be solved analytically.

```
[7]:
```

```
# You can see other attributes of the `qs.processes` folder with the `dir` function
set(dir(qs.processes)) - set(qs.processes.__all__)
```

While names with double “_” on both sides (e.g., `__path__`

) are built-in Python attributes, those with a single leading “_” (e.g., `_asm1`

) are individual scripts included in the folder.

```
[8]:
```

```
# For example, `_asm1.py` is a script in the `processes` folder.
qs.processes._asm1
```

For the rest of this section, we will use activated sludge model no. 1 (ASM1) as an example to learn about the core classes in the `_process`

module.

Back to top

#### ASM1¶

ASM1 is a mathematical model for the activated sludge processes in wastewater treatment systems. It uses **13** components to characterize wastewater and activated sludge, including biomass, organic substrates (COD), organic and inorganic nitrogens (N), dissolved oxygen (DO), alkalinity etc.

ASM1 includes **8** transformation processes of the components (e.g., hydrolysis, biomass growth and decay). The stoichiometry of each process considers conservations of COD, N, and electric charges. The process rates are dependent on environmental conditions and state variables (i.e., relevant component concentrations).

```
[9]:
```

```
# Before we get to the subsections, let's get ready by loading ASM1-related objects in qsdsan
from qsdsan.processes import create_asm1_cmps, ASM1
cmps = create_asm1_cmps()
cmps.show()
```

This means the `create_asm1_cmps`

function returns a `CompiledComponents`

object corresponding to the 13 components in ASM1 besides water.

Note:Classes in the`Process`

module cannot be used independent of the`Component`

module. Therefore, it is common practice to create corresponding`CompiledComponents`

objects and set thermo before creating`Process`

instances.

```
[10]:
```

```
# By default, the thermo is set with this `CompiledComponents` object upon its creation.
# We can verify that by calling the `get_thermo` function
qs.get_thermo()
```

```
[11]:
```

```
# Next we need to create an instance of the ASM1 model
# We can see that `ASM1` is a subclass of `CompiledProcesses`, so it can be used for demonstration
# of the capabilities of the `_process` module.
ASM1.__bases__
```

```
[12]:
```

```
# Let's start by reading the documentation of the `ASM1` class
?ASM1
```

```
[13]:
```

```
# Without getting into the details of ASM1, we will leave all parameters at their default values
asm1 = ASM1()
asm1.show()
```

This shows that `asm1`

is composed of the 8 transformation processes.

Back to top

### 1.2. `Process`

¶

The `Process`

class is developed to model individual *dynamic* processes. *Dynamic* means that the rate of reaction is finite (i.e., the reaction is not instantaneous) and often changes with the state of the system. Each of the eight transformation process in ASM1 is described by a `Process`

object.

```
[14]:
```

```
# Let's take the 0th process in `asm1` as an example to learn about `Process`:
# p1 = asm1.tuple[0]
# You can also access the process using the typical `getattr` method because
# individual processes in a `CompiledProcesses` object are stored as attributes
p1 = asm1.aero_growth_hetero
p1.show()
```

### Tip¶

`<Object>.<Attribute>`

is essentially identical to `getattr(<Object>, '<AttributeName>')`

.

`p1`

is an instance of `Process`

describing the growth of heterotrophic biomass under aerobic conditions. The `show`

method prints `p1`

’s defining characteristics: stoichiometry, reference, rate equation, parameters, and dynamic parameters. All of them are stored as the object’s attributes.

```
[15]:
```

```
# For example, we can retrieve information on the stoichiometry of this process
p1.stoichiometry
```

Unlike what’s printed by the `show`

method, a dictionary of evaluated stoichiometric coefficients is returned. The coefficients indicate that to produce 1g COD’s worth of `X_BH`

, approximately 1.49g COD’s worth of `S_S`

, 0.08g N’s worth of `S_NH`

, and 0.49g `S_O`

are consumed. The values of these stoichiometric coefficients are dependent on the parameter values in the formulas. We’ll see how stoichiometry will change in response to different parameter values later.

Note: The stoichiometry must be defined and interpreted in consistency with the componenets’ units of measure. You can check the units by calling`<Component>.measured_as`

.

Note that `X_BH`

has a stoichiometric coefficient of 1, which means `X_BH`

is a “reference component” (i.e., the process rate calculated by the rate equation is also `X_BH`

’s rate of production/consumption).

```
[16]:
```

```
# This information can also be accessed by calling the `ref_component` property.
p1.ref_component
```

```
[17]:
```

```
# Another defining characteristics of a process is its rate equation, which is stored as a
# property of the `Process` object
p1.rate_equation
```

The process rate here is a “Monod-type” function of the concentrations (in \(g/m^3\)) of the four components (\(S_{NH}\), \(S_O\), \(S_S\), and \(X_{BH}\)) involved in this process.

How do I know this function takes concentrations instead of, let’s say, mass flowrates as input? Answer is – we can tell by looking at the units of the parameters, which are defined upon initiation of the `Process`

object.

```
[18]:
```

```
# If we access the private attribute
p1._rate_equation
```

This shows \(\mu_H\), \(K_{NH}\), \(K_{OH}\), and \(K_S\) are parameters in this rate equation. Although the information on their units are not stored with the `Process`

object, we can see their definitions in the documentation of the `ASM1`

class.

```
[19]:
```

```
# ?ASM1
```

Since \(\mu_H\) are defined in \(d^{-1}\), and \(K_{NH}\), \(K_{OH}\), and \(K_S\) in \(g/m^3\), only when \(S_{NH}\), \(S_O\), \(S_S\), and \(X_{BH}\) are in \(g/m^3\) can this rate equation output process rate in the correct unit (i.e., \(d^{-1}\)). This is defined by the mathematical formulation of ASM1. It is also common for biokinetic models to define state variables as component concentrations.

```
[20]:
```

```
# Now that we've understood the required input to the rate equation, we can try evaluating the
# process rate. Let's try with all component concentrations equal to 1.
import numpy as np
state = np.ones(len(cmps)+1)
state
```

This is a numpy array with 15 elements. The first 14 elements are component concentrations in the exact same order as `cmps`

, while the last element indicates volumetric flowrate \(Q\). Although \(Q\) is not used in process rate calculation, it is often an important state variable for the modeling of reactor-scale dynamics. Therefore, it is common for a lot of rate functions to ask for an array of length `len(cmps)+1`

as input, but it is ultimately determined by how the rate function
is coded.

```
[21]:
```

```
# Note that the `rate_equation` attribute only stores the formula.
# The evaluation of process rate is done through the `rate_function` attribute, which is rendered
# automatically from the rate equation.
p1.rate_function(state)
```

This result indicates that at the component concentrations and flowrate defined by the `state`

array, this process happens at a rate of 0.152 \(d^{-1}\) (approximately), which translates to a heterotrophic biomass growth rate of 0.152 \(gCOD/m^{3}/d\).

Note: Both`stoichiometry`

and`rate_equation`

are properties of the`Process`

object. Changing the formulas of stoichiometric coefficients or rate equation are prohibited after initiation of the object.

As you may remember from the `ASM1`

documentation and the `asm1.show()`

printout above, parameters and/or dynamic parameters are also important attributes for a `Process`

object because they represent intrinsic characteristics of the process and affect the evaluation of stoichiometry and process rate.

```
[22]:
```

```
# For `asm1` and the individual processes within `asm1`, parameters are stored as a dictionary
p1.parameters
# asm1.parameters
# p1.parameters is asm1.parameters
```

You may notice `p1`

and `asm1`

share the dictionary of parameters, although not all parameters included in this dictionary are used in the stoichiometry or process rate calculations for `p1`

. This is to make sure the parameter values are consistent across different processes in `asm1`

.

```
[23]:
```

```
# If you update parameter values for `p1`, the same parameters will be updated accordingly for
# `asm1` and any other processes in `asm1`.
p1.set_parameters(mu_H=6.0, K_S=8.0, Y_H=0.8)
asm1.parameters
```

```
[24]:
```

```
# If you evaluate process rate or stoichiometry again with the same input,
# you should now expect different output.
# p1.stoichiometry
p1.rate_function(state)
```

In addition to the attributes above, `conserved_for`

is another useful attribute, expecially for biochemical processes. It dictates the conservation principles followed by the stoichiometric coefficients of the the process.

```
[25]:
```

```
p1.conserved_for
```

We can see `p1`

is subject to conservations of COD, N, and electric charge, as defined by ASM1. This information is useful for automatic derivation of unknown stoichiometric coefficients upon initiation of the `Process`

object. If numerical values are provided for all stoichiometric coefficients when creating the `Process`

object, this information can then be used by `check_conservation`

to check if these values satisfy the conservation principles.

```
[26]:
```

```
# No return indicates that all materials in `conserved_for` are conserved.
# Otherwise, a warning or a `RuntimeError` will be raised.
p1.check_conservation()
```

You may have noticed `p1`

’s “dynamic parameters” is empty. This is because all of the parameters for the calculation of stoichiometric coefficients are static values, i.e., they don’t change during simulation. If a process’s stoichiometric coefficient is a function of time or other state variables but doesn’t have a nice and simple mathematical form, we can create a `DynamicParameter`

object and add it to the process. We will learn how to do that in section 3.

Back to top

### 1.3. `Processes`

and `CompiledProcesses`

¶

It is more often than not you’re working with more than one dynamic process, which often makes it impossible to derive an analytical solution to the dynamics or the steady state of the system. This is when `Processes`

and `CompiledProcesses`

come in handy.

As is mentioned above, `Processes`

is simply a collection of multiple `Process`

objects. Compiling a `Processes`

object yields a `CompiledProcesses`

object. Note that `CompiledProcesses`

is also a subclass of `Processes`

, i.e., it has all of `Processes`

’s attributes and more.

```
[27]:
```

```
qs.CompiledProcesses.__bases__
```

We know that `ASM1`

is a subclass of `CompiledProcesses`

, which makes `asm1`

an instance of `CompiledProcesses`

.

```
[28]:
```

```
# Let's verify that
isinstance(asm1, qs.CompiledProcesses)
```

So, is there a way to “decompile” a `CompiledProcesses`

object? Yes!

```
[29]:
```

```
asm1 = qs.Processes(asm1.tuple)
asm1.show()
```

```
[30]:
```

```
# Let's see what difference "decompiling" made
asm1.__dict__
```

This shows that a `Processes`

object only stores individual `Process`

objects as its attributes.

```
[31]:
```

```
asm1.compile()
asm1.__dict__
```

In comparison, a `CompiledProcesses`

object has many more attributes that are convenient for the application of dynamic process models, although no additional information was provided for `compile`

.

```
[32]:
```

```
# For example, the stoichiometric coefficients of all processes are compiled into a table that
# is in consistent format as a Petersen matrix
asm1.stoichiometry
```

```
[33]:
```

```
# Similarly for the rate equations
asm1.rate_equations
```

```
[34]:
```

```
# More importantly, the `rate_function` attribute of a `CompiledProcesses` now outputs an array
# with each element corresponding orderly to the individual processes
asm1.rate_function(state)
```

With this, the evaluation of all process rates can be done with one command and the output data type is convenient for matrix operations.

To calculate the overall rate of production for a component, we need to first estimate its rates of production (positive means being produced, negative means being consumed) by each process and then sum them up. In mathematical forms, let’s say \(a_{ij}\) represents the stoichiometric coefficient of component \(j\) in process \(i\) (i.e., value on the \(i\)th row and \(j\)th column of the stoichiometry matrix). Denote \(\rho_i\) as process \(i\)’s reaction rate, now the overall production rate of component \(j\) should be

In matrix notation, this calculation can be neatly described as

where \(\mathbf{A}\) is the stoichiometry matrix and \(\mathbf{\rho}\) is the array of process rates.

```
[35]:
```

```
# This matrix operation is already streamlined for `CompiledProcesses` objects, you can see the
# mathematical form of the rates of production as a function of component concentrations.
asm1.production_rates
```

```
[36]:
```

```
# To evaluate the rates of production for all components, all you need to
# do is to call the `production_rates_eval` method.
asm1.production_rates_eval(state)
```

With these convenient evaluation methods, `asm1`

can now be used in `SanUnit`

to model system dynamics (e.g., reactor effluent concentrations over time). We will introduce how to integrate process models with `SanUnit`

subclasses and how to enable dynamic simulations in a separate tutorial.

Back to top

## 2. Making a simple biokinetic model¶

Congratulations! You’ve understood the structure and the core capacities of the `_process`

module. 🎉

In this section, we will get some hands-on experience with the `_process`

module. For demonstration purpose, we will implement a simple biokinetic model with the classes we learned in the previous section. It’s even better if you have your own mathematical model of dynamic processes in mind.

Let’s get started with defining the mathematical model.

### The model¶

Let’s consider a simplified activated sludge system with 5 components besides water: - particulate organic substrate (\(X_S\), measured as COD) - soluble organic substrate (\(S_S\), measured as COD) - dissolved oxygen (\(O_2\), measured as itself) - dissolved CO2 (\(CO_2\), measured as C) - biomass (\(X_B\), measured as COD)

This model uses 3 dynamic processes to describe the transformations among the components: 1. Hydrolysis (i.e., \(X_S\) transforming into \(S_S\)) 2. Growth of biomass (i.e., \(X_B\) utilizing \(S_S\) and \(O_2\) to grow) 3. Decay of biomass (i.e., \(X_B\) dying and transforming into \(X_S\))

To mathematically characterize these processes, we need to define the stoichiometry and the process rate equation of each process. Let’s say, for process 1 and 3:

where \(k_{hyd}\) is the hydrolysis rate constant, \(b\) is the biomass decay rate constant, both in \(d^{-1}\cdot m^3\cdot g^{-1}\).

For biomass growth, assuming the yield of biomass is \(y_B\) [gCOD \(X_B\)/gCOD \(S_S\)], then the stoichiometry can be defined as:

where \(a_1\) and \(a_2\) need to be determined by conservations of COD and carbon. Similar to other activated sludge models, we can use a Monod-type function to represent biomass growth rate, i.e.,

### 2.1. Defining the components¶

For convenience, we will find similar components that have been defined in the default set. Note that we will need to consider conservations of COD and C in the stoichiometry of process 2. Therefore, we should make sure the `Component`

objects we pick have appropriate values for COD content (`i_COD`

) and carbon content (`i_C`

).

```
[37]:
```

```
# Load the default set of components
cmps_all = qs.Components.load_default()
# Copy select components from the default set and give them new IDs
X_S = cmps_all.X_B_Subst.copy('X_S')
S_S = cmps_all.S_F.copy('S_S')
O2 = cmps_all.S_O2.copy('O2')
CO2 = cmps_all.S_CO3.copy('CO2')
X_B = cmps_all.X_OHO.copy('X_B')
# Compile them into a new set. Don't forget the water!
cmps_bkm = qs.Components([X_S, S_S, O2, CO2, X_B, cmps_all.H2O])
cmps_bkm.compile()
cmps_bkm.show()
```

Note: You might’ve noticed H2O isn’t an active component in any process but is always included in the`CompiledComponents`

set. This is because water is assumed to be the bulk liquid and is needed for the convergence of total volumetric flows in simulation regardless of the process model.

```
[38]:
```

```
# Now we can check if their `measured_as` attributes are correctly set
# X_S.measured_as == 'COD'
# CO2.measured_as == 'C'
O2.measured_as == None
```

```
[39]:
```

```
# Then you can check relevant `i_` properties of the components
# For example, O2 should have a negative COD content, or more specifically -1 gCOD/gO2
O2.i_COD
```

```
[40]:
```

```
# A final preparation step is to set thermo with the components
qs.set_thermo(cmps_bkm)
```

It is common practice to wrap up the creation of model-specific `CompiledComponents`

objects in a function for convenient reuse. The `create_<ModelName>_cmps`

functions we’ve seen in the previous section is a good example.

Back to top

### 2.2. Defining the processes¶

With the components ready, we can now define individual dynamic processes in the model as `Process`

objects.

```
[41]:
```

```
# As always, let's take a look at the documentation, especially the init signature.
?qs.Process
```

For process 1 and process 3, the stoichiometric coefficients are known and the rate equations take simple forms with just one parameter.

```
[42]:
```

```
# To define `Process` objects for these two processes, we just need to input their
# stoichiometries and rate equations
# Process 1: hydrolysis
pc1 = qs.Process(ID='hydrolysis',
reaction='X_S -> S_S',
ref_component='X_S',
rate_equation='k_hyd*X_S',
conserved_for=(),
parameters=('k_hyd',))
```

Note that the information on stoichiometry is reflected in the `reaction`

argument. `ref_component`

is set to be `X_S`

because its stoichiometric coefficient is -1.

No conservation principle is used in the initiation of `pc1`

(i.e., `conserved_for`

is an empty tuple). This is fine, because there’s no unknown stoichiometric coefficient that relies on conservation principles to solve. However, this doesn’t mean the process isn’t subject to any conservation principle. The 1 to 1 stoichiometry of this process reflects the conservation of COD, as both \(X_S\) and \(S_S\) are measured as COD.

```
[43]:
```

```
# We can verify this with the `check_conservation` method
# To do that, `conserved_for` needs to be set first
pc1.conserved_for = ('COD',)
pc1.check_conservation()
```

No return means `pc1`

conserves COD.

Note:`conserved_for`

is treated as an iterable in the`check_conservation`

method. So, when there’s only one conservation principle, we should make sure the entire string is recognized as one element in iteration.

```
[44]:
```

```
# For example, if we do this, we will get an error message because the function interprets it as
# checking conservations of 'C', 'O', and 'D', respectively
# pc1.conserved_for = ('COD')
# pc1.check_conservation()
```

Note that the initiation also asks for `parameters`

. This is to help identify any symbolic parameter in the stoichiometry and the rate equation.

```
[45]:
```

```
# Upon initiation, the parameters are stored as symbols. We still need to set values to them
# before we can evalute process rate.
pc1.show()
```

```
[46]:
```

```
# At this point, the initiation of process 3 should be quite straightforward.
# Here shows an alternative way to input stoichiometry
pc3 = qs.Process('decay',
reaction={'X_B':-1, 'X_S':1},
ref_component='X_B',
rate_equation='b*X_B',
conserved_for=('COD',),
parameters=('b',))
# But how do we initiate process 2, which has known stoichiometric coefficients?
pc2 = qs.Process('growth',
reaction='[1/y_B]S_S + [?]O2 -> X_B + [?]CO2',
ref_component='X_B',
rate_equation='mu_B * S_S/(K_S+S_S) * X_B',
conserved_for=('COD', 'C'),
parameters=('y_B', 'mu_B', 'K_S'))
pc2.show()
```

Because there’re two unknown stoichiometric coefficients, we need two conservation principles (i.e., two equations) to determine their values. It is convenient to solve for stoichiometric coefficients marked as `?`

by specifying `conserved_for`

in the initiation of `Process`

objects.

Back to top

### 2.3. Compiling the processes¶

```
[47]:
```

```
# Now the final step is to compile the individual processes into a biokinetic model
bkm = qs.Processes([pc1, pc2, pc3])
bkm.compile()
bkm.show()
```

```
[48]:
```

```
# Parameters in the stoichiometry and rate equations across all processes are compiled into
# a shared dictionary.
# This also means, we should absolutely avoid using the same symbol for different parameters in
# multiple processes.
bkm.parameters
```

```
[49]:
```

```
# After setting parameter values, the model will be ready
bkm.set_parameters(k_hyd=3.0, y_B=0.8, mu_B=4.0, K_S=9.0, b=0.4)
bkm.stoichiometry
```

```
[50]:
```

```
bkm.production_rates
```

Back to top

### 2.4. Batch mode for creating processes¶

It’d be very tedious to create individual `Process`

objects one by one and compile them, especially for models with a large number of parallel processes and/or components. Therefore, the `_process`

module includes convenient classmethods to create processes in batch.

```
[51]:
```

```
# Stoichiometry and rate equations are usually described in a table format
from qsdsan.utils import load_data
df_bkm = load_data('assets/_bkm.tsv', index_col=0)
df_bkm
```

This table contains the stoichiometry matrix and the array of rate equations. If the table is organized with process IDs as row indices and component IDs as column names, we can use it as direct input to the `Processes.load_from_file`

classmethod to create processes in batch.

```
[52]:
```

```
# The same amount of information is still required to create the model.
bkm_batch = qs.Processes.load_from_file(
path='assets/_bkm.tsv',
conserved_for=('COD', 'C'),
parameters=('k_hyd', 'y_B', 'mu_B', 'K_S', 'b'),
compile=True
)
bkm_batch.show()
```

```
[53]:
```

```
# You can see `bkm_batch` is equivalent to the `bkm` we created above
bkm_batch.stoichiometry
# bkm_batch.rate_equations
# bkm_batch.parameters
```

```
[54]:
```

```
# The reference component of each process is inferred from its stoichiometry
bkm_batch.decay.ref_component
```

```
[55]:
```

```
# `conserved_for` now applies to all processes,
# the following will trigger an error
# bkm_batch.decay.check_conservation()
```

This means the decay process defined above conserves COD but not carbon, suggesting that we could revisit the `i_C`

values for \(X_B\) and \(X_S\) or consider whether there are other carbon-containing products in this process.

```
[56]:
```

```
# `i_C` values for each component
dict(zip(cmps_bkm.IDs, cmps_bkm.i_C))
```

Back to top

## 3. Advanced features¶

`qsdsan`

aims to provide a flexible and user-friendly platform for process modeling. Therefore, many advanced features are included in the `_process`

module to streamline the implementation of process models with various mathematical complexity.

Let’s again use the model above as an example to learn about these features.

### 3.1. `DynamicParameter`

¶

What if the yield of biomass growth is dependent on the abundance of soluble organic substrate? Assume this relation can be described as

In this case, we will need to define \(y_B\) as a `DynamicParameter`

object because its value changes with system state during simulation.

```
[57]:
```

```
@pc2.dynamic_parameter(symbol='y_B', params={'y_max':0.8})
def yB_eval(state_arr, params):
X_S, S_S = state_arr[:2] # assuming the first 6 elements are component concentrations
y_max = params['y_max']
y_B = y_max * (S_S/(X_S+S_S))**(1/2)
return y_B
```

The function `yB_eval`

must take no more than two positional arguments, with the first one expecting an array of state variable values and the second one expecting a dictionary of additional inputs required for the evaluation of this dynamic parameter.

Previous tutorials have covered the use of `decorator`

. You can revisit this tutorial for more details.

```
[58]:
```

```
# At this point, the value for `y_B` is not updated yet, since we haven't evalutated it
# with input of component concentrations
bkm.parameters
```

```
[59]:
```

```
# But the list of dynamic parameters have been updated
pc2.show()
```

```
[60]:
```

```
# `y_B` is now a `DynamicParameter` object stored in the `_dyn_params` attribute of the process
type(pc2._dyn_params['y_B'])
```

```
[61]:
```

```
# Assuming component concentrations are all 1.
state_bkm = np.ones(6)
# This function evaluates the dynamic parameters
bkm.params_eval(state_bkm)
bkm.parameters
```

```
[62]:
```

```
# Then the stoichiometry is updated accordingly
bkm.stoichiometry
```

Note that in simulations, you don’t need to call the `params_eval`

method explicitly unless you would like to capture the dynamics of parameter values. The evaluation of dynamic parameters is performed within `production_rates_eval`

(i.e., the calculation of overall production rates of the components).

Back to top

### 3.2. `Kinetics`

and `MultiKinetics`

¶

`Kinetics`

and `MultiKinetics`

are classes that allow defining process rates as functions.

`Kinetics`

is a subclass of `DynamicParameter`

, because process rate in essence is also a function of state variables and additional static parameters.

```
[63]:
```

```
# The use of `Kinetics` to define process rate is also similar
# For example, for the decay process
@pc3.kinetics(parameters={'b':0.4})
def rho3_eval(state_arr, params):
X_B = state_arr[4]
b = params['b']
return b*X_B
pc3.rate_function
```

This may seem trivial, but it’s helpful when the form of process rate gets more complex. For example, the rate of substrate uptake in anaerobic digestion model no.1 is dependent on pH, which is dynamically determined by a linear equation system and thus cannot be expressed neatly as a function of state variables.

`MultiKinetics`

serves a similar purpose but calculates the rate of multiple processes at onces. It is designed to be used with `CompiledProcesses`

objects. This is convenient for models with similar process rate algorithms across multiple processes in the sense that vectorizing the shared algorithm can often improve computational speed in simulation. Another application of `MultiKinetics`

is models sharing intermediate variables (e.g., pH) across multiple process rate calculation.

Let’s assume the reaction rates of all three processes in the model above are now dependent on the abundance of oxygen. Here are the new rate equations:

How do we use `MultiKinetics`

to implement this change?

```
[64]:
```

```
# Define the function
def rhos_eval(state_arr, params):
X_S, S_S, O2, CO2, X_B = state_arr[:5]
I_O2 = (O2/(O2+params['K_O2']))**(1/2)
rho1 = params['k_hyd']*X_S*I_O2
rho2 = params['mu_B']*X_B*I_O2
rho3 = params['b']*X_B*(1.5-I_O2)
return np.array([rho1, rho2, rho3])
# Set it as the rate function of the `CompiledProcesses` object
bkm.set_rate_function(rhos_eval)
# Assign parameter values
bkm.rate_function.set_params(K_O2=2.0, **bkm.parameters)
# The evaluation of the rate function is the same as before
bkm.rate_function(state_bkm)
```

Back to top