Problem of the Day
May 3 (a) A(lambda) is an nxn matrix that depends on a real number lambda. Assume that the ij entry of this matrix is differentiable. How could fSolve be used to find lambda_star and a unit vector x_star so that
A(lambda_star)x_star = 0.
(b) N points on the surface of the earth have elevations h(1),...,h(N). A surveyor makes the rounds and assembles estimates of the elevation differentials in a matrix G. In particular, G(i,j) is an estimate of h(i)-h(j), obtained by standing at point j and measuring how much one must climb to reach point i. (Could be negative) We assume that G is very sparse. (When standing at point j you can only see a handful of other points.) Assuming that h(1) is known, how could LSQNONLIN be used to estimate h(2),...,h(N)? Solution
May 1 Suppose S is nxn, symmetric, and rank(S) = r. Explain why I + S has at most r+1 distinct eigenvalues? Solution
April 29 The method of steepest descent with exact line search is applied to phi(x) = .5*||Ax - b||^2 where it is assumed that A has full column rank. What is the gradient g_c of phi at x = x_c and give a recipe for computing x_new= x_c + mu_c*x_c where u_c is the exact line search step parameter. Solution
Apr 26 Look at ShowSteepDes and derive the recipe for the step-length parameter mu. Solution
Apr 24 Suppose A is mxn with m>n. Suppose rank(A) = r < n and that the SVD A = USV' has been computed. Compute a matrix Q with orthonormal columns that span the nullspace of B = [0 A; A' 0]. Solution
Apr 22 Look at the conic section fit script Apr22. What happens if the points (x(1),y(2)),...,(x(n),yn)) are colinear? Solution
Apr 19 To minimize
Q(s) = fc + gc'*s + (s'*Hc*s)/2
subject to the constraint that norm(s,2) = alfa ,we must choose mu so that the solution to
(Hc + mu*I)s = -gc
satisfies norm(s,2) = alfa. How could this be done? Here Hc is the current Hessian, and gc is the current gradient. Solution: Same as Apr 17!
Apr 17 Suppose H is positive definite and Q'*A*Q = D its schur decomposition. Let s(mu) be the solution to the linear system (H + mu*I)s = -g. How would you compute mu so that norm(s(mu) = delta where delta >0. Assume that norm(s(0)) > delta. Solution
Apr 15. Suppose t and y are column m-vectors. We wish to fit data points (t(i),y(i)), i=1:m with a function
f(t) = a*exp(bt).
Taking the least squares approach, the idea is to choose a and b so that norm(a*exp(bt) - y,2) is minimized. Show how this problem could be solved using the 1-D minimizer fminbnd. Hint: for a fixed b, what is the best a? Suggest a heuristic for picking the initial bracketing interval. Solution
Apr 12. Suppose f(x) = (1/2)x'Hx - x'*b where H is nxn and symmetric and b is nx1. Under what conditions does f have a finite minimum value? Solution
Apr 5 Suppose p(x)= a(1)+a(2)x + a(3)x^2 + a(4)x^3 + x^4 has 4 real roots a<b<c<d. We apply golden section search with initial interval [L,R] . Under what conditions will the process find a local minimum? Solution
Apr 3 How would you find the closet point on the ellipse (x/a)^2+(y/b)^2 = 1 to the point (u,v)? I.e., set up an objective function of a single variable that you want to minimize. Solution: This was an assignment
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Apr 1 Suppose F(x) is a function from Rn to Rn with the property that its Jacobian J(x) is tridiagonal for all x. Explain how one can generate a divided difference approximation of J with just four F-evaluations. Hint: for a cleverly chosen vector v, (F(x+hv)-F(x))/h will provide an approximation to J(:,1), J(:,4), J(:,7), etc. Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Mar 29 Consider applying the secant method to the function f(x) = (x^2 - A) where A>0. Show that
e(k+1) = e(k) * e(k -1) / (x(k) + x(k-1) ) where e(j) = x(j) - sqrt(A)
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Mar 27 Give a recipe for the line p(x) that interpolates F(x) = sqrt(x) at x = 1/4 and x = 1. How large can |F(x)-p(x)| be on the interval [ .25 , 1] ? Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Mar 25 Suppose B has eigenvalues lambda(1) > lambda(2) > ... > lambda(n). Assume that mu is a given scalar. What happens if we apply the power method with A = inv(B - mu*I)? Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Mar 15 Use the Schur decomposition to show that if A is symmetric with eigenvalues
lambda(1) > lambda(2) > ... > lambda(n)
then there exists a symmetric matrix B with repeated eigenvalues that satisfies norm(A-B) = del/2 where del is the gap between the two closest eigenvalues of A Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Mar 13 Suppose T is symmetric tridiagonal with diagonal entries a(1),..,a(n) and subdiagonal entries b(1),...,b(n-1).Let p_{r}(x) be det(A(1:r,1:r)-xI) for r=1:n. Let p_{0}(x) = 1 and note p_{1}(x) = a(1)-x. It can be shown that
p_{r}(x) = (a(r)-x)p_{r-1}(x) - b(r-1)^2p_{r-2}(x) r=2:n
SHow how to evaluate the derivative (d/dx)p_{n}(x). Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Mar 11 Using the Matlab functions hess and qr, show how to compute an orthogonal Q so that Q'(CC')Q = T is tridiagonal where C is mxn and m>>n. Hint. Do not form CC'. Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Mar 8 Suppose A is m-by-n and rank deficient. How would you compute the minimizer of norm(Ax-b) that is closest to a given vector d? Hint, Use the SVD and the fact that if xStar is a minimizer then so is any vector of the form xStar + z where z is a linear combo of singular vectors V(:,r+1),...,V(n). Solution.
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Mar 6 Suppose A is mxn and rank(A) = m < n. The system Ax = b where b is mx1 is underdetermined and has an infinite number of solutions. Show how to get a solution (any solution) by using [Q,R] = qr(A'). Solution.
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Mar 4 Recall that norm(M,'fro') returns the Frobenius norm, i.e., sqrt(sum(sum(A.*A))). By making effective use of the thin QR factorization show how to compute an nxp matrix X so that norm(AX-B,'fro') is minimized were A is mxn and B is mxp. Solution.
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Mar 1 Suppose u is a unit m-vector and u(1:k-1) = 0. Suppose A is an m-by-n matrix. Carefully show how you would overwrite A with H*A where H = I - 2*u*u'. Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Feb 27 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) Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Feb 18 Under what conditions will Gauss-Seidel converge when A = [a b ; c d]? Solution.
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Feb 15 If A is n-by-n, symmetric, and positive definite, explain why A(1:n-1,1:n-1) is symmetric and positive definite. Suppose b is a given column n-vector. How would you go about computing x and y so that Ax = b and A(1:n-1,1:n-1)y = b(1:n-1)? Solution.
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Feb 13 [L,D,piv] = ldl(A,'vector') computes P'*A*P = LDL' with P = I(piv,:). How would you solve Ax = b given L,D, and piv? Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Feb 11 (a) Vectorize the inner loop in MyChol. (b) A is symmetric positive definite and we compute A = LU. How would you compute lower triangular G so A = GG'? Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Feb 8 Modify the function luRecur so that it reports back the fraction of flops that are level-3 flops. Solution.
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Feb 6. This is about an output option associated with the Matlab LU function. [L,U,piv] = LU(A,'vector') does the same thing as [L,U,P] = lu(A) except that P is represented as an n-vector of integers. In particular, if I = eye(n,n). then P = I(piv,:). (a) How would you solve Ax = b given [L,U,piv] = lu(A,'vector'). (b) How would you solve A'x = b given [L,U,piv] = lu(A,'vector')? Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Feb 4. Assume that A is a given nxn matrix and that the comand [L,U,P] = lu(A) has been executed. (a) How would you solve Cx = b given that C = A*A' ? (b) Given the nxn matrix F and nx1 vectors f and g, how would you find nx1 vectors u and v so that [A F ; zeros(n,n) A] *[ u ; v ] = [f ; g]. Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Feb 1. What is the largest power of 10 that can be represented exactly in the IEEE format:
( +or-) (1 . b1 b2 ... b52)_{2} x 2^ ( (a1 a2 a11)_{2} - 1023)
Hint. 10 = 2 x 5 Solution
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Jan 30. Assume that the svd of an nxn A has been computed, i.e., U, S, and V are are available where [U,S,V] = svd(A). (a) How would you compute B = U*T*V' where T is the same as S except in the (n,n) position where we set T(n,n) = 0. (b) How would you solve A'x = b where b is a given nx1 vector? (c) How would you compute an orthogonal matrix Q and a symmetric matrix S so that A = QS? Solution.
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Jan 28. A 2x2 orthogonal matrix is either a rotation Q = [c s;-s c] or a reflection Q = [c s; s -c]. Here c^2+s^2=1. Note that a rotation has determinant 1 while a reflection has determinant -1. Complete the following Matlab function so that it performs as specified:
| function [U,S,V] = svd2(A) |
| % A is 2x2 and U*S*V' is its svd |
| % If det(A)>=0, then U and V are both rotations |
| % If det(A)<0, then U is a reflection and V is a rotation. |
Make use of the Matlab SVD function. Solution.
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Jan 25. Answer these questions by benchmarking with the Matlab tic and toc functions. (a) True or false, for large n and random nxn A, norm(A) takes about n times as long to execute as norm(A,1). (b) True or false, for large n, the product of an upper triangular matrix and a vector takes about half as long as the product of a full matrix and a vector. (c) True or false, for large n, x = A\b takes much longer to execute if A is full rather than triangular. How much longer as a function of n? Solution: Jan25.
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Jan 23. This problem gives you experience with Matlab's colon notation and Boolean capabilities. 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 distinct positive digits. 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. Solution: Jan23
+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
Jan 21. (a) Write a script that generates a random 2-by-2 matrix A using randn and then graphically displays all points of the form (y1,y2) where [y1 ; y2] = A*[cos(theta) ; sin(theta)] where 0 <=theta<=2*pi. (b) Write a script that generates a random 2-by-2 matrix A using randn and then graphically displays all points of the form (y1,y2) where [y1 ; y2] = A*[x1;x2] where |x1| + |x2| = 1. Solutions ideas Jan21, Jan21b