Low-Level Disinfection with Ozone

Note

This example is exactly the same as the corresponding the high level example. Thus the two interfaces can be compared.

While it does not make a huge difference for this simple example, the advantages grow already for the next one and especially for more complex systems.

As a first small example, we consider the following code, which implements a simple connection of two continuously stirred tank reactors (CSTR) as commonly used in water treatment. It serves as a quick example to start with. For more detailed explanations, we referred to the API of the ProcessSimulator.

This example implements two reactors used for ozonation of water. In the first reactor, the concentration of ozone is kept constant, to simulate a supply of ozone in the beginning. The stoichiometric matrix for the ozonation process is defined in the library and accessible under the name ozonation. This system is depicted as follows:

For building up this system, it is proceeded as follows:

Setup

First, the needed packages are imported and constants defined:

using BioChemicalTreatment # Reactors etc.
using ModelingToolkit # Modeling Framework
using DifferentialEquations # Solve the equations
using Plots # Plotting

# Define the needed constants
v = 1000.0/2    # Reactor volume
q = 10000.0/24  # Flow rate through the reactor
kO3 = 10.0/24  # Reaction rate for the decay of ozone
kd = 1500.0/24  # Reaction rate for the decay of the bacteria

System Creation

Then, all needed systems (here inflow and the reactors) are created.

We start with the reaction rates system for the disinfection dynamics. This system defines all compounds needed for this specific reaction, which are needed afterwards for defining the CSTRs. Thus it is convenient to start with this system. We need separate rates systems for each reactor, thus we create two identical ones, one for each CSTR. For the parameters, we set the decay of ozone in the first reactor to zero, to model that it is supplied in this tank.

@named disinfection_r1 = OzoneDisinfection(;kO3 = 0, kd)
@named disinfection_r2 = OzoneDisinfection(;kO3, kd)
nameof(disinfection_r2)

We continue with the two CSTRs. Each one is a system with an inflow and outflow port, and further an input with the provided rates (from the reaction rate system just created) and an output for the states, which are needed by the disinfection reaction rates system for computation. To create it we provide it with first the volume of the tank, and second with the compounds to have in the flow. The compounds in the flow can be extracted as states from the OzoneDisinfection systems above. Here we use the states of the corresponding rates, which is convenient. Note however, that this is not necessary and the states could as well be constructed independently of the rates.

@named CSTR1 = CSTR(v, unknowns(states(disinfection_r1)); initial_states = [0.5, 0])
@named CSTR2 = CSTR(v, unknowns(states(disinfection_r2)); initial_states = [0.5, 0])
nameof(CSTR2)

Finally, the influent is to be created: It needs the variables to provide an input for, and the corresponding values. We take the values directly from the variables in the inflow of the first CSTR, which are:

unknowns(inflows(CSTR1)[1])

3-element Vector{SymbolicUtils.BasicSymbolic{Real}}:
 q(t)
 S_O3(t)
 X_B(t)

To this, the value vector has to be aligned. The ozone concentration is set to 0.5 to keep this level in the first reactor (Consider the ozone to be supplied before the first reactor) and the bacteria concentration is set to 1, such that one can interpret the bacteria in the effluent as the percentage of bacteria that are remaining.

@named influent = Influent([q, 0.5, 1], unknowns(inflows(CSTR1)[1]))
nameof(influent)

Connection

Then, the connection equations for the systems are created. We need on one hand to connect the reaction rate systems to their reactors, and on the other the flows. For the reaction rates, the states of the reactor must be connected to the corresponding reaction rates systems, and the rates supplied by the reaction rate system to the reactor. For the flows, simply inflows and outflows can be connected as desired:

connection_eqns = [
	# Connect CSTR1 with its rates
	connect(states(CSTR1), states(disinfection_r1)),
	connect(rates(disinfection_r1), rates(CSTR1)),
	# Connect CSTR2 with its rates
	connect(states(CSTR2), states(disinfection_r2)),
	connect(rates(disinfection_r2), rates(CSTR2)),
	# Connect the flows
	connect(outflows(influent)[1], inflows(CSTR1)[1]),
	connect(outflows(CSTR1)[1], inflows(CSTR2)[1]),
]

With the connections, we can then build the overall model:

@named system = ODESystem(connection_eqns, t, systems = [
	influent,
	CSTR1, disinfection_r1,
	CSTR2, disinfection_r2
])

And simplify it to get the minimum set of equations:

system_simplified = structural_simplify(system)
full_equations(system_simplified) # Display the equations

\[ \begin{align} \frac{\mathrm{d} \mathtt{CSTR1.states.S\_O3}\left( t \right)}{\mathrm{d}t} &= \frac{ - \mathtt{influent.q\_const.k} \mathtt{CSTR1.states.S\_O3}\left( t \right)}{\mathtt{CSTR1.V}} + \frac{\mathtt{influent.S\_O3\_const.k} \mathtt{influent.q\_const.k}}{\mathtt{CSTR1.V}} - \mathtt{disinfection\_r1.kO3} \mathtt{CSTR1.states.S\_O3}\left( t \right) \\ \frac{\mathrm{d} \mathtt{CSTR1.states.X\_B}\left( t \right)}{\mathrm{d}t} &= \frac{\mathtt{influent.X\_B\_const.k} \mathtt{influent.q\_const.k}}{\mathtt{CSTR1.V}} + \frac{ - \mathtt{influent.q\_const.k} \mathtt{CSTR1.states.X\_B}\left( t \right)}{\mathtt{CSTR1.V}} - \mathtt{disinfection\_r1.kd} \mathtt{CSTR1.states.X\_B}\left( t \right) \mathtt{CSTR1.states.S\_O3}\left( t \right) \\ \frac{\mathrm{d} \mathtt{CSTR2.states.S\_O3}\left( t \right)}{\mathrm{d}t} &= \frac{\mathtt{influent.q\_const.k} \mathtt{CSTR1.states.S\_O3}\left( t \right)}{\mathtt{CSTR2.V}} + \frac{ - \mathtt{influent.q\_const.k} \mathtt{CSTR2.states.S\_O3}\left( t \right)}{\mathtt{CSTR2.V}} - \mathtt{disinfection\_r2.kO3} \mathtt{CSTR2.states.S\_O3}\left( t \right) \\ \frac{\mathrm{d} \mathtt{CSTR2.states.X\_B}\left( t \right)}{\mathrm{d}t} &= \frac{ - \mathtt{influent.q\_const.k} \mathtt{CSTR2.states.X\_B}\left( t \right)}{\mathtt{CSTR2.V}} + \frac{\mathtt{influent.q\_const.k} \mathtt{CSTR1.states.X\_B}\left( t \right)}{\mathtt{CSTR2.V}} - \mathtt{disinfection\_r2.kd} \mathtt{CSTR2.states.X\_B}\left( t \right) \mathtt{CSTR2.states.S\_O3}\left( t \right) \end{align} \]

Simulation and Plotting

Finally, we simulate and plot the result.

## Simulate the system (for 1 day)
# First create a ODEProblem for the time span (0, 1) and then solve it
prob = ODEProblem(system_simplified, [], (0, 1), [])
sol = solve(prob)

## Plot the output
plot(sol,
	idxs = [CSTR1.states.X_B ,CSTR2.states.X_B], # Plot bacteria in tank
	legend = :right,
	title = "Remaining fraction of bacteria in each tank")