Computational Mathematics and Applications
|
Thu, 29/04/2004 14:00 |
Dr Damien Jenkinson (University of Huddersfield) |
Computational Mathematics and Applications |
Rutherford Appleton Laboratory, nr Didcot |
Consider approximating a set of discretely defined values say at , with a chosen approximating form. Given prior knowledge that noise is present and that some might be outliers, a standard least squares approach based on norm of the error may well provide poor estimates. We instead consider a least squares approach based on a modified measure of the form , where is a constant to be fixed.
The choice of the constant in this estimator has a significant effect on the performance of the estimator both in terms of its algorithmic convergence to a solution and its ability to cope effectively with outliers. Given a prior estimate of the likely standard deviation of the noise in the data, we wish to determine a value of such that the estimator behaves like a robust estimator when outliers are present but like a least squares estimator otherwise.
We describe approaches to determining suitable values of and illustrate their effectiveness on approximation with polynomial and radial basis functions. We also describe algorithms for computing the estimates based on an iteratively weighted linear least squares scheme. |
|||
|
Thu, 06/05/2004 14:00 |
Dr Paul Dellar (University of Oxford) |
Computational Mathematics and Applications |
Comlab |
| The lattice Boltzmann equation has been used successfully used to simulate nearly incompressible flows using an isothermal equation of state, but much less work has been done to determine stable implementations for other equations of state. The commonly used nine velocity lattice Boltzmann equation supports three non-hydrodynamic or "ghost” modes in addition to the macroscopic density, momentum, and stress modes. The equilibrium value of one non-hydrodynamic mode is not constrained by the continuum equations at Navier-Stokes order in the Chapman-Enskog expansion. Instead, we show that it must be chosen to eliminate a high wavenumber instability. For general barotropic equations of state the resulting stable equilibria do not coincide with a truncated expansion in Hermite polynomials, and need not be positive or even sign-definite as one would expect from arguments based on entropy extremisation. An alternative approach tries to suppress the instability by enhancing the damping the non-hydrodynamic modes using a collision operator with multiple relaxation times instead of the common single relaxation time BGK collision operator. However, the resulting scheme fails to converge to the correct incompressible limit if the non-hydrodynamic relaxation times are fixed in lattice units. Instead we show that they must scale with the Mach number in the same way as the stress relaxation time. | |||
|
Thu, 13/05/2004 14:00 |
Dr Paul Matthews (University of Nottingham) |
Computational Mathematics and Applications |
Comlab |
|
The formation of steady patterns in one space dimension is generically
governed, at small amplitude, by the Ginzburg-Landau equation.
But in systems with a conserved quantity, there is a large-scale neutral
mode that must be included in the asymptotic analysis for pattern
formation near onset. The usual Ginzburg-Landau equation for the amplitude
of the pattern is then coupled to an equation for the large-scale mode.
These amplitude equations show that for certain parameters all regular periodic patterns are unstable. Beyond the stability boundary, there exist stable stationary solutions in the form of spatially modulated patterns or localised patterns. Many more exotic localised states are found for patterns in two dimensions. Applications of the theory include convection in a magnetic field, providing an understanding of localised states seen in numerical simulations. |
|||
|
Thu, 20/05/2004 14:00 |
Dr Brad Baxter (London) |
Computational Mathematics and Applications |
Comlab |
| We calculate an analytic value for the correlation coefficient between a geometric, or exponential, Brownian motion and its time-average, a novelty being our use of divided differences to elucidate formulae. This provides a simple approximation for the value of certain Asian options regarding them as exchange options. We also illustrate that the higher moments of the time-average can be expressed neatly as divided differences of the exponential function via the Hermite-Genocchi integral relation, as well as demonstrating that these expressions agree with those obtained by Oshanin and Yor when the drift term vanishes. | |||
|
Thu, 27/05/2004 14:00 |
Dr Miguel Anjos (University of Southampton) |
Computational Mathematics and Applications |
Comlab |
| The satisfiability (SAT) problem is a central problem in mathematical logic, computing theory, and artificial intelligence. We consider instances of SAT specified by a set of boolean variables and a propositional formula in conjunctive normal form. Given such an instance, the SAT problem asks whether there is a truth assignment to the variables such that the formula is satisfied. It is well known that SAT is in general NP-complete, although several important special cases can be solved in polynomial time. Extending the work of de Klerk, Warners and van Maaren, we present new linearly sized semidefinite programming (SDP) relaxations arising from a recently introduced paradigm of higher semidefinite liftings for discrete optimization problems. These relaxations yield truth assignments satisfying the SAT instance if a feasible matrix of sufficiently low rank is computed. The sufficient rank values differ between relaxations and can be viewed as a measure of the relative strength of each relaxation. The SDP relaxations also have the ability to prove that a given SAT formula is unsatisfiable. Computational results on hard instances from the SAT Competition 2003 show that the SDP approach has the potential to complement existing techniques for SAT. | |||
|
Thu, 03/06/2004 14:00 |
Prof Paul Houston (University of Leicester) |
Computational Mathematics and Applications |
Comlab |
|
In recent years, there has been considerable interest, especially in the context of
fluid-dynamics, in nonconforming finite element methods that are based on discontinuous
piecewise polynomial approximation spaces; such approaches are referred to as discontinuous
Galerkin (DG) methods. The main advantages of these methods lie in their conservation properties, their ability to treat a wide range of problems within the same unified framework, and their great flexibility in the mesh-design. Indeed, DG methods can easily handle non-matching grids and non-uniform, even anisotropic, polynomial approximation degrees. Moreover, orthogonal bases can easily be constructed which lead to diagonal mass matrices; this is particularly advantageous in unsteady problems. Finally, in combination with block-type preconditioners, DG methods can easily be parallelized.
In this talk, we introduce DG discretizations of mixed field and potential-based formulations of eddy current problems in the time-harmonic regime. For the electric field formulation, the divergence-free constraint within non-conductive regions is imposed by means of a Lagrange multiplier. This allows for the correct capturing of edge and corner singularities in polyhedral domains; in contrast, additional Sobolev regularity must be assumed in the DG formulation, and their conforming counterparts, when regularization techniques are employed. In particular, we present a mixed method involving discontinuous elements, which includes a normal jump stabilization term, and a non-stabilized variant employing discontinuous elements.The first formulation delivers optimal convergence rates for the vector-valued unknowns in a suitable energy norm, while the second (non-stabilized) formulation is designed to yield optimal convergence rates in both the –norm, as well as in a suitable energy norm. For this latter method, we also develop the a posteriori error estimation of the mixed DG approximation of the Maxwell operator. Indeed, by employing suitable Helmholtz decompositions of the error, together with the conservation properties of the underlying method, computable upper bounds on the error, measured in terms of the energy norm, are derived.
Numerical examples illustrating the performance of the proposed methods will be presented; here, both conforming and non-conforming (irregular) meshes will be employed. Our theoretical and numerical results indicate that the proposed DG methods provide promising alternatives to standard conforming schemes based on edge finite elements. |
|||
|
Thu, 10/06/2004 14:00 |
Dr Serge Gratton (CERFAX Toulouse) |
Computational Mathematics and Applications |
Rutherford Appleton Laboratory, nr Didcot |
|
We consider the solution of a linear system of equations using the GMRES iterative method. In some applications, performing inexact matrix-vector products in this method may be interesting, provided that a reasonable convergence of GMRES is achieved. A GMRES algorithm where the matrix vector product is performed inexactly is termed ”inexact GMRES algorithm”. An application of this idea occurs in computational electromagnetics, where the fast multipole method provides approximations of the matrix-vector product within a user-defined precision, and where these inaccurate matrix-vector products are all the more cheaper (in terms of CPU time) as the user-defined precision is low. The key point is then to design a control strategy of the accuracy of the matrix-vector product so that the GMRES converges better in the sense that 1) the inexact method achieves a satisfactory limiting accuracy, 2) within a reasonable number of steps.
In [1], a relaxation strategy is proposed for general systems and validated on a large set of numerical experiments. This work is based on heuristic considerations and proposes a strategy that enables a convergence of the GMRES iterates within a relative normwise backward error less than a prescribed quantity > 0, on a significant number of numerical experiments. Similar strategies have been applied to the solution of device simulation problems using domain decomposition [2] and to the preconditioning of a radiation diffusion problem in [5].
A step toward a theoretical explanation of the observed behaviour of the inexact GMRES is proposed in [3, 4]. In this talk, we show that in spite of this considerable theoretical study, the experimental work of [1] is not fully understood yet. We give an overview of the questions that still remains open both in exact arithmetic and in floating-point arithmetic, and we provide some insights into the solution of some of them. Possible applications of this work for the preconditioned GMRES method, when the matrix-vector product is accurate but the preconditioning operation is approximate, are also investigated, based on [3]. References [1] A. Bouras and V. Frayss´e. Inexact matrix-vector products in Krylov methods for solving linear systems: a relaxation strategy. SIAM Journal on Matrix Analysis and Applications, 2004. To appear. [2] A. Bouras, V. Frayss´e, and L. Giraud. A relaxation strategy for inner-outer linear solvers in domain decomposition methods. Technical Report TR/PA/00/17, CERFACS, Toulouse, France, 2000. [3] V. Simoncini and D. B. Szyld. Theory of inexact Krylov subspace methods and applications to scientific computing. SIAM Journal Scientific Computing, 25:454–477, 2003. [4] J. van den Eshof and G. L. G. Sleijpen. Inexact Krylov subspace methods for linear systems. SIAM Journal on Matrix Analysis and Applications, February 2004. To appear. [5] J. S. Warsa, M. Benzi, T. A. Warein, and J. E. Morel. Preconditioning a mixed discontinuous finite element method for radiation diffusion. Numerical Linear Algebra with Applications, 2004. To appear. |
|||
|
Tue, 15/06/2004 14:00 |
Dr Yifan Hu (Wolfram Research) |
Computational Mathematics and Applications |
Rutherford Appleton Laboratory, nr Didcot |
| The talk will start with an introduction to recent development in Mathematica, with emphasis on numerical computing. This will be followed by a discussion of graph drawing algorithms for the display of relational information, in particular force directed algorithms. The talk will show that by employing multilevel approach and octree data structure, it is possible to achieve fast display of very large relational information, without compromising the quality. | |||
|
Thu, 17/06/2004 14:00 |
Prof Gilbert Strang (MIT) |
Computational Mathematics and Applications |
Comlab |
An essential first step in many problems of numerical analysis and
computer graphics is to cover a region with a reasonably regular mesh.
We describe a short MATLAB code that begins with a "distance function"
to describe the region: is the distance to the boundary
(with d < 0 inside the region). The algorithm decides the location
of nodes and the choice of edges. At each iteration the mesh becomes a
truss (with forces in the edges that move the nodes). Then the Delaunay
algorithm resets the topology (the choice of edges). The input includes
a size function to vary the desired meshlengths over the region.
The code is adaptable and public (math.mit.edu/~persson/mesh). It extends directly to dimensions. The 2D meshes are remarkably regular. In higher
dimensions it is difficult to avoid "slivers". The theory is not complete.
We look for the "right proof" of a known result about the inverse of a tridiagonal matrix : Every 2 by 2 submatrix on and above the diagonal
of inv( ), or on and below the diagonal, has determinant = zero. This is
the case of a general result for (band + low rank) matrices:
All submatrices above the th superdiagonal of have rank( ) < if
and only if all submatrices above the th subdiagonal of inv( ) have
rank( ) < .
A quick proof uses the Nullity Theorem. The complementary subspaces and have the same nullity. These full matrices inv( ) appear in
integral equations and in the boundary element method. The large low rank
submatrices allow fast multiplication. The theory of decay away from
the main diagonal could extend to include the property of low rank. |
|||

