Scientific Computing

This module uses the DFT to help in understanding the convergence
rates of iterative methods for solving a system of linear equations.
Iterative methods are often used for solving systems of linear
equations *A*** x** =

The user begins by selecting the number of interior mesh points to
use in the numerical solution of the Laplace equation. The upper graph
displays an initial guess for the solution (blue dots), including the
boundary conditions *u*(0) = 0*u*(1)
= 0*Specific Initial
Guess*, which cycles through a set of seven different initial
guesses that illustrate properties of the different iterative methods,
or by clicking *Random Initial Guess*, which generates a random
initial guess for the solution, choosing a value between −1.5 and
1.5 at each interior mesh point. Finally, the user selects an
iterative method from the menu provided. For the SOR method, the user
can adjust the relaxation parameter *ω* using the slider.
The theoretically optimal value of *ω* for the given number
of mesh points is shown below the slider. The SOR method is
implemented here using the “red-black” ordering, meaning in
this case that on each iteration the solution components corresponding
to the odd-numbered mesh points are updated before those for the
even-numbered mesh points. Clicking *Iterate* carries out one
iteration of the chosen iterative method for solving the linear
system. After each iteration, the approximate solution in the upper
graph is updated to show the new values, which should eventually
converge to zero. The user can execute as many iterations of the
method as desired; an iteration counter keeps track of the number of
iterations, and the norm of the residual is displayed, so the user can
compare different methods after the same number of iterations.
Clicking *Reset* allows the user to select a different method for
the same initial guess or to change the number of mesh points or
starting values.

Throughout the iterations, the bottom graph shows the power spectrum
(red dots) of the first half of the DFT of the values displayed in the
top graph (interior mesh points only). The power spectrum of a
sequence of complex values is computed by multiplying each value by its
complex conjugate. For a sequence of real values, the power spectrum
of the DFT is symmetric, so only the first half of the power spectrum
need be displayed. The power spectrum represents the amount of energy
present at each frequency. The tick labels on the horizontal axis show
the wave numbers of selected frequencies; higher wave numbers
correspond to higher frequencies. As the tick labels on the vertical
axis suggest, the values of the power spectrum are shown on a log base
2 scale, with all values less than 2^{−5} lying on the
horizontal axis.

As watching the power spectrum during iterations reveals, some of the iterative methods quickly reduce some frequencies in the residual, but then convergence bogs down as other frequencies decay very slowly. This observation, combined with the fact that the apparent frequency of a component of the residual depends on the coarseness of the mesh, motivates multigrid methods, in which the residual is transferred between grids of varying refinement in order to obtain maximum benefit from an iterative method's abililty to suppress various frequencies most quickly.

**Reference:** Michael T. Heath, *Scientific Computing,
An Introductory Survey*, 2nd edition, McGraw-Hill, New York,
2002. See Section 11.5 and Computer Problem 12.14 on page 510.

**Developers:** Evan VanderZee and Michael Heath