Practice Set B: Problems 10-13

10.

(a)

A1 = [3 4 5; 2 -3 7; 1 -6 1]; b = [2; -1; 3];

format short; x = A1\b

x =

2.6196

-0.2283

-0.9891

A1*x

ans =

2.0000

-1.0000

3.0000

(b)

A2 = [3 -9 8; 2 -3 7; 1 -6 1]; b = [2; -1; 3];

x = A2\b

Warning: Matrix is close to singular or badly scaled.

Results may be inaccurate. RCOND = 4.189521e-018.

x =

-6.0000

-1.3333

1.0000

The matrix A2 is singular. In fact

det(A2)

ans =

0

(c)

A3 = [1 3 -2 4; -2 3 4 -1; -4 -3 1 2; 2 3 -4 1]; b3 = [1; 1; 1; 1];

x = A3\b3

x =

-0.5714

0.3333

-0.2857

0

A3*x

ans =

1.0000

1.0000

1.0000

1.0000

(d)

syms a b c d x y u v;

A4 = [a b; c d]; A4\[u; v]

ans =

[ -(-d*u+b*v)/(d*a-b*c)]

[ (-c*u+v*a)/(d*a-b*c)]

det(A4)

ans =

d*a-b*c

The determinant of the coefficient matrix is the denominator in the answer. So one gets an answer only if the coefficient matrix is non-singular.

11.

(a)

rank(A1)

ans =

3

rank(A2)

ans =

2

rank(A3)

ans =

4

rank(A4)

ans =

2

MATLAB implicitly assumes ad - bc ¹ 0 here.

(b)

Only the second one computed is singular.

(c)

det(A1)

inv(A1)

ans =

92

ans =

0.4239 -0.3696 0.4674

0.0543 -0.0217 -0.1196

-0.0978 0.2391 -0.1848

det(A2)

ans =

0

The matrix A2 does not have an inverse.

det(A3)

inv(A3)

ans =

294

ans =

0.1837 -0.1531 -0.2857 -0.3163

0 0.1667 0 0.1667

0.1633 0.0306 -0.1429 -0.3367

0.2857 -0.0714 0 -0.2143

det(A4)

inv(A4)

ans =

d*a-b*c

ans =

[ d/(d*a-b*c), -b/(d*a-b*c)]

[ -c/(d*a-b*c), a/(d*a-b*c)]

12.

(a)

[U1, R1] = eig(A1)

U1 =

Columns 1 through 2

-0.9749 0.6036

-0.2003 0.0624 + 0.5401i

0.0977 -0.5522 + 0.1877i

Column 3

0.6036

0.0624 - 0.5401i

-0.5522 - 0.1877i

R1 =

Columns 1 through 2

3.3206 0

0 -1.1603 + 5.1342i

0 0

Column 3

0

0

-1.1603 - 5.1342i

A1*U1 - U1*R1

ans =

1.0e-014 *

Columns 1 through 2

0.3109 0.2776 - 0.3997i

-0.0444 0 - 0.0777i

-0.0833 0.0999 - 0.1776i

Column 3

0.2776 + 0.3997i

0 + 0.0777i

0.0999 + 0.1776i

This is zero. Notice the "e-014".

[U2, R2] = eig(A2)

U2 =

Columns 1 through 2

0.9669 0.7405

0.1240 0.4574 - 0.2848i

-0.2231 0.2831 + 0.2848i

Column 3

0.7405

0.4574 + 0.2848i

0.2831 - 0.2848i

R2 =

Columns 1 through 2

-0.0000 0

0 0.5000 + 6.5383i

0 0

Column 3

0

0

0.5000 - 6.5383i

A2*U2 - U2*R2

ans =

1.0e-014 *

Columns 1 through 2

-0.2224 -0.2554 - 0.2665i

-0.1498 0.0888

-0.1156 -0.0222 + 0.0666i

Column 3

-0.2554 + 0.2665i

0.0888

-0.0222 - 0.0666i

Same comment as in (a).

[U3, R3] = eig(A3)

U3 =

Columns 1 through 2

-0.2446 - 0.4647i -0.2446 + 0.4647i

0.6254 0.6254

0.0025 + 0.3017i 0.0025 - 0.3017i

-0.1736 - 0.4603i -0.1736 + 0.4603i

Columns 3 through 4

-0.5621 + 0.1062i -0.5621 - 0.1062i

0.1982 + 0.0654i 0.1982 - 0.0654i

-0.5833 -0.5833

-0.2215 - 0.4898i -0.2215 + 0.4898i

R3 =

Columns 1 through 2

4.0755 + 4.1517i 0

0 4.0755 - 4.1517i

0 0

0 0

Columns 3 through 4

0 0

0 0

-1.0755 + 2.7440i 0

0 -1.0755 - 2.7440i

A3*U3 - U3*R3

ans =

1.0e-014 *