say at
, with a chosen approximating form. Given prior knowledge that noise is present and that some might be outliers, a standard least squares approach based on
norm of the error
may well provide poor estimates. We instead consider a least squares approach based on a modified measure of the form
, where
is a constant to be fixed.
elements, which includes a normal jump stabilization term, and a non-stabilized variant employing discontinuous
elements.The first formulation delivers optimal convergence rates for the vector-valued unknowns in a suitable energy norm, while the second (non-stabilized) formulation is designed to yield optimal convergence rates in both the
–norm, as well as in a suitable energy norm. For this latter method, we also develop the a posteriori error estimation of the mixed DG approximation of the Maxwell operator. Indeed, by employing suitable Helmholtz decompositions of the error, together with the conservation properties of the underlying method, computable upper bounds on the error, measured in terms of the energy norm, are derived.
within a relative normwise backward error
less than a prescribed quantity
> 0, on a significant number of numerical experiments. Similar strategies have been applied to the solution of device simulation problems using domain decomposition [2] and to the preconditioning of a radiation diffusion problem in [5].
is the distance to the boundary
(with d < 0 inside the region). The algorithm decides the location
of nodes and the choice of edges. At each iteration the mesh becomes a
truss (with forces in the edges that move the nodes). Then the Delaunay
algorithm resets the topology (the choice of edges). The input includes
a size function
to vary the desired meshlengths over the region.
dimensions. The 2D meshes are remarkably regular. In higher
dimensions it is difficult to avoid "slivers". The theory is not complete.
: Every 2 by 2 submatrix on and above the diagonal
of inv(
of a general result for (band + low rank) matrices:
All submatrices
above the
th superdiagonal of
if
and only if all submatrices
above the
.
and