Problem of the Day

Aug 26  Here is a talk and a paper that discusses the polygon averaging problem. The matrix M_{n} plays a critical role, e.g.,

M_{4} = .5*[1 1 0 0; 0 1 1 0; 0 0 1 1; 1 0 0 1]

What would the "M" matrix look like if we took an arbitrary convex combination of vertices instead of their average, i.e.,

xNew(k) = mu* x(k) + (1-mu)*x(k+1)               0 < mu < 1

instead of xNew(k) = (x(k)+x(k+1))/2. Feel free to play around with PolyAve.m to see what the impact would be with this generalization. (Visit the computation of  xNew and yNew  in the subfunction Smooth.)

--------------------

Aug 28   This will get you thinking about rows, columns, and blocks and how to "grab 'em" in Matlab. Write a Matlab function Sudoku(A) that returns "1" if A is a valid Sudoku solution and zero otherwise. Thus, A must be a 9x9 matrix with the property that every column, row, and 3x3 block is made up of the integers 1 through 9.. Hints on checking. If sort(A(5,:)) is the same as 1:9, then row 5 is "OK".  If B= A(4:6,7:9) and sort(B(:)) is the same as (1:9)', then the  (2,3)  block is OK. Some built-in functions that you may find useful: sort, sum, any, all

--------------------

Aug 31   This problem is about "spotting" matrix-vector operations given some awful-looking multiple summation with subscripts all over the place. As a warm-up, take a look at problem 1.1.6 in GVL4. Now here is the problem. Assume that A is nxn and that  u, x, y, z, and d are column n-vectors of reals. Assume that v is a permutation of 1:n. Rewrite the following O(n^4)  algorithm so that it is O(n^2):

 alfa = 0 for i = 1:n for j = 1:n for p = 1:n for q = 1:n alfa = alfa + A(i,p)*A(v(j),q)*u(i)*y(j)*x(q)*z(p)*d(i) end end end end

A handy test script SpotMatVec is available.

-------------------

Sep 2   Without using Matlab, what is the SVD of the 2x2 matrix  A =  [ 3   4 ; 8   -6] ?

--------------------

Sep 4  The Haar matrix W_n is defined in terms of W_m where n = 2m and W_{1} = [1]:

W_n = [ kron(W_m , [1;1])  kron(I_{m},[1;-1]) ]

Can you develop a closed form expression for the number of nonzeros in W_n ?

--------------------

Sep 9 Get solid with the floating representation. Understand fpFacts.m. Give a recipe for the largest positive integer k such that 10^k can be represented EXACTLY in IEEE double precision. Hint: How many bits are in the binary expansion of 5^k. Give a recipe for the largest positive integer k such that k! can be represented EXACTLY in IEEE double precision. Hint: 6! = 1*2*3*4*5*6 = 45*(2^4)

--------------------

Sep 11    Assume L1 and L2  are nxn, lower triangular, and nonsingular. Assume that B is nxn. How would you compute an nxn matrix X so that L1*X*L2' = B using the \ operator? Hint. How would you do it if L2 = I? (How would you do it if L1 = 0?) Now assume the availability of a function x = LTriSol(L,b) that returns the solution to the lower triangular system Lx = b where b is a column vector. How would you compute X using LTriSol?

--------------------

