Department
of
Mathematics

University of

Lincoln

 CSCE/Math 447

Welcome to the Math 447, Numerical Analysis II: Numerical Linear Algebra home page. You're probably here for information, so let's start with the vital statistics of the course.

### Essential Information

• Announcements Announcements about current class activities.
• Class Policy Statement
• Course Syllabus
• Assignments Assignments and due dates.
• Class Materials Some useful course materials for this course, including a printable copy of the Class Policy Statement and Syllabus are here in our home directory.
• Notes and FAQ Notes and answers to questions about class matters that students have asked. Check this periodically for the most current information. Also, check the class discussion board for information.
• Lecture Notes You'll find an outline of our activities and some extra notes here.

### Numerical Linear Algebra Resources

Lots of information can be found on the WWW. Go to your favorite search engine (like the Google or Yahoo sites listed on my home page) and try searching on "numerical linear algebra and optimization". See how many web pages you hit and visit a few interesting looking sites. As we cover specialized topics, try searching on them with keywords such as "Householder triangularization".

## Notes and FAQ

I've been asked several questions. Here they are, along with my answers:

On the revised question #5b, you've written "b=Ax, where x is a vector of ones. �nbsp;Solve the system Ax=b using the P L U from your program.." �nbsp;Do you really mean to give us the x and have us solve for b? �nbsp;Don't you rather mean let b be a vector of ones and then we solve for x?

Answer: No, I mean what I wrote. As we have done in other experiments, you are to generate the right hand side vector b with the assumption that the solution is actually a vector of ones. So when you solve for x, you expect to approximate a vector of ones. Just keeps the numbers nicer.

What is meant by a unit lower triangular matrix?

Answer: This means a lower triangular matrix whose diagonal entries are all ones...like the L that you get in an LU factorizaton.

08/08/06: Programming odds & ends...
Here's a question I've answered for other students in different classes that might be helpful to you. It identifies a problem with understanding of how MATLAB functions work:

I have a question about the Euler function for matlab. One of it's arguments is function name. What exactly would the function name be? We have an ODE we're trying to solve. So what is our function? Also if we define a function, we have to give it points to evaluate at.

