Past Computational Mathematics and Applications Seminar

No seminar this week

Let the thin plate spline radial basis function method be applied to

interpolate values of a smooth function $f(x)$, $x \!\in\! {\cal R}^d$.

It is known that, if the data are the values $f(jh)$, $j \in {\cal Z}^d$,

where $h$ is the spacing between data points and ${\cal Z}^d$ is the

set of points in $d$ dimensions with integer coordinates, then the

accuracy of the interpolant is of magnitude $h^{d+2}$. This beautiful

result, due to Buhmann, will be explained briefly. We will also survey

some recent findings of Bejancu on Lagrange functions in two dimensions

when interpolating at the integer points of the half-plane ${\cal Z}^2

\cap \{ x : x_2 \!\geq\! 0 \}$. Most of our attention, however, will

be given to the current research of the author on interpolation in one

dimension at the points $h {\cal Z} \cap [0,1]$, the purpose of the work

being to establish theoretically the apparent deterioration in accuracy

at the ends of the range from ${\cal O} ( h^3 )$ to ${\cal O} ( h^{3/2}

)$ that has been observed in practice. The analysis includes a study of

the Lagrange functions of the semi-infinite grid ${\cal Z} \cap \{ x :

x \!\geq\! 0 \}$ in one dimension.

In this talk we review experiences of using the Harwell Subroutine

Library and other numerical software codes in implementing large scale

solvers for commercial industrial process simulation packages. Such

packages are required to solve problems in an efficient and robust

manner. A core requirement is the solution of sparse systems of linear

equations; various HSL routines have been used and are compared.

Additionally, the requirement for fast small dense matrix solvers is

examined.

The sphere is an important setting for applied mathematics, yet the underlying approximation theory and numerical analysis needed for serious applications (such as, for example, global weather models) is much less developed than, for example, for the cube.

\\

\\

This lecture will apply recent developments in approximation theory on the sphere to two different problems in scientific computing.

\\

\\

First, in geodesy there is often the need to evaluate integrals using data selected from the vast amount collected by orbiting satellites. Sometimes the need is for quadrature rules that integrate exactly all spherical polynomials up to a specified degree $n$ (or equivalently, that integrate exactly all spherical harmonies $Y_{\ell ,k}(\theta ,\phi)$ with $\ell \le n).$ We shall demonstrate (using results of M. Reimer, I. Sloan and R. Womersley in collaboration with

W. Freeden) that excellent quadrature rules of this kind can be obtained from recent results on polynomial interpolation on the sphere, if the interpolation points (and thus the quadrature points) are chosen to be points of a so-called extremal fundamental system.

\\

\\

The second application is to the scattering of sound by smooth three-dimensional objects, and to the inverse problem of finding the shape of a scattering object by observing the pattern of the scattered sound waves. For these problems a methodology has been developed, in joint work with I.G. Graham, M. Ganesh and R. Womersley, by applying recent results on constructive polynomial approximation on the sphere. (The scattering object is treated as a deformed sphere.)

Continuing advances in computing technology provide the power not only to solve

increasingly large and complex process modeling and optimization problems, but also

to address issues concerning the reliability with which such problems can be solved.

For example, in solving process optimization problems, a persistent issue

concerning reliability is whether or not a global, as opposed to local,

optimum has been achieved. In modeling problems, especially with the

use of complex nonlinear models, the issue of whether a solution is unique

is of concern, and if no solution is found numerically, of whether there

actually exists a solution to the posed problem. This presentation

focuses on an approach, based on interval mathematics,

that is capable of dealing with these issues, and which

can provide mathematical and computational guarantees of reliability.

That is, the technique is guaranteed to find all solutions to nonlinear

equation solving problems and to find the global optimum in nonlinear

optimization problems. The methodology is demonstrated using several

examples, drawn primarily from the modeling of phase behavior, the

estimation of parameters in models, and the modeling, using lattice

density-functional theory, of phase transitions in nanoporous materials.

This talk reviews some recent joint work with Michael Eiermann and Olaf

Schneider which introduced a framework for analyzing some popular

techniques for accelerating restarted Krylov subspace methods for

solving linear systems of equations. Such techniques attempt to compensate

for the loss of information due to restarting methods like GMRES, the

memory demands of which are usually too high for it to be applied to

large problems in unmodified form. We summarize the basic strategies which

have been proposed and present both theoretical and numerical comparisons.

Support Vector Machines are a new and very promising approach to

machine learning. They can be applied to a wide range of tasks such as

classification, regression, novelty detection, density estimation,

etc. The approach is motivated by statistical learning theory and the

algorithms have performed well in practice on important applications

such as handwritten character recognition (where they currently give

state-of-the-art performance), bioinformatics and machine vision. The

learning task typically involves optimisation theory (linear, quadratic

and general nonlinear programming, depending on the algorithm used).

In fact, the approach has stimulated new questions in optimisation

theory, principally concerned with the issue of how to handle problems

with a large numbers of variables. In the first part of the talk I will

overview this subject, in the second part I will describe some of the

speaker's contributions to this subject (principally, novelty

detection, query learning and new algorithms) and in the third part I

will outline future directions and new questions stimulated by this

research.

The flow in a cylinder with a rotating endwall has continued to