Sep 14  A and B are nxn matrices and A is nonsingular. Assume that b and c are column n-vectors. How would you use the Matlab lu function and \ operator to solve the 2n-by-2n linear system [A  B ; zeros(n,n) A' ] [y ; z] = [b ; c]?

--------------------

Sep 16 (a) Assume nxn matrices A1 and A2 are differ only in their last columns. Assume that lu(A1) produces P1*A1 = L1*U1 and that lu(A2) produces P2*A2 = L2*U2. Describe the differences between P1 and P2, L1 and L2, and U1 and U2. (b) Assume  that A is a singular nxn matrix and that PA = LU. Assuming exact arithmetic, how could a unit 2-norm vector z be found so that Az = 0.

--------------------

Sep 18   Modify the function BlockLU so that it reports back the fraction of flops that are level-3 flops. Modify it again so that it is nonrecursive.

--------------------

Sep 21  Assume that we have computed the Cholesky factorization A = GG' of a symmetric positive definite nxn matrix A.  Let k be an integer that satisfies 1<=k<n. Recall that A(k,k) is positive and that A(1:k,1:k) is also sym pos def. What is the smallest smallest positive tau for which A - tau*ek*ek' is not positive definite? Here, ek = I(:,k). In plain English, we want to subtract "just enough" stuff from the (k,k) entry so that A loses definiteness.

--------------------

Sep 23 A is symmetric positive definite and tridiagonal. If A = LDL' and D = diag(d), how would you compute a vector x so that x'*A*x = d(k)?

--------------------

Sep 25  Suppose Tn = Gn*Gn' is the Cholesky factorization of Tn, the -1,2,-1 tridiagonal matrix. Gn is lower bidiagonal. What can you say about the entries Gn(n,n) and Gn(n,n-1) as n-->infty?

--------------------

Sep 28  Under what conditions does the Gauss-Seidel iteration converge when A = [a b; c d]? Repeat for the Jacobi iteration.

--------------------

Sep 30

--------------------

Oct 2  Suppose A = [ d1 e; e d2] is positive definite. What is the condition of this matrix? Now let Atilde = inv(D)Ainv(D) where D = diag(sqrt(d1),sqrt(d2)). What is the condition of Atilde? This is an issue because the latter condition nember determines how quickly PCG converges with D as a preconditioner.

--------------------

Oct 5 (a) For the 3-by-2 case, verify that the method of normal equations is equivalent to zeroing the gradient of

f(x) = (A*x - b)'*(A*x - b)

(b) Given matrices A1 (m1xn), A2 (m2xn) and A3 (m3xn) and vectors b1 (m1x1), b2 (m2x1) and b3 (m3x1), how would you compute a vector x that minimizes || A1*x - b1 ||^2  + || A2*x - b2 ||^2  + ||A3*x - b3 ||^2? (c) Given [Q,R] = qr(A,0), how would you compute the minimum value of || Ax -b|| ?

---------------------

Oct 7   Given A (mxn) and b (mx1), how would you compute the closest vector to b that is orthogonal to every column of A? (b) Given that Q = [c s;s -c] is orthogonal, how would you compute a unit vector u (2x1) so that Q = I-2u*u'?

--------------------

Oct 9  Suppose f and g are column n-vectors. Show how to compute an orthogonal Q as a product of Givens rotations so that Q' * (f * g') * Q  is zero everywhere except in the (1,1) and (1,2) entry. Hint. Put a lot of zeros into f and then put a lot of zeros in g.

--------------------

Oct 14 (a) What can you say about Q, R, and P if [Q,R,P] = qr(A) and A has more columns than rows. (b) Assuming that A has full row rank, how could you use Q, R, and P to compute a solution to  Ax = b?

--------------------

Oct 16  Suppose A is mxn with full column rank and that b is mx1. How would you minimize || Ax - b|| subject to the constraint that x(1) +...+ x(n) = 0?

---------------------

Oct 19 Suppose A and B are are 2x2 and that det(A)det(B) < 0. How would you compute Q = [c s;-s c] so that norm(A - BQ,'fro') is minimum?

--------------------

Oct 21 This is about the case when the TLS problem has no solution. Show that if the smallest singular value of [A b] is unique and equals the smallest singular value of A, then the TLS problem has no solution. Use this result to show that if A'*b = 0 and norm(b) > sigma_min(A), then the TLS problem has no solution.

