# Course 18.086: Computational Science and Engineering II (Spring 2008)

Department of Mathematics
Massachusetts Institute of Technology

## Policy

The grading will be based on homework exercises and a project. There will be no exams.
• Homework exercises: Most homework problems will involve programs. The focus will be on Matlab. It is illegal to consult solutions from previous years.
• Computational project: Over the whole lecture time, each participant works on an extended problem related to the content of the course. Project may relate to your thesis work, but the project work must be unique to this course. It is illegal to recycle old projects. Please consult the course information for more details.
• Till 02/19/2008 (second week): Submit project proposal
• Till 03/21/2008 (before spring break): Submit midterm report on project
• On 05/10/2008: Project presentations (presentation program)
• Till 05/12/2008: Submit final report on project

## Outline

This course addresses students in engineering and related fields with interest in computational science. Presented will be fundamental principles and equations, thus equipping the audience with the ability to address specific problems in their future career.
The course consists of three main topics: initial value problems, solving large systems, and optimization. The goal of the course is to provide a good start into each of these fields, focusing more on fundamental ideas than on involved details. Focus will be given on the mathematical understanding as well as on applying the presented concepts. Homework will include programming exercises, using Matlab. Outline of the topics:
1. Initial value problems
Linear initial value problems such as the wave equation and the heat equation admit closed form solutions in simple geometries. In a more complex setup they have to be solved numerically. Many important aspects seen for the linear problems regarding stability and covergence transfer to nonlinear initial value problems. Important examples for these are flow problems (fluid flow, traffic flow, etc.). Solutions to nonlinear problems may become non-smooth, and numerical methods have to consider this fact.
Topics: Stiff problems, wave equation, heat equation, Airy equation, convection-diffusion, conservations laws, front propagation, Navier-Stokes equations, Fourier methods, finite differences, consistency, stability, covergence order, Lax equivalence theorem, CFL-condition, leapfrog method, staggered grids, shocks, upwind, Lax-Wendroff, finite volume methods, level set method
2. Solving large systems
The discretization of partial differential equations by finite difference or finite element methods leads to large sparse linear systems, either directly for linear problems or as an auxiliary subproblem for many nonlinear problems. Gaussian elimination destroys the sparse structure, so solvers are required which make use of the specific sparse matrix structure.
Topics: Applications yielding sparse matrices, elimination with reordering, iterative methods, preconditioning, incomplete LU, multigrid, Krylov subspaces, conjugate gradient method
3. Optimization and minimum principles
Optimization problems search for the minimizer of some quantity (cost function), possibly given constraints. Quadratic cost functions lead to linear systems using Lagrange multipliers and Kuhn-Tucker conditions. Saddle point problems, regularization and calculus of variations will be presented as fundamental concepts. A different world in encountered in the case of linear cost functions. Applications are operations research and network problems. Solution algorithms are the simplex method or interior point methods. The underlying principle in all approaches is the concept of duality.
Topics: Weighted least squares, duality, constrained minimization, inverse problems, calculus of variations, saddle point problems, linear programming, simplex method, interior point methods, adjoint methods
What you get when you take this course:
• An introduction into and overview over numerical methods for initial values problems, large systems and optimization
• Focus on finite difference methods, efficient Krylov methods and multigrid, and various optimization techniques
• Solid knowledge about fundamental mathematical concepts and principles in computational science, such as stability, convergence study, model equations, numerical diffusion, staggered grids, preconditioning, Krylov subspaces, duality, regularization
• Knowledge about modern problems and techniques, such as CFD, level set methods, particle methods, multigrid, adjoint methods
• Hints for and practice in good techniques to do scientific programming in matlab.
• Useful matlab codes, to a large extend made by yourselves. You will get your fingers dirty in this course.
• The opportunity to work on a computational project related to your own research
What this course does not offer:
• Finite element methods, Galerkin methods, special matrices, heuristics for global optimization (e.g. simulated annealing)
• Sophisticated aspects (e.g. ENO/WENO schemes for advection) can at most be touched
• Involved coding in C++ or fortran, parallelization
• Simple and straightworward homework problems
• Exams

## Schedule

