Project 1.7. Great Lakes Pollution

Understanding

Most of the water flowing into Lake Ontario is from Lake Erie. Suppose that the pollution of  these two lakes ceased suddenly.  How long would it take for the pollution level in eaach lake to be reduced to 10% of its present level?

Plan

We will model the system using two discrete deterministic recurrence relations, which are coupled. We will make a number of simplifying assumptions:
	1) all water that leaves Erie goes into Ontario, and Ontario gets water only from Erie;
	2) water leaving Erie is replaced, so its water level is constant (as is Ontario's);
	3) pollution is well-mixed in the lake (i.e. uniform pollution levels);
	4) pollution is not longer entering Erie, nor Ontario (except via the polluted waters of Erie).
	
Furthermore, we are told that
	1) each year 38% of Erie's water flows into Ontario;
	2) each year 13% of Ontario's water flows out and is replaced by the water from Erie;
	3) a(n) and b(n) are the total amounts of pollution in Ontario and Erie, respectively. (This is not 	     the same as concentration!)
	4) we are to assume that the initial pollution contents are a(0)=b(0)=1.

Part 1: From these assumptions we can formulate the following model for the pollution levels as a function of time given by the variable:

	a(n+1) = 0.87a(n) + .38b(n)
	b(n+1)=0.62b(n)

for Ontario and Erie, respectively. 
	The latter says that the pollution content of Erie at time n+1 is 62% of what it was in the previous year   (b(n)),  since 38% of the pollution  moved into Ontario.
	The former says that Ontario's pollution content is increased each year by 38% of Erie's total from the previous year, and that it loses 13% of its previous year's pollution load (a(n)).

Part 2: 
	a(n+2) = .87a(n+1) + .38b(n+1)
	a(n+1) = .87a(n) + .38b(n)
but since there is a relationship between b(n+1) and b(n),
	b(n+1)=.62b(n)
we deduce that
	a(n+2) = .87a(n+1) + .38*.62b(n)
	a(n+1) = .87a(n) + .38b(n)
and since we have two equations involving b(n), we can eliminate b(n) by subtracting 0.62 times the second from the first:
	a(n+2) = .87a(n+1) + .38*0.62b(n)
	.62*a(n+1) = .62*.87a(n) + .62*.38b(n)
to give
	a(n+2) - .62*a(n+1) = 0.87a(n+1) + .38*0.62b(n) - [.62*0.87a(n) + .62*.38b(n)]
which yields
	a(n+2) = .87a(n+1)  + .62*a(n+1) - .62*.87a(n) 
	            = (.87  + .62)a(n+1) - .62*0.87a(n)

Part 3: 
	Since a(0)=1 and b(0)=1, we can find a(1):
	a(1) = .87a(0) + .38b(0) = 1.25
	
	a(n+2) =  (.87  + .62)a(n+1) - .62*0.87a(n)  
is a linear equations of second order, so we can solve it using "The Method" of page 26:
We assume a solution of the form
	a(n)=Ck^n ,
plug it into the linear equation, and solve for values of  k that make it work. This gives
k^2 -(.87  + .62)k+.62*0.87=0,
which has solutions visible to the naked eye k=.87 and k=.62. The general solution is a sum of these two solutions, i.e.
	a(n)=C_1 .87^n + C_2 .62^n
where the constants are determined by meeting the initial conditions (i.e. the values of a(0) and a(1)):
	a(0)=C_1 + C_2=1
	a(1)= .87C_1 +  .62C_2=1.25
When we solve these two equations for C_1 and C_2, we find that
	a(n)=2.52(.87)^n - 1.52(.62)^n
Checking these values below, we confirm that they are correct.

Carry Out the Plan

We set out the equations and the initial conditions:

Clear[Erie,Ontario,n,Years];

Years=30;
Erie[0]=1.0;
Ontario[0]=1.0;

Erie[n_]:=  Erie[n] =  Erie[n-1]*0.62;
Ontario[n_]:=Ontario[n] = .87*Ontario[n-1] + .38*Erie[n-1];

Here's our simulation, including both the recurrence relation equations and our closed-form solutions:

ErieVals=Table[{n,Erie[n]},{n,0,Years}];
OntarioVals=Table[{n,Ontario[n]},{n,0,Years}];
EventualVals=Table[{n,2.52(.87)^n},{n,0,Years}];
BothVals=Table[{n,
				Erie[n],
				Ontario[n],
				(.62)^n,
				2.52(.87)^n - 1.52(.62)^n},
				{n,0,Years}];

A comparison of the values obtained from the analytic solutions and the recurrence relations shows that they're the same:

BothVals//TableForm

Let's plot the result:

pplot=ListPlot[ErieVals,PlotJoined->True,AxesLabel->{"Years",  "Total Pollution"}];
dplot=ListPlot[OntarioVals,PlotJoined->True,AxesLabel->{"Years",  "Total Pollution"}];
eplot=ListPlot[EventualVals,PlotJoined->True,AxesLabel->{"Years",  "Eventual Pollution"}];

\(Show[pplot, dplot]; 
Show[pplot, dplot, \ eplot]; \)

Evaluate

Erie falls to 90% of its initial pollution level within the first year. Ontario takes several years to do so, in year 7.
	Erie reaches 10% of its initial pollution level by the end of year 5. Ontario takes 24 years to do so.
	While this model is rather simplistic, it still shows a number of interesting features. For one thing, it suggests that, given a cessation of source pollution, a lake might recover relatively quickly (Erie). Erie's water source is clean, and Erie recovers at an exponential rate (the level of pollution declines exponentially).
	The model as it stands starts with lake Erie with a higher concentration of pollution than Ontario (their gross quantities of pollution are both 1, but Erie contains much less water than Ontario, so its concentration is much greater). Since the only source to Ontario is taken to be Erie, the water pouring into Ontario is initially more polluted than that exiting, so the pollution level of Ontario climbs initially.
	Then Ontario too declines, and does so exponentially. This is the result of the "competition" between the eigenvalues .87 and .62: eventually powers of .62 are much smaller than powers of .87. So eventually the solution becomes indistinguishable from the values of 2.52(.87)^n (shown in the plots above as EventualVals). As the curve rises to meet its "eventual curve", it is experiencing "transient behavior".
	
	The model is obviously neglecting a lot, including the myriad sources of pollution that afflict a huge body of water with many bordering cities, passing through two countries, lots of agriculture bordering its shores, etc. The book mentions neglecting evaporation and rainfall (e.g. acid rainfall in Canada and the northeast). The assumption that Erie is the only source of water for Ontario is highly suspect, of course. The mixture assumption is absurd. Even so, these sorts of assumptions are the kind that permit us to get a first approximation to a solution which is perhaps more reasonable; and in the meantime, the first approximation may be somewhat meaningful.
	
	For example, this model would  also serve for the assumption of a sudden source of pollution entering Ontario with an equilibrium pollution profile (i.e., on top of a steady state of pollution). Rather than tailing off to zero pollution, Ontario would recover to its steady state in the same fashion as we see here (exponentially).