kill(all)$

peter.vlasschaert@gmail.com,19/12/2016

 general equation (second-order linear equations)

discrete : second-order linear equation( => alternate approach)
-------------------------------------------------------------
p1*a(n+2)+p2*a(n+1)+p3*a(n)=0


amount pollution => simplification to understand principles.

Lake Ontari:a(n)
Lake Eri   :b(n)
Lake Ontario connected with Lake Eri

system of equations :difference equations

           a(n+1)=0.87*a(n)+0.38*b(n)    '1'
           b(n+1)=0.62*b(n)              '2'

eqn1:a(n+1)=a1*a(n)+b1*b(n)$
eqn2:b(n+1)=a2*a(n)+b2*b(n)$

b(n):each year: Lake Eri 38%(0.38) 'replaced' by rain          
      

a2:0;
b2:1-0.38;

eqn1:a(n+1)=a1*a(n)+b1*b(n);
eqn2:b(n+1)=a2*a(n)+b2*b(n);

rem:replaced : 1-..

a1:1-0.13;
b1:0.38;

eqn1:a(n+1)=a1*a(n)+b1*b(n);
eqn2:b(n+1)=a2*a(n)+b2*b(n);

eqn_1:part(solve(eqn1,b(n)),1);

eqn_2:ev(eqn_1,n=n+1);

eqn_3:ev(eqn2,eqn_1,eqn_2);

p1*a(n+2)+p2*a(n+1)+p3*a(n)=0

eqn_4:lhs(eqn_3)-rhs(eqn_3)=0;

(0.62/0.38)/100;

a(n+2) = x^2
a(n+1) = x
a(n)   = 1
p1*x^2+p2*x+p3 = 0, solve ?

eqn_5:ev(eqn_4,a(n)=1,a(n+1)=x,a(n+2)=x^2);

eqs:solve(eqn_5,x);

31.0/50.0;

pc1:float(part(eqs,1,2));
pc2:float(part(eqs,2,2));

eq_model:a(k)=c1*(pc1)^k+c2*(pc2)^k;

given : a(0),a(1) to find c1,c2

eq_model1:ev(eq_model,k=0);
eq_model2:ev(eq_model,k=1);

solve the system : a(0),a(1)

cf:solve([eq_model1,eq_model2],[c1,c2]);

cf1:part(cf,1,1);
cf2:part(cf,1,2);

model_lake:ev(eq_model,cf1,cf2);