Archives

Nov 30. Clarifciation for the case c = .05 in Assignment 6. Either find d_opt(.05) or show that it doesn't exist. "Show" can be by pencil and paper and/or by experimental evidence.

Nov 28 Take-Home Prelim test scripts and solutions. (P1, P2, P3, P4, QuadModelCV, SkewSchur3CV, expofitCV,EnclosingCV.

Nov 28. A6 is available.

Nov 17. The A5 solutions are available.

Nov 17. The Take-Home Prelim is due Tuesday, Nov 21, 11:59PM

Nov 17. The course-evaluation system will be open for students to complete evaluations from Monday, 20 November to Sunday 3 December (inclusive). You can see what the course evaluation system looks like at this URL:
http://www.engineering.cornell.edu/courseeval/

Nov 16. P2Hint. Find a unit-2norm null vector v for A, i.e., Av = 0. You then have a Householder problem: find unit 2-norm u so (I-2uu')*e1 = v . Note that the given matrix is skew-symmetric: A' = -A, i.e., A(j,i) = -A(i,j). Skew-symmetric matrices have zeros on the diagonal. Also, Q'*A*Q is skew-symmetric if A is skew-symmetric. This means that if C = Q'*A*Q has a zero first column (and row) then C(2:3,2:3) will automatically have the form [0 tau;-tau 0]for some tau since B is skew-symmetric.

Recall how we construct Householder matrices and what they do. In particular P = I - 2uu' is orthogonal and symmetric.
so Px = y implies Py = P*P*x = P'*P*x = x.

Nov 16. P1Hint. A statement of the form "x is a linear combination of certain vectors" means x = My where M is a matrix whose columns are the vectors and where y is a column vector.

Nov 13. Sample test scripts and templates. P1, QuadModel, P2, SkewSchur3, P3, ExpoFit, P4, EncloseE. Note: your implementations will be tested against scripts that may be more "hostile" than P1, P2, P3, and P4.

Nov 13. The take home prelim is available. Email for a pdf.

Nov 6 Today's 2:30-3:30 office hours are shifted to Tuesday 11-12.

Nov 2. Extra credit for P3. As defined in P3Hint, the objective function B(theta) = f(theta)/d(theta)^2 uses somewhat faulty "fraction" function f(theta). In particular, the line thru points E and C should be tangent to the "Venus Circle" in the diagram. Likewise, the line thru points E and B should also be tangent line. Rework the problem with B and C chosen in this way. The resulting value of f(theta) will be smaller. To do this you need to know the diameter of the inner planet. Use .007564 for Venus and .007973 for Earth.

Nov 2. In P4, your implementation of PolyInEllipse will have two subfunctions, one that computes the negative of the area and one that computes the negative of the perimeter. They will have the form MyF(t,b) where t is an n-vector and b is the semiaxis parameter.

You will have to modify the default parameters for the optimizer fminsearch. This is generally a slow method. However, it is attractive because it requires no gradient or Hessian encoding. You'll have to include something like this in your implementation of PolyInEllipse:

    N = 20000;
    options = optimset('MaxFunEvals',N,'MaxIter',N);
    topt = fminsearch(@(t) MyF(t,b),tInitial,options);

A special feature of our problem is that the t-vectors are ordered; 0<=t(1)<=t(2)<=...<=t(n)<=2*pi. fminsearch will not automatically enforce this property as it iterates. So it is critical that you make this the first statement in your perimeter and area objective functions: t = sort(t)

Nov 2. In P2, Broyden's method (described on p. 241) "updates" an approximate Jacobian each step. The question arises, what is the initial approximate Jacobian? (B_0 in the book notation). You should set B_0 to be the exact Jacobian evaluated at xInitial.

A cool feature of Broyden's method is that the linear system solving overhead can be reduced to O(n^2) by clever updating of the QR factorizaion. As shown in class, if you have QR, = A then you can get Qnew*Rnew = A +rank-1 in O(n^2). Anyway, you are not expected to include this feature in your implementation of TwoWay.

Oct 31. In P1, note that if L and R are adjacent roots of f ' (x) and f(L)*f(R)< 0 , then f has a root

in [L,R]. You also have to check if [-1,xmin] and [xmax,1] are bracketing where

xmin and xmax are the smallest and largest zeros of f ' in [-1,1]. See the hint.

Oct 31 The only trig you need for P3 is summed up in tau = acos(u'*v/(norm(u)*norm(v)) which is the angle between the vectors u and v. See the hint.

Oct 23. To be eligible for the extra credit problems you must email me (pdf's please') the corresponding "regular" problem by Tuesday, 11:59PM.

Oct 23. Clarifications for Problem 2. Here is a theorem which you can use to figure out the tolerance that you hand over to JacSVD:

Theorem. If sigma(n) is the smallest singular value of A and sigma_tilde(n) is the smallest singular value of A_tilde, then |sigma(n) - sigma_tilde(n) | <= norm(A - Atilde,'fro').

Another point. [U,S,V] = svd(A) always returns a real S even if A is complex. Fact: Because JacSVD uses Matlab's SVD to solve the 2-by-2 problems, A(i,i) stays real once we do a subproblem that involves index i. However, if the incoming A is "too diagonal" and JacSVD doesn't do anything, you'll exit with a complex S. What are you going to do about that?

Oct 23. There was a missing absolute values in the JacSVD subfunction function off. It has been corrected. Since A is complex you want to sum all the | A(i,j) |^2 not all the A(i,j)^2 .

Oct 18. Some extra credit problems associated with the current assignment are available.

Oct 16 Take a look at this subset of the complex plane. Each x+iy is colored (red, green, blue,orrange) according to whether the Newton iteration zNext = zCurrent - f(zCurrent)/f'(zCurrent) converges to 1, i, -1, or -i where f(z) = z^4 - 1. You can see that if you draw a small enough circle around one of these roots, then the iteration converges to that root. Part of P3 is to figure out (approximately) the largest radius of these neighborhoods.

In problem 3 it is fine to use the matlab det function for computing the determinant of a 3-by-3 orthogonal matrix. Usually, det has limited value, but in our setting it is fine.

As we presented the min norm(A - BQ,'fro') problem in class, we spoke of the need for Q to be a rotation, i.e., to have determinant equal to +1. In problem 3, we do NOT add this constraint to the problem. Just go about finding the best orthogonal Q regardless of its determinant.

Note: feel free to use norm( <matrix> ,'fro') when assessing how close you can align via rotation and translation

Oct 2 Another typo in Prob 3. In (b) should be "For any orthogonal Q, show that v_opt = (AQ' - B)'e/m. The " Q' " is missing. You do NOT need any special 'trace" properties of the Frobenius norm to prove the optimal selections for v and Q.

Sept 29. Typo in Problem 3 of Assignment 3. The last line under item "(b)" should end with ... " = USV^{T} is the SVD." The middle matrix (the diagonal matrix of singular values) is missing.

Sept 29. Some relevant textbook Review questions for Prelim 1 (scheduled for Friday, October 6, during class): 1.1-1.10, 1.32, 1.42, 1.48, 1.49, 2.1-2.5, 2.7, 2.9, 2.13-2.17, 2.21-2.28, 2.35, 2.42,2.44, 2.46, 2.49, 2.53, 2.57, 2.61, 2.65, 2.76, 3.6, 3.9, 3.10-3..20 , 3.29, 3.38, 3.41, 3.42, 3.43. Study the 'Problems of the Day" too!

Suppose f and g are column n-vectors and that you know Alfa = f'*g. Assume that P is an m-by-m orthogonal matrix and that u = Pf and v = Pg. Show how to compute NewAlfa = u(2:m)'*v(2:m) in O(1) flops.

Sept 18. The benchmarking script for P1 in Assignment 2 has the correct placement of tic and toc. For the Cholesky method it makes sense when benchmarking to include the set-up-A portion of the calculation.

Sept 18 Solutions to Assignment 1 are available. Answers to the questions appear in the comments of the solution codes. The mean scores on the 5 problems were [ 3.1 3.5 4.4 4.3 3.3].

Sept 14 There were some typos. In P2, replace all instances of M1/M2 with M1*M2. Reason: Since M1 approximates norm(D+CC',1) and M2 approximates norm(inv(D+CC'),1) we see that M1*M2 approximates cond(D+CC',1). Incidentally, in P2A you are free to use Matlab's cond function to get the "exact 1-norm condition.

In P3 there are 10 data points so the "12" should be changed to "10" in the paragraph that talks about the introduction of noise.

Regarding the plotting of your computed ellipses in P3, I have written a demo script ShowContour that shows how you can use Matlab's contour function to accomplish this task easily. Submit only one plot. It should depict the graph of the "A\b ellipse" based on the original data, the 'A\b ellipse" based on the noisey data, AND the 10 data points. This will nicely communicate one of the main points of the problem: how the fitting ellipse changes when the data changes

Sept. 8. Assignmennt 1 comments. P1. 10^k = 5^k*2^k. The mantissa must hold 5^k exactly. P2. Take a look at the partial sums. P3. Typing the command "format long" tells Matlab to display results to full precision.. Review Matlab's pointwise vector operations. P4 Review the effect of condition number on Ax = b accuracy. P5 It is fine to write a script that sheds light on the two questions.