Well, functions are variables, and like all variables they have names and content. In Matlab names are just strings of characters (strings for short). Now notice that the comments of the function Euler say that the function has to have a certain format, namely it has to have arguments t,y, in that order. If you are solving a vector system, then y must be a vector; otherwise it is a scalar. In our exercise, it is just a scalar. For example, if I wanted to create a function represents the rhs of the ode dy/dt = t*sin(y), I would create a "function file" by using my favorite editor (Matlab has one, just type "edit" while you are in the command mode and the editor may pop up with a prompt asking if you want to create a new file. If so, hit the OK button, and you're off and running. Another option is to use your own favorite plain text editor.) In any case, the next thing that I would do is save the empty file and give the future function file a name like "fcn.m". You must use the extension ".m". Here would be the contents of the file, the first line of which must start with "function...":

function retval = fcn(t,y)
% description: FYI, the percent sign is a comment symbol
% in both matlab and octave. and in fact whatever I put
% here will appear right up to the first blank line
% whenever I type "help fcn". Matlab and octave first see
% if the function is already in memory. They then search
% the current directory, then their own special function file
% directories for a name that matches "fcn.m". Then they open
% the file to see if its first line is "function x = fcn(xx,xxx)".
% the arguments of the function are xx and xxx, and the return value
% is whatever name follows "function". If so it is a function file
% and the programs will then print out all comment lines up to the
% first blank.

% This line won't appear in your help request, and in fact you don't
% need to type any comment lines at all. I like to remind myself of
% syntax by putting a brief help description in my functions.
% "retval" is just my favorite name for a returned value. You can use
% anything else. Also, the parameter list for the function, "xx,xxx"
% is important only inside the function. These are dummy variables or
% placeholders, just like the dummy variable of a definite integral.

retval = t*sin(y);

end % this last end is optional in matlab (or octave in a function file.)

Whew! That's all there is to it. Notice that the only mandatory lines in this file are the first line and the line "retval = t*sin(y)". Now if I invoke the function interactively in Matlab, say by typing "fcn(2,4*pi)", I'll get the return value of 2*sin(4*pi), which will be a specific number (in this case, 0.) BTW, Matlab has lots of builtin functions like sin, cos, ln, exp, etc, and builtin variable names like pi, eps, etc.
And finally, for the answer you've been waiting for, the name of the function I just defined is 'fcn'. NB: don't use double quotes in Matlab, only single. This is how you pass a function as a parameter in Matlab/Octave. Just say, for example, to solve the ode with IC y(0)=1 from 0 to 3 in steps of size 0.1, just type

Here's another question I got from a student visit to my office: "Did you have a copy of an explicit trapezoidal method somewhere?"
Well, actually, yes, but I forgot that it's somewhat hidden. Look at the Theta function that I defined. It implements the theta method, but it uses fixed point iteration to make the method explicit. I set the iteration number, "iterlim", to 2 which, if you pass the parameter theta=0.5, gives you the explicit trapezoidal method.

Finally, I'm going to offer a few tips to novice programmers that might make their life a little easier. Feel free to use or not use them or misuse them.

A few practical programming tips...

• Use the Matlab help files and the command line "help." Sometimes, they're pretty useful.
• Have a look at my tutorial files, whose URL is listed above under Numerical PDE Resources, and see if they help you.
• Create and use at least one "standard directory". For example, I have a subdirectory in my Math942 directory called ODEMethods. I'm keeping all the .m files for this course in that directory. When I run Octave from a terminal, I first change to that directory. You can run Matlab by icon and and use menus to set the Matlab path to this directory. In this way, when you issue a command like "fcn(2,3) and you have put the file fcn.m in this directory, Matlab will find the function.
• Recycle, recycle, recycle. When I'm going to create a new function, say "gcn.m" and it's somewhat like the function file "fcn.m", I load up the file fcn.m and immediately Save As to "gcn.m". I then edit the file to make the necessary changes that I need for gcn. You should do the same thing with function files you create or I've given you.
• Use script files to automate your work. A script file is a .m file that is not a function file, but has commands for Matlab to carry out. The first line of a script file should not start with "function ...", in contrast to a function file. Commands are then executed as though you typed them at the command line prompt. For instance, suppose I want to run Euler repeatedly using step sizes 1, 1/2, 1/4, 1/8, ..., 1/128, record the values of y(3) in an array, and then do a loglog plot of stepsize versus absolute value of y(3) (can't do loglog of negative numbers). I would write a script to do the work for me called, say, "EulerTest.m". The file would look like this:

% script: EulerTest.m
% description: does some calculations for an exercise.
% This script file *must* be run from a directory where the
% function files "Euler.m" and "fcn.m" live.
disp('Running calculations for Euler');
% start the experiment
numdivides = 7; % remove the semicolon if you want to see a result
h = 0.5.^(0:numdivides); % array of stepsizes
results = zeros(1,numdivides+1); % initialize an array to hold results
for kk = (1:numdivides+1)
msteps = round(3/h(kk)); % number of steps to get to 3
temp = Euler('fcn',0,1,h(kk),msteps);
results(kk) = temp(msteps+1);
end;
loglog(h,results); % the loglog plot

Run the script by just saying its name to Matlab "EulerTest" and pressing Enter. When it's done, you'll have a few new variables in your workspace, one of which is "results". (Type "who" and see its name.) Now do with it what you will.
BTW, Octave fans, here's an advantage of Octave over Matlab: you can include function definitions in Octave files, just like you can do in C or Fortan, etc. Then you could have a script portion of the file at the bottom that runs like a "main". Unfortunately, Matlab requires separate files for each function except for subfunctions (private functions) that get used by the main function sharing the name of the file and are invisible to everything else.
• About debugging something you've written...well, this is a big subject. Just for the record, Matlab has a pretty fancy debugger built in. For a quick and dirty solution (sometimes!) to what ails your program, try removing semicolons to actually see intermediate output, and also try sprinkling temporary extra display comments in the code you're debugging at key points, so you can see where your code is going. Use "disp()"
• One of the concepts you should make sure you understand is the idea of a "local" variable versus a "global" variable. Global variables are just that: visible to everyone and everything. Functions are global variables. So, if in my function definition of a function, say gcn, I make a reference to an already created function "fcn(t,x)" then this reference will work correctly. However, all other user variables are local by default (it is possible to declare them as global, but it takes an extra step.) This means that they are recognized only in the environment in which they were created. If I create a variable "xtmp" in a script file, then reference "xtmp" inside a function I create, then I really have two separate variables because the function environment is different from the user environment, which is where the script variables live -- remember, running a script is like typing in the same commands at the keyboard.

I've been asked to explain what significant digits of an approximation to a number mean. Eventually, we'll need to know, so here goes: to get the number of significant digits, first *subtract* (rather than just looking at the numbers) the two (may as well be larger - smaller,) then find the position of the leading digit of the error relative to the position of leading digit of the exact answer. (We're thinking in fixed point representation in this discussion.) If the difference in that position is less than 5, then number of significant digits is one less than that position, else two less.

For example if 1.006 and .996 are used to approximate 1, calculate 1.006 - 1 = 0.006. Notice I put the zero in front of the decimal to start counting from the right position. There is nonzero digit at the 4th position, counting from the (base 10) position of the leading digit of 1, and the size of this digit is greater than 5, so this approximation has 2 significant digits. On the other hand, 1 - .995 = 0.005, which again has a nonzero digit in the 4th position, and the size of the digit is at most 5, so .995 has 3 significant digits as an approximation to 1. BTW, it's also perfectly correct to say that each answer has one significant digit, though this doesn't give all the available information. Hope this helps.

## Course Description

So what is numerical linear algebra? Math 314 II? An analysis of rounding errors? Not exactly. The textbook for this course is a gem, and its preface contains one of the most elegant descriptions that I know:

"The field of numerical linear algebra is more beautiful, and more fundamental, than its rather dull name may suggest. More beautiful, because it is full of powerful ideas that are quite unlike those normally emphasized in a linear algebra course in a mathematics department. (At the end of the semester, students invariably comment that there is more to this subject than they ever imagined.) More fundamental, because, thanks to a trick of history, 'numerical' linear algebra is really applied linear algebra. It is here that one finds the essential ideas that every mathematical scientist needs to work effectively with vectors and matrices. In fact, our subject is more than just vectors and matrices, for virtually everything we do carries over to functions and operators. Numerical linear algebra is really functional analysis, but with the emphasis always on practical algorithmic ideas..."

Here is a list of topics we will cover: fundamentals (review of some linear algebra with a smattering of new ideas), QR Factorization and least squares, Gram-Schmidt orthogonalization, conditioning and stability, systems of equations, eigenvalues and iterative methods.

There will be a midterm, a final and homework assignments. Some Matlab programming will be required, but all of it can be done in the Mathematics Computer Lab and no prior experience with Matlab is required. Most of the programming can also be done in Octave, an open source software project that runs under most flavors of linux and Windows.

## Class Policy Statement

Place/Time: AvH 12, 3:30-4:45 MWF, Fall 2006

Preq: CSE/Math 340, Math 221 and 314 or equivalent, or permission.

Objectives: To help students achieve competence in the following areas:

• Understanding of basic theory of numerical linear algebra.
• Knowledge of the most common algorithms and their applicability.
• Stability and efficiency (complexity) of these algorithms.
• Implementation of and experiments with the algorithms (mainly via Matlab platform).

Instructor: Dr. Thomas Shores

Telephone: Office 472-7233 Home 489-0560

Email: tshores1@math.unl.edu

Office Hours: Monday 10:00-11:30, Tuesday 11:00-13:00, Thursday 12:00-13:30, Friday 8:30-10:30, and by appointment. Office: 229 AvH

Class Attendance: Is required. If absent, it is incumbent upon the student to determine what has been missed as soon as possible. It is advisable to consult with the instructor.

Homework/Projects: Homework will be assigned in class and collected according to the syllabus collection dates. All homework assigned one week or more before the collection date is due on that date. Most of the homework assignment problems will be graded in detail for homework points. It is strictly forbidden to copy someone else's homework, though some discussion with other students is allowed. For some of the more substantial exercises I will allow teams of two to turn in a single solution with clear attribution of partnership. The official programming language for this course is Matlab. Prior experience in Matlab is not required. Students will be given an account in the Mathematics Computer Lab for computer related exercises and can obtain written lab instructions in the lab itself. Current information about the course will be available on the web (via the 447 homepage or my home page). Using the web is strongly recommended for keeping track of due dates for homework collections and other current activities in the course.

Reading Assignment: Read the sections of the texts as, or before, they are covered in class lectures. This is a standing assignment throughout the semester.

Grade: A midterm will be given outside of class meeting time and will account for 150 points. The final exam will count 200 points. Exams may have a take home component. All in-class exams are closed book. Homework will count 150 points. The final grade will be based on these 500 points.

Final Exam: Will be comprehensive. To be given on Tuesday, December 12, 8:30 - 10:30 pm in AvH 12.

Grades of "I", "W" or "P": These grades will be given in strict accordance with University policy. (See any Schedule of Classes for the relevant information and dates.)

Keep This Information!!!

## Syllabus for CSE/Math 4/847, Spring 2006

• TEXT: Numerical Linear Algebra, Lloyd Trefethen and David Bau, SIAM, Philadelphia, 1997.

The times listed below are approximate, and may be adjusted as the semester progresses. Text sections are identified by lecture number. All assignments are due on Tuesday of the weeks specified below.

 WEEK DATES LECTURE TOPICS

 1 Aug 21-25 1 Matrix-Vector Multiplication 2 Orthogonal Vectors and Matrices Matlab Orientation

 2 Aug 28-Sep 1 3 Norms 4 Singular Value Decomposition 5 More on the SVD

Friday, September 1, is the last day to withdraw from the course and not have it appear on your transcript.

 3 Sep 4 Labor Day Sep 5-8 6 Projectors 7 QR Factorization

 4 Sep 11-15 Assignment 1 Due 8 Gram-Schmidt Orthogonalization 9 Mathematical Software 10 Householder Triangularization

 5 Sep 18-22 11 Least Squares Problems 12 Conditioning and Condition Numbers 13 Floating Point Arithmetic

 6 Sep 25-29 14 Stability 15 More on Stability 16 Stability of Householder Triangularization

 7 Oct 2-6 Assignment 2 Due 17 Stability of Back Substitution 18 Conditioning of Least Squares Problems

 8 Oct 9-13 19 Stability of Least Squares Algorithms 20 Gaussian Elimination 21 Pivoting

Friday, October 13, is the last day to change your grade option to or from ``Pass/No Pass''.

 9 Oct 16-17 (no class) Fall Break Oct 18-20 22 Stability of Gaussian Elimination 23 Cholesky Factorizations

 10 Oct 23-27 Assignment 3 Due REVIEW MIDTERM

 11 Oct 30-Nov 3 24 Eigenvalue Problems 25 Overview of Eigenvalue Algorithms 26 Reduction to Hessenberg or Tridiagonal Form

 12 Nov 6-10 27 Rayleigh Quotient and Inverse Iteration 28 QR Algorithm without Shifts 29 QR Algorithm with Shifts

Friday, November 10, is the last day to withdraw from the course and receive a grade of W.

 13 Nov 13-17 Assignment 4 Due 31 Computing the SVD 32 Overview of Iterative Methods

 14 Nov 20-21 33 Arnoldi Iteration Nov 22-24 Thanksgiving Vacation

 15 Nov 27-Dec 1 34 How Arnoldi Locates Eigenvalues 35 GMRES 38 Conjugate Gradients

 16 Dec 4-8 Assignment 5 Due 40 Preconditioning REVIEW

Final Exam: The final exam is a comprehensive exam to be given on Tuesday, December 12, 8:30 - 10:30 pm in AvH 12.

Department Grading Appeals Policy: The Department of Mathematics and Statistics does not tolerate discrimination or harassment on the basis of race, gender, religion or sexual orientation. If you believe you have been subject to such discrimination or harassment, in this or any math course, please contact the department. If, for this or any other reason, you believe that your grade was assigned incorrectly or capriciously, appeals may be made to (in order) the instructor, the department chair, the departmental grading appeals committee, the college grading appeals committee and the university grading appeals committee.

## Assignments

1. (Points: 30, Due: September 14, finished)
Textbook numbers 1.1, 1.3, 2.1, 2.5, 2.6, 3.2, 3.3, 4.1, 4.4, 5.3.

2. (Points: 30, Due: October 3, finished)
Textbook numbers 6.4, 6.5, 7.1, 7.5, 9.2(only math grads justify (c)), 10.2 (also (c): repeat Experiment 2 of Lecture 9 as done in class with the Householder QR included and compare results.), 11.3, 12.1.

3. (Points: 30, Due October 24, finished)
Textbook numbers 12.3ab, 13.3, 14.1, 15.1abcd, 15.2, 16.1a, 17.1, 18.1.

4. (Points: 30, Due November 16, finished)
Textbook numbers 20.1, 20.3, 21.4, 22.1, 22.2, 23.3 (do the experiment with m=1000), 24.2, 25.1.

5. (Points: 30, Due December 5, finished)
Textbook numbers 26.3, 27.5(ugrads need only exhibit a nontrivial example that illustrates the claim of the exercise), 28.2, 29.1ab, 33.2, 35.4, N1.