Columns 1 through 2

-0.1554 - 0.2665i -0.1554 + 0.2665i

0.4441 + 0.0888i 0.4441 - 0.0888i

-0.0888 + 0.2220i -0.0888 - 0.2220i

0.0888 - 0.4441i 0.0888 + 0.4441i

Columns 3 through 4

-0.1887 - 0.0444i -0.1887 + 0.0444i

0.3608 + 0.3331i 0.3608 - 0.3331i

0.2442 - 0.0444i 0.2442 + 0.0444i

-0.3775 - 0.0666i -0.3775 + 0.0666i

And ditto yet again.

[U4, R4] = eig(A4)

U4 =

[ 1, 1]

[ -(-1/2*d+1/2*a-1/2*(d^2-2*d*a+a^2+4*b*c)^(1/2))/b, -(-1/2*d+1/2*a+1/2*(d^2-2*d*a+a^2+4*b*c)^(1/2))/b]

R4 =

[ 1/2*d+1/2*a+1/2*(d^2-2*d*a+a^2+4*b*c)^(1/2), 0]

[ 0, 1/2*d+1/2*a-1/2*(d^2-2*d*a+a^2+4*b*c)^(1/2)]

A4*U4 - U4*R4

ans =

[ 0, 0]

[ c-d*(-1/2*d+1/2*a-1/2*(d^2-2*d*a+a^2+4*b*c)^(1/2))/b+(-1/2*d+1/2*a-1/2*(d^2-2*d*a+a^2+4*b*c)^(1/2))/b*(1/2*d+1/2*a+1/2*(d^2-2*d*a+a^2+4*b*c)^(1/2)), c-d*(-1/2*d+1/2*a+1/2*(d^2-2*d*a+a^2+4*b*c)^(1/2))/b+(-1/2*d+1/2*a+1/2*(d^2-2*d*a+a^2+4*b*c)^(1/2))/b*(1/2*d+1/2*a-1/2*(d^2-2*d*a+a^2+4*b*c)^(1/2))]

simplify(ans)

ans =

[ 0, 0]

[ 0, 0]

(b)

A = [1 0 2; -1 0 4; -1 -1 5];

clear U1 U2 R1 R2

[U1, R1] = eig(A)

U1 =

-0.8165 0.5774 0.7071

-0.4082 0.5774 -0.7071

-0.4082 0.5774 0.0000

R1 =

2.0000 0 0

0 3.0000 0

0 0 1.0000

B = [5 2 -8; 3 6 -10; 3 3 -7];

[U2, R2] = eig(B)

U2 =

0.8165 -0.5774 0.7071

0.4082 -0.5774 -0.7071

0.4082 -0.5774 0.0000

R2 =

2.0000 0 0

0 -1.0000 0

0 0 3.0000

We observe that the columns of U1 are negatives of the corresponding columns of U2. Finally,

A*B - B*A

ans =

0 0 0

0 0 0

0 0 0

13.

(a)

If we set X_n to be the column matrix with entries x_n, y_n, z_n, and M the square matrix with entries 1, 1/4, 0; 0, 1/2, 0; 0, 1/4, 1 then X_n+1 = MX_n.

(b)

We have X_n = MX_n-1 = M²X_n-2 = … = MⁿX₀.

(c)

M = [1, 1/4, 0; 0, 1/2, 0; 0, 1/4, 1];

[U,R] = eig(M)

U =

1.0000 0 -0.4082

0 0 0.8165

0 1.0000 -0.4082

R =

1.0000 0 0

0 1.0000 0

0 0 0.5000

(d)

M should be URU ^-1. Let's check:

M - U*R*inv(U)

ans =

0 0 0

0 0 0

0 0 0

It is totally evident that R_¥ is the diagonal matrix with entries 1, 1, 0. Therefore M_¥ = UR_¥U ^-1. So

Minf = U*diag([1, 1, 0])*inv(U)

Minf =

1.0000 0.5000 0

0 0 0

0 0.5000 1.0000

(e)

syms x0 y0 z0; X0 = [x0; y0; z0]; Minf*X0

ans =

[ x0+1/2*y0]

[ 0]

[ 1/2*y0+z0]

Half of the mixed genotype migrates to the dominant genotype and the other half of the mixed genotype migrates to the recessive genotype. These are added to the two original pure types, whose proportions are preserved.

(f)

M^5*X0

ans =

[ x0+31/64*y0]

[ 1/32*y0]

[ 31/64*y0+z0]

M^10*X0

ans =

[ x0+1023/2048*y0]

[ 1/1024*y0]

[ 1023/2048*y0+z0]

(g)

If you use the suggested alternate model, then only the first three columns of the table are relevant, the transition matrix M becomes M = [1 1/2 0; 0 1/2 1; 0 0 0], and we leave it to you to compute that the eventual population distribution is [1 0 0], independent of the initial population.