HomeworkDayDateLectures
1 Wed02/06 1Introduction & motivation
Fri02/08 2Fourier methods for linear initial value problems6.1
Mon02/11 3Finite difference methods & ordinary differential equations6.2
Wed02/13 4Von Neumann stability analysis6.1 & 6.2
Fri02/15 5Finite difference methods for the one-way wave equation6.3
project proposalTue02/19 6Modified equationLeVeque
Wed02/20 7Lax equivalence theorem6.3
21Fri02/22 8Numerical error analysis & convection-diffusion problems6.5
Mon02/25 9Wave equation, leapfrog, staggered grids6.4
Wed02/2710Conservation laws: examples, method of characteristics6.6
Fri02/2911Conservation laws: weak solutions6.6
Mon03/0312Conservation laws: high resolution shock capturing6.6
Wed03/0513Conservation laws: finite volume & particle methods6.6
32Fri03/0714Higher space dimensions
Mon03/1015Systems of initial value problems
Wed03/1216Level set method6.8
Fri03/1417Navier-Stokes equation6.7
Mon03/1718Navier-Stokes equation6.7
Wed03/1919Generating sparse matrices
3, midterm reportFri03/2120Elimination with reordering7.1
03/24-03/28Spring break
4 Mon03/3121ILU, preconditioning, iterative methods7.1 & 7.2
Wed04/0222Iterative methods7.2
Fri04/0423Multigrid methods7.3
Mon04/0724Multigrid methods7.3
Wed04/0925Krylov subspace methods7.4
Fri04/1126Krylov subspace methods7.4
Mon04/1427Least squares problems8.1
Wed04/1628Weighted least squares & duality8.1
54Fri04/1829Minimization with constraints8.1
04/21-04/22Patriot's vacation
Wed04/2330Regularization & pseudoinverse8.2
Fri04/2531Inverse problems8.2
Mon04/2832Tychonov regularization & ill posed operators in function spaces8.2
Wed04/3033Calculus of variations8.3
5Fri05/0234Calculus of variations8.3
Wed05/0736Linear programming8.6
Fri05/0937Linear programming8.6
Wed05/1439Project presentations

## Homework Exercises

