A Potion Master's Guide
Welcome, apprentice! This is not just any simulator; it is a magical vessel where the arcane arts of biotechnology and fermentation come to life. Use this guide to master the fundamental spells and advanced enchantments needed to brew the perfect potion.
Step 1: The Recipe's Foundation
Every great potion begins with a foundational plan. Your first two choices—the Fermentation Mode and the Product Induction Strategy—define the entire course of the brew.
Fermentation Mode
- Batch: Simple and direct. Add all ingredients at the start and let the magic unfold.
- Fed-Batch: Strategic and powerful. Continuously add ingredients over time to reach higher yields.
Induction Strategy
This is the core "spell" that tells your microbes when to create the product, based on triggers like adding a catalyst (Inducer), removing a blocker (Repressor), or even near-starvation (Carbon Limitation).
Step 2: Fine-Tuning the Spell
True mastery comes from understanding the subtle forces at play. These parameters allow you to manipulate the environment and the very nature of your microbes.
- The Cell's Spellbook: Define your microbe's inherent abilities.
µ_max sets its maximum growth speed, while Yield Coefficients (Y_...) determine its efficiency.
- Environmental Stress: Microbes are sensitive to Temperature and pH. Push them too far from their optimal conditions, and their magic will wane.
- Product Inhibition (
P_crit): Beware a potion too potent! As your product accumulates, it can become toxic to its creators, halting all growth even if food remains.
- Catalyst Effects (
m_I): The inducer that starts your reaction might not last forever. If it gets consumed, your production will mysteriously cease.
- Life and Death (
k_d): Microbes aren't immortal. The death rate (k_d) causes the population to decline over time, especially during starvation. When cells lyse, their nutrients can be recycled (Y_S1/lysis), providing a last bit of sustenance for the survivors.
Step 3: The Automated Cauldron
The best Potion Masters don't stir by hand. The cauldron is enchanted with powerful automation charms that mimic real-world control systems.
- DO-Stat Feeding: Instead of guessing, use this "smart" feed strategy. The cauldron senses when the microbes are hungry and adds food automatically, preventing waste and stress.
- PID Controllers: Sophisticated algorithms constantly adjust Agitation Speed and pH Correction Feeds to keep the environment perfectly stable, allowing your microbes to focus on their magical work.
Step 4: Taming the Unpredictable
Even the most skilled Potion Master knows that magic is not always predictable. The real world is full of randomness. This simulator allows you to embrace that chaos to create more realistic brews.
- Biological Variability: No two microbes are identical twins. This setting introduces random fluctuations to the growth and production rates, simulating the inherent stochasticity of life. This is why no two experiments are ever exactly the same!
- Process & Sensor Noise: Your cauldron and its tools are not perfect. Pumps may fluctuate slightly, and magical sensors can be affected by bubbles or electrical spirits. This noise adds a layer of instrumental realism, forcing your control algorithms to react to an imperfect world.
- Parameter Uncertainty: What if your ancient grimoire is slightly inaccurate? This allows you to simulate having imperfect knowledge of the true biological parameters, testing how robust your recipe is against the unknown.
If the cauldron spills over on you, contact Anupam Chowdhury.
The Grimoire of Equations
The cauldron's behavior is governed by a system of ordinary differential equations (ODEs). These formulae dictate the kinetics of growth, the magic of production, and the balance of ingredients within the brew. The model is structured in three layers: the per-cell kinetic rates, the overall mass balances in the cauldron, and the environmental control systems.
Part 1: The Cell's Engine - Microbial Kinetics
This layer describes the fundamental rates at which individual cells grow and create products based on their immediate surroundings. These are the "per-cell" instructions.
A. Specific Growth Rate ($\mu$)
The specific growth rate, $\mu$ (1/hr), is the cornerstone of the simulation. It's determined by the availability of substrates, the concentration of dissolved oxygen, and other environmental conditions like pH and temperature.
Oxygen's Influence: Metabolic Shifting
Cells dynamically shift their metabolism based on dissolved oxygen ($C_{O2}$). We model this with two competing factors:
Aerobic Factor ($f_{O2,aerobic}$): $$f_{O2,aerobic} = \frac{C_{O2}}{K_{O2} + C_{O2}}$$
Anaerobic Factor ($f_{O2,anaerobic}$): $$f_{O2,anaerobic} = \frac{K_{O2,ana\_inhib}}{K_{O2,ana\_inhib} + C_{O2}}$$
Parameter Definitions:
- $K_{O2}$: Oxygen affinity constant for aerobic metabolism.
- $K_{O2,ana\_inhib}$: Oxygen concentration that causes 50% inhibition of anaerobic pathways. A lower value means anaerobic metabolism is more sensitive to oxygen.
Calculating the Total Growth Rate
These oxygen factors scale the maximum potential growth rates on each substrate. These are then combined and adjusted by environmental factors to get the final growth rate.
Growth on Primary Carbon (S1): $$\mu_{S1} = \left( \mu_{max,S1,aerobic} \cdot f_{O2,aerobic} + \mu_{max,S1,anaerobic} \cdot f_{O2,anaerobic} \right) \cdot \frac{S_1}{K_{S1} + S_1}$$
Growth on Secondary Carbon (S2): $$\mu_{S2} = \left( \mu_{max,S2,aerobic} \cdot f_{O2,aerobic} \right) \cdot \frac{S_2}{K_{S2} + S_2}$$
Environmental Factors ($f(pH), f(T)$): These are Cardinal models that return a value between 0 and 1.
$$f(pH) = \frac{(pH - pH_{min})(pH - pH_{max})}{(pH - pH_{min})(pH - pH_{max}) - (pH - pH_{opt})^2}$$
$$f(T) = \frac{(T - T_{max})(T - T_{min})^2}{(T_{opt} - T_{min})[(T_{opt} - T_{min})(T - T_{opt}) - (T_{opt} - T_{max})(T_{opt} + T_{min} - 2T)]}$$
Product Inhibition ($f(P)$): High product concentrations can be toxic and inhibit growth. This is modeled with an additional term:
$$f(P) = \left(1 - \frac{P}{P_{crit}}\right)^{n_P}$$
Total Specific Growth Rate:
$$\mu = (\mu_{S1} + \mu_{S2}) \cdot f(I) \cdot f(pH) \cdot f(T) \cdot f(P)$$
Parameter Definitions:
- $\mu_{max,S1,aerobic}, \mu_{max,S2,aerobic}$: Maximum specific growth rate on S1 and S2 under optimal aerobic conditions.
- $\mu_{max,S1,anaerobic}$: Maximum specific growth rate on S1 under optimal anaerobic conditions.
- $K_{S1}, K_{S2}$: Substrate affinity constants for growth.
- $f(I)$: Inducer inhibition factor for growth, typically $\frac{K_{Ii}}{K_{Ii} + I}$.
- $P_{crit}$: The critical product concentration at which growth stops.
- $n_P$: The exponent for product inhibition, controlling the sharpness of the effect.
B. Specific Production Rate ($q_p$)
Similar to growth, the specific production rate ($q_p$) depends on substrates, oxygen, and its own unique set of environmental conditions for temperature and pH. These factors are combined to determine the final rate.
Blended Maximum Production Rate ($q_{p,max,blended}$):
$$q_{p,max,blended} = q_{p,max,aerobic} \cdot f_{O2,aerobic} + q_{p,max,anaerobic} \cdot f_{O2,anaerobic}$$
Environmental Factors for Production ($f_p(pH), f_p(T)$):
These use the same Cardinal model structure as growth, but with parameters specific to production (denoted with a subscript 'p').
$$f_p(pH) = \frac{(pH - pH_{min,p})(pH - pH_{max,p})}{(pH - pH_{min,p})(pH - pH_{max,p}) - (pH - pH_{opt,p})^2}$$
$$f_p(T) = \frac{(T - T_{max,p})(T - T_{min,p})^2}{(T_{opt,p} - T_{min,p})[(T_{opt,p} - T_{min,p})(T - T_{opt,p}) - (T_{opt,p} - T_{max,p})(T_{opt,p} + T_{min,p} - 2T)]}$$
This blended rate is then multiplied by terms specific to the induction strategy and the environmental factors:
Constitutive: $$q_p = q_{p,max,blended} \left( \frac{S_1}{K_{p,S1} + S_1} \right) \cdot f_p(pH) \cdot f_p(T)$$
Induction Only: $$q_p = q_{p,max,blended} \left( \frac{S_1}{K_{p,S1} + S_1} \right) \left( \frac{I}{K_{p,I} + I} \right) \cdot f_p(pH) \cdot f_p(T)$$
Second Substrate (S2) Only: $$q_p = q_{p,max,blended} \left( \frac{S_2}{K_{p,S2} + S_2} \right) \cdot f_p(pH) \cdot f_p(T)$$
Co-feed: $$q_p = q_{p,max,blended} \left( \frac{S_1}{K_{p,S1} + S_1} \right) \left( \frac{S_2}{K_{p,S2} + S_2} \right) \cdot f_p(pH) \cdot f_p(T)$$
Induction + Co-feed: $$q_p = q_{p,max,blended} \left( \frac{S_1}{K_{p,S1} + S_1} \right) \left( \frac{S_2}{K_{p,S2} + S_2} \right) \left( \frac{I}{K_{p,I} + I} \right) \cdot f_p(pH) \cdot f_p(T)$$
Repressor-based: $$q_p = q_{p,max,blended} \left( \frac{S_1}{K_{p,S1} + S_1} \right) \left( \frac{K_{pR}^{n_R}}{K_{pR}^{n_R} + R^{n_R}} \right) \cdot f_p(pH) \cdot f_p(T)$$
Carbon Limitation: $$q_p = q_{p,max,blended} \left( \frac{S_1}{K_{p,S1} + S_1} \right) \left( \frac{K_{pS1,lim}^{n_{S1,lim}}}{K_{pS1,lim}^{n_{S1,lim}} + S_1^{n_{S1,lim}}} \right) \cdot f_p(pH) \cdot f_p(T)$$
Parameter Definitions:
- $q_{p,max,aerobic}, q_{p,max,anaerobic}$: The maximum specific production rate under aerobic/anaerobic conditions.
- $K_{p,S1}, K_{p,S2}, K_{p,I}$: Affinity constants for production related to S1, S2, and Inducer (I).
- $K_{pR}, n_R$: The Repressor constant and Hill Coefficient for repression.
- $K_{pS1,lim}, n_{S1,lim}$: The S1 limitation constant and Hill Coefficient for induction by limitation.
- $pH_{opt,p}, pH_{min,p}, pH_{max,p}$: The optimal, minimum, and maximum pH for product formation.
- $T_{opt,p}, T_{min,p}, T_{max,p}$: The optimal, minimum, and maximum temperature for product formation.
Part 2: The Cauldron's Contents - Mass Balances
This layer scales up the per-cell rates to describe the change in total amounts of biomass, substrates, and product within the entire volume of the cauldron over time.
Blended Yield Coefficients
The efficiency of converting substrate to useful products changes with oxygen availability. We calculate effective, blended yields that are weighted averages of the aerobic and anaerobic values.
Effective Biomass Yield ($Y_{X/S1,eff}$): $$Y_{X/S1,eff} = \frac{Y_{X/S1,aerobic} f_{O2,aerobic} + Y_{X/S1,anaerobic} f_{O2,anaerobic}}{f_{O2,aerobic} + f_{O2,anaerobic} + \epsilon}$$
Effective Product Yield ($Y_{P/S1,eff}$): $$Y_{P/S1,eff} = \frac{Y_{P/S1,aerobic} f_{O2,aerobic} + Y_{P/S1,anaerobic} f_{O2,anaerobic}}{f_{O2,aerobic} + f_{O2,anaerobic} + \epsilon}$$
Effective Maintenance ($m_{S1,eff}$): $$m_{S1,eff} = \frac{m_{S1,aerobic} f_{O2,aerobic} + m_{S1,anaerobic} f_{O2,anaerobic}}{f_{O2,aerobic} + f_{O2,anaerobic} + \epsilon}$$
Core Balance Equations
The change in each component is its inflow minus outflow plus what's generated minus what's consumed.
Biomass (X): $$\frac{dX}{dt} = (\mu - k_d) \cdot X - \frac{X}{V} (F_{total})$$
Primary Substrate ($S_1$): $$\frac{dS_1}{dt} = \frac{F_1 S_{1,in}}{V} - \left( \frac{\mu_{S1}}{Y_{X/S1,eff}} + m_{S1,eff} \right) X - \frac{q_p}{Y_{P/S1,eff}} X + (k_d \cdot X \cdot Y_{S1/lysis}) - \frac{S_1}{V} (F_{total})$$
Secondary Substrate ($S_2$): $$\frac{dS_2}{dt} = \frac{F_2 S_{2,in}}{V} - \frac{\mu_{S2}}{Y_{X/S2,aerobic}} X - \frac{q_p}{Y_{P/S2,aerobic}} X - \frac{S_2}{V} (F_{total})$$
Product (P): $$\frac{dP}{dt} = q_p \cdot X - \frac{P}{V} (F_{total})$$
Parameter Definitions:
- $F_{total}$: Sum of all incoming feed rates ($F_1+F_2+F_I+F_{base}+F_{acid}$).
- $Y_{...}, m_{...}$: Various yield and maintenance coefficients for aerobic and anaerobic conditions.
- $k_d$: The specific death rate (1/hr), a first-order constant representing the rate of cell decay.
- $Y_{S1/lysis}$: The yield of substrate from lysis (g/g). Represents the fraction of dead biomass that is recycled back into a usable primary substrate (S1).
Part 2A: Genetic Drift & Population Dynamics
To model genetic drift, the total biomass (X) is split into two competing populations: Producers ($X_p$) and Non-producers ($X_n$). Non-producers have a growth advantage ($\alpha_n$) but do not create product. They are generated from producers at a specific mutation rate ($k_m$).
Total Biomass: $$X_{total} = X_p + X_n$$
Growth Rate of Non-Producers: $$\mu_n = \mu \cdot \alpha_n$$
Producer Biomass Balance ($X_p$):
$$\frac{dX_p}{dt} = (\mu \cdot X_p) - (k_d \cdot X_p) - (k_m \cdot X_p) - \left(\frac{X_p}{V} \cdot F_{total}\right)$$
Non-Producer Biomass Balance ($X_n$):
$$\frac{dX_n}{dt} = (\mu_n \cdot X_n) + (k_m \cdot X_p) - (k_d \cdot X_n) - \left(\frac{X_n}{V} \cdot F_{total}\right)$$
Part 3: Environmental Control Systems
This final layer models the physical systems and control algorithms that maintain the fermentation environment.
A. Dissolved Oxygen & Agitation Control
The DO concentration is a dynamic balance between the Oxygen Transfer Rate (OTR) from air and the Oxygen Uptake Rate (OUR) by the cells.
DO Mass Balance:
$$\frac{dC_{O2}}{dt} = \underbrace{k_L a (C^{*}_{O2} - C_{O2})}_{\text{OTR}} - \underbrace{OUR_{current}}_{\text{OUR}}$$
Mass Transfer Coefficient ($k_L a$): $$k_L a = k_{La,c} \cdot N^{k_{La,a}} \cdot F_{air}^{k_{La,b}}$$
Target OUR ($OUR_{target}$): Calculated from the actual aerobic components of metabolism.
$$OUR_{target} = \underbrace{\left[ \frac{\mu_{aerobic\_actual}}{Y_{X/O2}} + m_{O2} \right] X_{total}}_{\text{Growth & Maintenance}} + \underbrace{\left[ \frac{q_{p,aerobic\_actual}}{Y_{P/O2}} \right] X_p}_{\text{Production}} $$
Dynamic OUR ($OUR_{current}$): The actual OUR smoothly approaches the target, governed by a metabolic time constant, $\tau_{OUR}$.
$$\frac{d(OUR_{current})}{dt} = \frac{OUR_{target} - OUR_{current}}{\tau_{OUR}}$$
PID Controller for Agitation
To maintain the DO at its setpoint, a Proportional-Integral-Derivative (PID) controller adjusts the agitation speed ($N$). A deadband is used to prevent the controller from reacting to minor noise.
Error: $$e(t) = DO_{setpoint} - DO_{measured}(t)$$
Controller Output (if $|e(t)| > \text{Deadband}$): $$ N(t) = K_{c,DO} e(t) + \frac{K_{c,DO}}{\tau_{i,DO}} \int_{0}^{t} e(t') dt' + K_{c,DO} \tau_{d,DO} \frac{de(t)}{dt} $$
Parameter Definitions:
- $\tau_{OUR}$: The metabolic response time constant for OUR.
- $C^{*}_{O2}$: The saturation concentration of oxygen in the liquid.
- $Y_{X/O2}, Y_{P/O2}, m_{O2}$: Yield and maintenance coefficients for oxygen consumption.
- $K_{c,DO}$: The controller gain (proportional term).
- $\tau_{i,DO}$: The integral time constant.
- $\tau_{d,DO}$: The derivative time constant.
B. pH Dynamics & Control
As cells grow, they produce acidic byproducts. A control system adds base or acid to maintain the pH at a desired setpoint.
Effective Acid Yield ($Y_{H/S1,eff}$): This is a blended yield based on oxygen availability.
$$Y_{H/S1,eff} = \frac{Y_{H/S1,aerobic} f_{O2,aerobic} + Y_{H/S1,anaerobic} f_{O2,anaerobic}}{f_{O2,aerobic} + f_{O2,anaerobic} + \epsilon}$$
Net Acid Accumulation: $$\frac{d(\text{Total Acid})}{dt} = \left( \frac{\mu_{S1}}{Y_{X/S1,eff}} \cdot Y_{H/S1,eff} \cdot X \cdot V \right) - (F_{base} C_{base,in}) + (F_{acid} C_{acid,in})$$
pH Calculation: $$pH_t = pH_{initial} - \frac{(\text{Total Acid}_t / V_t)}{\text{Buffer Capacity}}$$
PID Controller for pH
A PID controller calculates a required correction. If the output is positive, base is added; if negative, acid is added. Action is only taken if the error exceeds the deadband.
Error: $$e(t) = pH_{setpoint} - pH_{measured}(t)$$
Controller Output ($P_{out}$ if $|e(t)| > \text{Deadband}$): $$ P_{out}(t) = K_{c,pH} e(t) + \frac{K_{c,pH}}{\tau_{i,pH}} \int_{0}^{t} e(t') dt' + K_{c,pH} \tau_{d,pH} \frac{de(t)}{dt} $$
IF $P_{out} > 0$, THEN $F_{base} \propto P_{out}$
IF $P_{out} < 0$, THEN $F_{acid} \propto -P_{out}$
C. Advanced Feed Strategies: DO-Stat Control
Instead of a pre-programmed feed, the DO-Stat strategy feeds the culture intelligently based on its metabolic state.
Trigger Condition: A feed "shot" is delivered when the rate of change of DO exceeds a threshold, indicating substrate exhaustion.
IF $\frac{d(DO)}{dt} > \text{Threshold}$ AND $S_1$ is low THEN trigger a feed shot.
Feed Rate Equation: $$ F_{1,rate} = F_{basal} + F_{shot} $$
Parameter Definitions:
- DO Spike Threshold: The minimum rate of DO change (%/hr) required to trigger a feed shot.
- Feed Shot Size: The volume of feed (mL/L) delivered per shot.
- Basal Feed Rate: A low, constant feed to prevent complete starvation.
- Lockout Period: A brief duration after a shot where no new shots can be triggered.
Part 4: Introducing Stochasticity (Noise)
The previous equations describe a perfect, deterministic world. To model a more realistic process, we introduce random fluctuations, transforming our model from a system of Ordinary Differential Equations (ODEs) to Stochastic Differential Equations (SDEs). An SDE has two parts: a deterministic "drift" (the original ODE) and a random "diffusion" term.
General SDE Form: $$dX_t = \underbrace{f(X_t, t)dt}_{\text{Drift (deterministic)}} + \underbrace{g(X_t, t)dW_t}_{\text{Diffusion (stochastic)}}$$
Here, $dW_t$ represents a Wiener process or "Brownian motion," which is the mathematical formalization of random noise that accumulates over time.
A. Biological Variability
We model biological variability by adding a diffusion term to the mass balance equations for biomass and product. The magnitude of the noise is proportional to the current amount of biomass or product.
Biomass (X): $$dX = \left[ (\mu - k_d)X - \frac{X}{V}F_{total} \right]dt + \sigma_b X dW_t$$
Product (P): $$dP = \left[ q_p X - \frac{P}{V}F_{total} \right]dt + \sigma_p P dW_t$$
Parameter Definitions:
- $\sigma_b$: The intensity of the noise for biomass, corresponding to Growth Variability.
- $\sigma_p$: The intensity of the noise for the product, corresponding to Production Variability.
B. Process and Sensor Noise
Other sources of noise are modeled more directly. Sensor noise is added to the "measured" value that the PID controllers see, while process noise is added to the final output of an actuator (like a feed pump).
Example (Sensor Noise): $$DO_{measured} = DO_{actual} + \mathcal{N}(0, \sigma^2_{sensor})$$
Example (Process Noise): $$F_{1,actual} = F_{1,setpoint} \times (1 + \mathcal{N}(0, \sigma^2_{process}))$$
Where $\mathcal{N}(0, \sigma^2)$ is a random number drawn from a Gaussian (normal) distribution with a mean of 0 and a variance of $\sigma^2$.