Ch 2.6 Case Study: Lake Burley Griffin

Lake Burley Griffin

title

Background

  • Lake Burley Griffin is located in Canberra, Australia's capital.

title

Background

  • Lake Burley Griffin was created artificially in 1962

title

Background

  • Serves recreational and aesthetic purposes.

title

Background

  • Queanbeyan is in the state of New South Wales (NSW)
  • Lake Burley Griffin is in the Australian Capital Territory
  • Safety levels/standards are different between territories.

title

Background

  • Extensive measurements of pollution levels (1970s)
  • Upstream sewage plants (Queanbeyan) contribute pollution

title

Background

  • Significant also was rural and urban runoff during summer rainstorms.
  • These contributed to dramatic increases in pollution levels
  • Responsible for lifting concentration levels above safety limits.

title

Background

  • In 1974 the mean concentration of the bacteria faecal coliform count was approximately \( 10^7 \) bacteria per \( \mathrm{m}^3 \) at the point where the river feeds into the lake.

title

Background

Safety threshold for faecal coliform count in water for contact recreational sports:

  • No more than 10% of total samples over a 30-day period should exceed \( 4 \times 10^6 \) bacteria per \( \mathrm{m}^3 \).

title

Use Model to Investigate

  • If sewage management were improved, would lake would flush out?
  • Would pollution levels would drop below the safety threshold?

title

General Compartmental Model & Word Equations

  • The lake is the only compartment.
  • Applying the balance law, we have the following.

title title

Assumptions, Model & Solution

  • Lake is continuously well mixed
  • Pollution is uniform throughout
  • Lake has a constant volume \( V = 28 \times 10^6 \mathrm{m}^3 \).
  • Lake has constant flow F, with \( F_{in} = F_{out} \).
  • IVP and solution:

\[ \begin{align*} \frac{dC}{dt} & = \frac{F}{V}c_{in} - \frac{F}{V}C(t) , \,\, C(0)= c_0 \\ \\ C(t) & = c_{in} - (c_{in} - c_0) e^{-Ft/V } \end{align*} \]

Lake Pollution Model

  • How long will it take for the lake's pollution level to reach a specified level \( C_L \) if only fresh water flows into the lake?
  • Only fresh water coming in: \( c_{in} = 0 \)

\[ \begin{align*} C(t) & = c_{in} - (c_{in} - c_0) e^{-Ft/V } \\ & = c_0 e^{-Ft/V } \end{align*} \]

  • Solve for \( t \), using \( C = C_L \):

    \[ t = - \frac{V}{F} \ln \left( \frac{C_L}{c_0} \right) \]

Lake Pollution Model: Predictions

title

  • Fresh water entering the lake: \( c_{in} = 0 \)
  • Volume: \( V = 28 \times 10^6 \mathrm{m}^3 \).
  • Mean monthly flow: \( F = 4 \times 10^6 \mathrm{m}^3/ \mathrm{month} \)
  • Initial faecal coliform count: \( c_0 = 10^7\mathrm{bacteria}/ \mathrm{m}^3 \)
  • \( C_L = ? \)

Lake Pollution Model: Predictions

title

  • Fresh water entering the lake: \( c_{in} = 0 \)
  • Volume: \( V = 28 \times 10^6 \mathrm{m}^3 \).
  • Mean monthly flow: \( F = 4 \times 10^6 \mathrm{m}^3/ \mathrm{month} \)
  • Initial faecal coliform count: \( c_0 = 10^7\mathrm{bacteria}/ \mathrm{m}^3 \)
  • \( C_L = 4 \times 10^6 \mathrm{bacteria}/ \mathrm{m}^3 \)

Lake Pollution Model: Predictions

  • Parameter Values
  • \( V = 28 \times 10^6 \mathrm{m}^3 \)
  • \( F = 4 \times 10^6 \mathrm{m}^3/ \mathrm{month} \)
  • \( C_L = 4 \times 10^6\mathrm{bacteria}/ \mathrm{m}^3 \)
  • \( c_0 = 10^7\mathrm{bacteria}/ \mathrm{m}^3 \)

  • Formula \[ t = - \frac{V}{F} \ln \left( \frac{C_L}{c_0} \right) \]

  • Use R for calculations:
V <- 28*10^6 
F <- 4*10^6 
CL <- 4*10^6
c0 <- 10^7
t <- - V/F*log(CL/c0)

t

[1] 6.414035

Lake Pollution Model

  • However, pure water entering the lake is not a very realistic scenario with three sewage plants and much farmland upstream.
  • Thus including the entrance of polluted river water into the lake model is essential.

\[ C(t) = c_{in} - (c_{in} - c_0) e^{-Ft/V } \]

title

Numerical Solutions: Maple, MATLAB, R

  • Including the entrance of polluted river water into the lake model is essential.