• Theoretical parts: On paper. Submit in the lecture or in office hours.
• Programming parts: Matlab m-files. Submit by email, in the following way:
• Name your matlab file: ex{#}_part{#}_{yourname}[_{more}].m.
Examples: My program, solving u_t=-u_xxxx would be named ex1_part2_seibold.m.
• Email your files as attachments to both the lecturer (seibold) and the grader. Use subject: 18.086 homework {#}
• Write clean code. Comment your programs appropriately. Programs which do not run on our matlab are not graded.

## Matlab Programs

Please consult the Computational Science and Engineering web page for matlab programs and background material for the course.
• Finite differences for the heat equation: mit18086_fd_heateqn.m (CSE)
Example uses homogeneous Dirichlet b.c. on the left, and homogeneous Neumann b.c. on the right, and explicit Euler in time, which can easily be changed to implicit Euler.
• Four linear PDE solved by Fourier series: mit18086_linpde_fourier.m
Shows the solution to the IVPs u_t=u_x, u_t=u_xx, u_t=u_xxx, and u_t=u_xxxx, with periodic b.c., computed using Fourier series. The initial condition is given by its Fourier coefficients. In the example a box function is approximated.
• Finite differences for the one-way wave equation, additionally plots von Neumann growth factor: mit18086_fd_transport_growth.m (CSE)
Approximates solution to u_t=u_x, which is a pulse travelling to the left. The methods of choice are upwind, downwind, centered, Lax-Friedrichs, Lax-Wendroff, and Crank-Nicolson. For each method, the corresponding growth factor for von Neumann stability analysis is shown.
• Compute stencil approximating a derivative given a set of points and plot von Neumann growth factor: mit18086_stencil_stability.m
Computes the stencil weights which approximate the n-th derivative for a given set of points. Also plots the von Neumann growth factor of an explicit time step method (with Courant number r), solving the initial value problem u_t = u_nx.
Example for third derivative of four points to the left:
>> mit18086_stencil_stability(-3:0,3,.1)
• Numerical error analysis for upwind, Lax-Friedrichs, and Lax-Wendroff for transport equation: mit18086_error_analysis.m
Applies the upwind, Lax-Friedrichs, and Lax-Wendroff method to the transport equation u_t=u_x for a sequence of step sizes for final time and Courant number fixed. The norm of the errors is plotted and compared to the theoretical expected results.
• Finite differences for the wave equation: mit18086_fd_waveeqn.m (CSE)
Solves the wave equation u_tt=u_xx by the Leapfrog method. The example has a fixed end on the left, and a loose end on the right.
• Traffic flow solved by a characteristic particle method: mit18086_traffic.m
Shows the approximate weak entropy solution plus a set of cars moving on the road. Observe that the velocity of the cars differs from the velocity of information. Note that the computation is done by a particle method, which is a different approach as a fixed grid finite difference method.
• Nonlinear finite differences for the one-way wave equation with discontinuous initial conditions: mit18086_fd_transport_limiter.m (CSE)
Solves u_t+cu_x=0 by finite difference methods. Of interest are discontinuous initial conditions. The methods of choice are upwind, Lax-Friedrichs and Lax-Wendroff as linear methods, and as a nonlinear method Lax-Wendroff-upwind with van Leer and Superbee flux limiter.
• Shallow water equations solved by a particle method: mit18086_shallowwater.m
Shows the solution of the nonlinear shallow water equations with particles moving with the flow. As disturbances get small, the solution behaves like the linear wave equation with a base flow velocity.
• Level set method for front propagation under a given front velocity field: mit18086_levelset_front.m (CSE)
Uses the level set method with reinitialization to compute the movement of fronts under a given velocity field.
• Finite differences for the incompressible Navier-Stokes equations in a box: mit18086_navierstokes.m (CSE)
Solves the incompressible Navier-Stokes equations in a rectangular domain with prescribed velocities along the boundary. The standard setup solves a lid driven cavity problem.
This Matlab code is compact and fast, and can be modified for more general fluid computations. You can download a Documentation for the program.
• Sets up and solves a sparse system for the 1d, 2d and 3d Poisson equation: mit18086_poisson.m (CSE)
Sets up a sparse system by finite differences for the 1d Poisson equation, and uses Kronecker products to set up 2d and 3d Poisson matrices from it. The systems are solved by the backslash operator, and the solutions plotted for 1d and 2d.
• Computes the LU decomposition of a 2d Poisson matrix with different node ordering: mit18086_fillin.m (CSE)
Sets up a 2d Poisson problem and computes the LU decomposition of the system matrix, firstly with lexicographic ordering, then with red-black ordering and then with symamd ordering.
• Demo for elimination with reordering: moe.m
Visulizes elimination with various reordering methods for a 2d Poisson problem. Reordering methods are reverse Cuthill-McKee, minimum degree, symamd, and nested dissection.
• Computes the incomplete LU factorization of a 2d Poisson matrix for different tolerances: mit18086_ilu.m (CSE)
Sets up a 2d Poisson problem and computes the LU factorization and the incomplete LU factorization for different tolerances. Run times, condition numbers and nonzero pattern are plotted for comparison.
• Smoothing and convergence for Jacobi and Gauss-Seidel iteration: mit18086_smoothing.m (CSE)
Sets up a 1d Poisson problem with an oscillatory right hand side and applies Jacobi and Gauss-Seidel iteration to it. One can observe how the error gets smoothed very quickly, but the actual convergence to the correct solution is very slow.
• Multigrid solver for 1d Poisson problem: mit18086_multigrid.m (CSE)
Sets up a 1d Poisson test problem and solves it by multigrid. The method uses twogrid recursively using Gauss-Seidel for smoothing and elimination to solve at coarsest level. The number of pre- and postsmoothing and coarse grid iteration steps can be prescribed.
• Conjugate gradient method for 2d Poisson problem: mit18086_cg.m (CSE)
Sets up a 2d Poisson problem and solves the arising linear system by a conjugate gradient method. The error over iteration steps is plotted.

## Course Projects

The following project have been done in this course:
• Dongfang Bai: Solving piston secondary motion of internal combustion engines -> report as pdf
• Sudhish Kumar Bakku: Modeling elastic wave propagation in a layered medium with surface topography -> report as pdf
• Ishan Barman: Understanding constrained regularization to develop accurate multivariate calibration schemes -> report as pdf
• Hussam Busfar: Locating earthquakes and minimizing the error
• Haijie Chen.: Micro fluid surface tension driven flow
• Goh Chun Fan: Fluid flow in a changing domain
• Legena Henry: A simulation of a viscous boundary layer upon the introduction of vortices
• Tian Fook Kong: Modeling a microfluidic switch based on hydrodynamic spreading -> report as pdf
• Robert Legge: Fluid flow through porous media with different porosity gradients -> report as pdf
• Junlun Li: Coupling of focal mechanism with finite difference wave propagation in a rotated staggered grid -> report as pdf
• Robert Panish: Exploration of shock formation for a reentering space shuttle -> report as pdf
• Vernella Vickerman: Membrane dynamics during stimulation with chemotatic gradients
• Yinchun Wang: Simulation of droplets on surfaces using a level set method -> report as pdf

## Project Presentations

The projects for the courses 18.086 and 18.336 will be presented on the Computational Saturday, May 10th, 2008, from 9:40am to 3:15pm, in 2-132 and 2-136.