attract much attention since Vogel (1968) first observed the vortex

breakdown of the central core vortex that forms. Recent experiments

have observed a multiplicity of unsteady states that coexist over a

range of the governing parameters. In spite of numerous numerical and

experimental studies, there continues to be considerable controversy

with fundamental aspects of this flow, particularly with regards to

symmetry breaking. Also, it is not well understood where these

oscillatory states originate from, how they are interrelated, nor how

they are related to the steady, axisymmetric basic state.

\\

\\

In the aspect ratio (height/radius) range 1.6 2.8. An efficient and

accurate numerical scheme is presented for the three-dimensional

Navier-Stokes equations in primitive variables in a cylinder. Using

these code, primary and secondary bifurcations breaking the SO(2)

symmetry are analyzed.

\\

\\

We have located a double Hopf bifurcation, where an axisymmetric limit

cycle and a rotating wave bifurcate simultaneously. This codimension-2

bifurcation is very rich, allowing for several different scenarios. By

a comprehensive two-parameter exploration about this point we have

identified precisely to which scenario this case corresponds. The mode

interaction generates an unstable two-torus modulate rotating wave

solution and gives a wedge-shaped region in parameter space where the

two periodic solutions are both stable.

\\

\\

For aspect ratios around three, experimental observations suggest that

the first mode of instability is a precession of the central vortex

core, whereas recent linear stability analysis suggest a Hopf

bifurcation to a rotating wave at lower rotation rates. This apparent

discrepancy is resolved with the aid of the 3D Navier-Stokes

solver. The primary bifurcation to an m=4 traveling wave, detected by

the linear stability analysis, is located away from the axis, and a

secondary bifurcation to a modulated rotating wave with dominant modes

m=1 and 4, is seen mainly on the axis as a precessing vortex breakdown

bubble. Experiments and the linear stability analysis detected

different aspects of the same flow, that take place in different

spatial locations.

We combine linear algebra techniques with finite element techniques to obtain a reliable stopping criterion for the Conjugate Gradient algorithm. The finite element method approximates the weak form of an elliptic partial differential equation defined within a Hilbert space by a linear system of equations A x = b, where A is a real N by N symmetric and positive definite matrix. The conjugate gradient method is a very effective iterative algorithm for solving such systems. Nevertheless, our experiments provide very good evidence that the usual stopping criterion based on the Euclidean norm of the residual b - Ax can be totally unsatisfactory and frequently misleading. Owing to the close relationship between the conjugate gradient behaviour and the variational properties of finite element methods, we shall first summarize the principal properties of the latter. Then, we will use the recent results of [1,2,3,4]. In particular, using the conjugate gradient, we will compute the information which is necessary to evaluate the energy norm of the difference between the solution of the continuous problem, and the approximate solution obtained when we stop the iterations by our criterion.

Finally, we will present the numerical experiments we performed on a selected ill-conditioned problem.

References

• [1] M. Arioli, E. Noulard, and A. Russo, Vector Stopping Criteria for Iterative Methods: Applications to PDE's, IAN Tech. Rep. N.967, 1995.
• [2] G.H. Golub and G. Meurant, Matrices, moments and quadrature II; how to compute the norm of the error in iterative methods, BIT., 37 (1997), pp.687-705.
• [3] G.H. Golub and Z. Strakos, Estimates in quadratic formulas, Numerical Algorithms, 8, (1994), pp.~241--268.
• [4] G. Meurant, The computation of bounds for the norm of the error in the conjugate gradient algorithm, Numerical Algorithms, 16, (1997), pp.~77--87.

Boundary layers are often studied with no pressure gradient

or with an imposed pressure gradient. Either of these assumptions

can lead to difficulty in obtaining solutions. A major advance in fluid

dynamics last century (1969) was the development of a triple deck

formulation for boundary layers where the pressure is not

specified but emerges through an interaction between

boundary layer and the inviscid outer flow. This has given rise to

new computational problems and computations have in turn

fed ideas back into theoretical developments. In this survey talk

based on my new book, I will look at three problems:

flow past a plate, flow separation and flow in channels

and discuss the interaction between theory and computation

in advancing boundary layer theory.

This talk reviews some theoretical and practical aspects

of incompressible flow modelling using finite element approximations

of the (Navier-) Stokes equations.

The infamous Q1-P0 velocity/pressure mixed finite element approximation

method is discussed. Two practical ramifications of the inherent

instability are focused on, namely; the convergence of the approximation

with and without regularisation, and the behaviour of fast iterative

solvers (of multigrid type) applied to the pressure Poisson system

that arises when solving time-dependent Navier-Stokes equations

using classical projection methods.

\\

\\

This is joint work with David Griffiths from the University of Dundee.

We examine the convergence characteristics of iterative methods based

on a new preconditioning operator for solving the linear systems

arising from discretization and linearization of the Navier-Stokes

equations. With a combination of analytic and empirical results, we

study the effects of fundamental parameters on convergence. We

demonstrate that the preconditioned problem has an eigenvalue

distribution consisting of a tightly clustered set together with a

small number of outliers. The structure of these distributions is

independent of the discretization mesh size, but the cardinality of

the set of outliers increases slowly as the viscosity becomes smaller.

These characteristics are directly correlated with the convergence

properties of iterative solvers.