Section 4. Matrix Methods

3.* Define the three column vectors b1, b2, b3 as follows

b1,b2,b3 := <0,1,-1>,<1,1,0>,<-1,0,1>;

> with(LinearAlgebra):

> b1,b2,b3 := <0,1,-1>,<1,1,0>,<-1,0,1>;

b1, b2, b3 := Vector[column]([[0], [1], [-1]]), Vector[column]([[1], [1], [0]]), Vector[column]([[-1], [0], [1]])

a. Define the matrix B having b1, b2, b3 as its columns with the entry B := <b1|b2|,b3>.

> B := <b1|b2|b3>;

B := Matrix([[0, 1, -1], [1, 1, 0], [-1, 0, 1]])

b. Find the characteristic polynomial of B and factor it with factor(%).

> CharacteristicPolynomial(B,t);
factor(%);

t^3-2*t^2-t+2

(t-1)*(t-2)*(t+1)

c. Find B's eigenvalues with Eigenvalues(B).

> Eigenvalues(B);

Vector[column]([[-1], [1], [2]])

d. Find B's eigenvectors using Eigenvectors(B). Name the input lambda,V.

> lambda,V := Eigenvectors(B);

lambda, V := Vector[column]([[1], [2], [-1]]), Matrix([[0, -1, 2], [1, -1, -1], [1, 1, 1]])

e. Calculate MatrixInverse(V).B.V.

> MatrixInverse(V).B.V;

Matrix([[1, 0, 0], [0, 2, 0], [0, 0, -1]])

4.* Continuing 3. Obtain the solution to v' = Bv satisfying v(0) = <1,2,3>.

a. Do if first using dsolve applied to the system of linear differential equations defined by v' = Bv with the appropriate initial conditions.

> u := <x(t),y(t),z(t)>:
DEsystem := 'diff(u[k],t) = (B.u)[k]' $ k=1..3;

DEsystem := diff(x(t), t) = y(t)-z(t), diff(y(t), t) = x(t)+y(t), diff(z(t), t) = -x(t)+z(t)

> soln := dsolve( {DEsystem,x(0)=1,y(0)=2,z(0)=3} );

soln := {z(t) = 1/2*exp(-t)+5/2*exp(t), x(t) = exp(-t), y(t) = -1/2*exp(-t)+5/2*exp(t)}

> VectorForm = subs(soln,u);

VectorForm = Vector[column]([[exp(-t)], [-1/2*exp(-t)+5/2*exp(t)], [1/2*exp(-t)+5/2*exp(t)]])

> Check;
eval(%%,t=0);

eval(DEsystem,soln);

Check

VectorForm = Vector[column]([[1], [2], [3]])

-exp(-t) = -exp(-t), 1/2*exp(-t)+5/2*exp(t) = 1/2*exp(-t)+5/2*exp(t), -1/2*exp(-t)+5/2*exp(t) = -1/2*exp(-t)+5/2*exp(t)

b. Do it second by making a fundamental matrix solution X(t) defined as the matrix with the eigenvector solutions in the columns then computing

v(t) = X(t)*X(0)^(-1)*v(0)

> for k from 1 to 3 do W[k] := exp(lambda[k]*t)*Column(V,k): end do: unassign('k');

> <W[1]|W[2]|W[3]>;
X := unapply(%,t):

Matrix([[0, -exp(2*t), 2*exp(-t)], [exp(t), -exp(2*t), -exp(-t)], [exp(t), exp(2*t), exp(-t)]])

> v := X(t).X(0)^(-1).<1,2,3>;

v := Vector[column]([[exp(-t)], [-1/2*exp(-t)+5/2*exp(t)], [1/2*exp(-t)+5/2*exp(t)]])

> Check;
eval(v,t=0);

map(diff,v,t) = B.v;

Check

Vector[column]([[1], [2], [3]])

Vector[column]([[-exp(-t)], [1/2*exp(-t)+5/2*exp(t)], [-1/2*exp(-t)+5/2*exp(t)]]) = Vector[column]([[-exp(-t)], [1/2*exp(-t)+5/2*exp(t)], [-1/2*exp(-t)+5/2*exp(t)]])

c. To it third by using the Matrix Exponential, Exp_At. Once you have it, the solution is

v(t) = Exp_At*v(0)

> Exp_At := MatrixExponential(B,t);

Exp_At := Matrix([[2/3*exp(-t)+1/3*exp(2*t), 1/3*exp(2*t)-1/3*exp(-t), -1/3*exp(2*t)+1/3*exp(-t)], [1/3*exp(2*t)-1/3*exp(-t), 1/6*exp(-t)+1/3*exp(2*t)+1/2*exp(t), -1/3*exp(2*t)+1/2*exp(t)-1/6*exp(-t)]...

> v := Exp_At.<1,2,3>;

v := Vector[column]([[exp(-t)], [-1/2*exp(-t)+5/2*exp(t)], [1/2*exp(-t)+5/2*exp(t)]])

> Check;
See_above;

Check

See_above

d. Check the solution in each case.