\[ \frac{dC}{dt} = \frac{F}{V} \left[ c_{in} - C(t) \right] , \,\, C(0)= c_0 \]

  • Case 1: \[ \begin{align*} c_{in} & = 3 \times 10^6 \mathrm{bacteria}/ \mathrm{m}^3 \\ F & = (4*12) \times 10^6 \mathrm{m}^3/ \mathrm{year} \end{align*} \]
  • Case 2: \[ \begin{align*} c_{in} & = 10^6 \left[ 10 + 10 \cos(2 \pi t) \right] \mathrm{bacteria}/ \mathrm{m}^3 \\ F(t) & = 10^6 \left[ 10 + 6 \sin(2 \pi t) \right] \mathrm{m}^3/ \mathrm{year} \end{align*} \]

Case 1: Maple (Ch2.5 Listing 2.7, Part 1)

title

Case 1: Maple (Ch2.5 Listing 2.7, Part 2)

title

Case 1: Maple (Ch2.5 Listing 2.7)

\[ \lim_{t \rightarrow \infty} C(t) = \lim_{t \rightarrow \infty} \left( c_{in} - (c_{in} - c_0) e^{-Ft/V }\right) = c_{in} \]

title

Case 1: MATLAB (Ch2.5 Listing 2.8)

title

\[ \begin{align*} c_{in} & = 3 \times 10^6 \mathrm{bacteria}/ \mathrm{m}^3 \\ F & = (4*12) \times 10^6 \mathrm{m}^3/ \mathrm{year} \end{align*} \]

Case 1: MATLAB (Ch2.5 Listing 2.8)

  • ch25Lake

title

Case 1: R (Fourth-Order Runge-Kutta)

rk4plot <- function (f, x0 , y0 , h, n) {
  x <- rep(0, n)
  y <- rep(0, n)
  x[1] <- x0
  y[1] <- y0
  for(i in 1:n) {
    s1 <- h*f(x[i], y[i])
    s2 <- h*f(x[i] + h/2, y[i] + s1/2)
    s3 <- h*f(x[i] + h/2, y[i] + s2/2)
    s4 <- h*f(x[i] + h, y[i] + s3)
    y[i+1] <- y[i] + s1/6 + s2/3 + s3/3 + s4/6
    x[i+1] <- x[i] + h
  }
  return(plot(x,y,type="l", col="blue"))
}

Case 1: R (Fourth-Order Runge-Kutta)

f <- function(x, C) { 48/28*(3 - C) }
rk4plot(f, 0, 10, 0.1, 40)

plot of chunk unnamed-chunk-4

Case 1: R (Fourth-Order Runge-Kutta)

rk4data <- function (f, x0 , y0 , h, n) {
  x <- rep(0, n)
  y <- rep(0, n)
  x[1] <- x0
  y[1] <- y0
  for(i in 1:n) {
    s1 <- h*f(x[i], y[i])
    s2 <- h*f(x[i] + h/2, y[i] + s1/2)
    s3 <- h*f(x[i] + h/2, y[i] + s2/2)
    s4 <- h*f(x[i] + h, y[i] + s3)
    y[i+1] <- y[i] + s1/6 + s2/3 + s3/3 + s4/6
    x[i+1] <- x[i] + h
  }
  return(data.frame (x = x, C = y))
}

Case 1: R (Fourth-Order Runge-Kutta)

f <- function(x, C) { 48/28*(3 - C) }
rk4data(f, 0, 10, 0.2, 10)

     x         C
1  0.0 10.000000
2  0.2  7.968438
3  0.4  6.526483
4  0.6  5.503016
5  0.8  4.776583
6  1.0  4.260978
7  1.2  3.895013
8  1.4  3.635260
9  1.6  3.450893
10 1.8  3.320033
11 2.0  3.227152

Case 2: Maple, Part 1 (Ch2.6 DrG)

title

Case 2: Maple, Part 2 (Ch2.6 DrG)

title

Case 2: MATLAB (Ch2.6 DrG)

title

Case 2: MATLAB (Ch2.6 DrG)

  • ch26Lake2

title

Case 2: MATLAB (Ch2.6 DrG)

title

Case 2: MATLAB (Ch2.6 DrG)

  • ch26Lake(6)

title

Case 2: R (Fourth-Order Runge-Kutta)

F1 <- function(x) {(10 + 6*sin(2*pi*x))}
c1 <- function(x) {(10 + 10*cos(2*pi*x))}
f <- function(x, C) { F1(x)/28*(c1(x) - C) }
rk4plot(f, 0, 6, 0.1, 80)

plot of chunk unnamed-chunk-7

Case 2: R (Fourth-Order Runge-Kutta)

F1 <- function(x) {(10 + 6*sin(2*pi*x))}
c1 <- function(x) {(10 + 10*cos(2*pi*x))}
f <- function(x, C) { F1(x)/28*(c1(x) - C) }
rk4data(f, 0, 6, 0.1, 10) 

     x        C
1  0.0 6.000000
2  0.1 6.550003
3  0.2 7.020280
4  0.3 7.183187
5  0.4 7.032754
6  0.5 6.770244
7  0.6 6.592136
8  0.7 6.545226
9  0.8 6.595532
10 0.9 6.768673
11 1.0 7.132961

Case 2: R (Fourth-Order Runge-Kutta)

title

title