Load the Linear Algebra package as usual.
Additionally we will need to load two new packages for today:
> with(plots); with(ArrayTools);

Q1
We are going to be looking at the recurrence such that
a(k+1)=2*a(k)+13*a(k-1)-38*a(k-2)+24*a(k-3) and a(0)=10, a(1)=79, a(2)=309 and a(3)=835.
We set up the matrix representing this recurrence as follows:
 > R:=<<2,1,0,0>|<13,0,1,0>|<-38,0,0,1>|<24,0,0,0>>;

Verify that the eigenvalues and eigenvectors have the properties discussed in class.
That is, the eigenvalues are roots of the polynomial p defined below 
and the eigenvectors are made from the powers of each eigenvalue. 
 > (e,v):=Eigenvectors(R);
 > p:=x -> x^4 - 2*x^3 - 13*x^2 + 38*x - 24;

You can evaluate and/or plot the polynomial as follows but then use words to explain.
 > p(5);
 > plot(p(x),x=-5..5,-150..150);
[1 mark]

We can get the general product of the power of the matrix and the initial values as follows:
> Multiply(MatrixPower(R,n),<835,309,79,10>);
All four lines give formulae for a term of the recurrence. Identify which row gives
the correct values for a(0), a(1), a(2) and a(3) and check it.
[1 mark]

Since the dominant eigenvalue is negative, this sequence will eventually be negative for some
values of n. Increase the value that is 8 in this plot command until you see that first happen.
> pointplot([seq( [n,(-4)^n+31*3^n+12*2^n- 34],n=0..8)],symbol=solidcircle,color=red,symbolsize=20);

Maple can verify when (4/3)^n > 31 to see when the two largest eigenvalue powers will cancel out:
> evalf(log(31)/log(4/3));
Why does this value not match that found using the plot?
[1 mark]

The eigenvector matrix is of a type called a Vandermonde Matrix, Maple can produce one as follows:
> B:=Reverse(Transpose(VandermondeMatrix(<1,2,3>,3,4)));
Use B to find the value of a(3) (which is not 835) which will guarantee that the sequence increases
for ever. To do this you will need to find t using Pivot and/or RowOperation on
> Bv:=<B | <t,309,79,10>>;
Check your answer using MatrixPower as before.
[2 marks]


Q2
Let us pick some random data to fit lines and curves to:
> xv:=<-3,2,4>; randomize(Value(Time())): yv:=RandomVector(3,generator=rand(-2..6));
xv are the x values and yv are the y values.

We can now visualise the data using this command:
> p1:=pointplot(xv,yv,symbol=solidcircle,symbolsize=20,colour=[green],view=[-5..6, -4..9]);

Now make the matrix that we will use for the best fit straight line:
> A1:=VandermondeMatrix(xv,3,2);
If you prefer to use the ones on the left format you can use:
> A1:=FlipDimension(VandermondeMatrix(xv,3,2),2);
Whichever A1 you choose to use you can multiply the transpose as we have talked about:
> A2:=Multiply(Transpose(A1),A1);
> B2:=Multiply(Transpose(A1),yv);

Use MatrixInverse to find the matrix which when multiplied by A2 gives B2 
and hence find the coefficients of the best fit straight line giving new variables m and b. 
[1 mark]

Having assigned m and b we can plot the line using:
> f:=x -> m*x+b;
> p2:=plot(f(x),x=-5..6,color=blue);
and display both data and straight line using
> display([p1,p2]);
Based on the picture, do you think you have correctly found the line of best fit?

Check that the differences as follows have a sum of zero:
> d1:=f(xv[1])-yv[1]; d2:=f(xv[2])-yv[2]; d3:=f(xv[3])-yv[3];
[2 mark]

Now repeat this process to find and plot the exact fit quadratic q using
> A3:=VandermondeMatrix(xv,3,3);
Find the coefficients r, s and t to make the quadratic q 
> q:=x -> r*x^2 + s*x + t;
and then create a plot p3:
> p3:=plot(q(x),x=-5..6,color=red);

You can display all three plots as follows:
> display([p1,p2,p3]);
How are you convinced you have the exact fit quadratic?
[2 marks]

Once you have completed all of these activities, please save your worksheet then use
File / Export As and choose pdf. Now upload the file to Crowdmark.