Chemostat modell

Implement a chemostat model with one resource and Monod kinetics.

library("deSolve")
chemostat <- function(time, y, parms) {
  with(as.list(c(y, parms)), {
    mu   <- mumax * P/(kp + P)
    dAlg <- mu * Alg - D * Alg
    dP   <-  D *(P0 - P) - 1/Y * mu * Alg
    list(c(dAlg, dP), mu=mu)
   })
}
parms <- c(
  mumax = 0.5,    # 1/d
  kp    = 0.01,   # half saturation constant (mg/L)
  Y     = 41,     # Yield coefficient
  D     = 0.1,    # Dilution rate (1/d)
  P0    = 0.05    # Phosphorus in the medium (mg/L)
)
times <- seq(0, 40, 0.1)  # (d)
y0    <- c(Alg=2, P=0.05) # Phytoplankton C and Phosphorus P (mg/L)
out   <- ode(y0, times, chemostat, parms)
plot(out)

Exercises

  1. Which measurement unit has \(Y\)?
  2. Is it possible to turn the chemostat model into a batch, by changing only the parameters, not the equations?
  3. Does the equilibrium after a reasonably long time depend on the selection of the initial conditions of Alg and P? Are there exceptions?

Dependence of the equilibrium on the dilution rate

Run simulations with different dilution rates D ans plot \(P\) and \(Alg\) as a function of \(D\).

Protocol

  1. Find a reasonable range for D by experimenting
  2. Subdivide this range in 10 equidistant steps.
  3. Run the model, write the last values down
  4. Enter the values of \(D\), \(P\) and \(Alg\) in a data frame or 3 vectors (or into Excel).
  5. Calculate \(D * Alg\)
  6. Plot \(P\), \(Alg\) and \(D * Alg\) as a function of \(D\)

As an alternative, you may also consider to write a loop in R with more than 10 intermediate steps.

Hints

  • make sure that the model is always run to an equilibrium
  • the final state of a simulation can be retrieved with:
out[nrow(out),]
  • as an alternative, package rootSolve may be used:
library(rootSolve)
times <- c(0, Inf)
out   <- stode(y0, times, chemostat, parms)
out
  • rootSolve has two solvers, stode (that solves the system with a Newton-Raphson method by setting the derivatives to zero) and runsteady that runs the model numerically until an equilibrium is reached.

For the chemostat model, an analytical solution would also be possible, i.e. mathematics with pencil and paper instead of the computer. Though this is in gereral the preferred method, we use simulation here as a didactical reasons, because an analytical solution is ofteh imposible for complex ecological models.