Mar 13  Let A = [lambda_1  m ; 0 lambda_2] and assume that the eigenvalues lambda_1 and lambda_2 are distinct. Determine X so inv(A)*A*X = diag(lambda_1,lambda_2). Recall that if Ax = lambda* x and y'*A = lambda*y' then the condition of lambda is  1/abs(x'*y) assuming that x and y have unit 2-norm. What is the condition of lambda_1 and lambda_2 in the above 2-by-2 example? More generally, how would you compute the condition number of all the eigenvalues using [X,D] = eig(A) assuming that the eigenvalues of A are distinct?

Mar 9  If A-sI = QR, why is B = RQ +sI  similar to A?

Mar 6  Write a fragment to handle the update A = H'*A*H where H = I - 2*u*u' with u(1:r) = 0. Assume A is symmetric.

Mar 4

Mar  2  Show that after a Jacobi update that zeros A(p,q), then off(Anew)^2 = off(A)^2 - 2A(p,q)^2. Here, off(A)^2 is the sum of the squares of the off-diagonal entries.

--------------------

Oct 23  Suppose Q'*A*Q = diag(d1,d2) is the schur decomposition of a symmetric 2x2 matrix A. Why is it always possible to choose Q so that det(Q) = -1?

--------------------

Oct 26 Suppose Ax = lambda*x is the dominant eigenpair for a symmetric A. What happens if we apply the power method with the matrix A - lambda*x*x' ?

--------------------

Oct 28  If you have a function [c,s] = Schur2(A) that computes a rotation that diagonalizes a symmetric 2-by-2 A, show how it can be used to compute a 2-by-2 svd. Hint: How would you determine c and s so that [c s;-s c]'*[w x ; y z] is symmetric.?

--------------------

Oct 30  Suppose A is nxn and symmetric. How would you compute an nxk matrix Q with orthonormal columns so that trace(Q'*A*Q) is maximized.?

--------------------

Nov 2 Suppose Ax = mu*x + r where x has unit 2-norm. Use the Wielandt-Hoffman theorem to prove that mu is within norm(r,2) of an exact eigenvalue.

--------------------

Nov 4   In the Householder tridiagonalization Q'*A*Q = T, the first column of Q is the first column of I. Explain. Now suppose v is a unit vector.How would you compute an orthogonal U so that (a) Uv = + or - I(:,1) and (b) U'*A*U is tridiagonal.

--------------------

Nov 6  Suppose A = I + UDU' where U is n-by-k and D is kxk. If q is a unit vector, what can you say about the rank of [q, Aq, A^2*,..,A^(n-1)q]?

--------------------

Nov 9  Let A = [lambda_1  m ; 0 lambda_2] and assume that the eigenvalues lambda_1 and lambda_2 are distinct. Determine X so inv(A)*A*X = diag(lambda_1,lambda_2). Recall that if Ax = lambda* x and y'*A = lambda*y' then the condition of lambda is  1/abs(x'*y) assuming that x and y have unit 2-norm. What is the condition of lambda_1 and lambda_2 in the above 2-by-2 example? More generally, how would you compute the condition number of all the eigenvalues using [X,D] = eig(A) assuming that the eigenvalues of A are distinct?

--------------------

Nov 11   A is a real matrix and mu is a complex eigenvalue. Suppose Az = mu*z. Show that the real and imaginary parts of the (complex) eigenvector z define a real, 2-dimensional invariant subspace. Suppose the dominant eigenvalue lambda_1 of a real matrix is complex. Note that the conjugate of lambda_1 is also an eigenvalue so the dominant eigenvalue is not unique. What happens if we apply the power method with a real starting vector? Assume that all the other eigenvalues are smaller in maginitude.

--------------------

Nov 13   Suppose x and y are unit 2-norm n-vectors. How would you choose mu so that ||A*x-m*x||^2 + ||A'*y - mu*y||^2  is minimized?

--------------------

Nov 16

--------------------

Nov 18

--------------------

Nov 20

------------------

Nov 23

--------------------

Nov 30

--------------------

Dec 2

--------------------

Dec 4

--------------------