Norman Bleistein <norm@dix.mines.edu> Colorado School of Mines, Golden, CO
Abstract
In an article in 1954, Hagedoorn described a construction technique that created an image of a reflector (in the earth) from acoustic reflection data observed at the boundary (surface of the earth). It took twenty years and more for us analyst-types to derive a mathematical formulation that translated that constructive technique into mathematical formulas that could be used to process data in a computer to produce an image of reflectors and to characterize the reflection strength across those reflectors. In this talk, I present a version of Hagedoorn's construction that I called modeling and migration (reflector imaging) on a shoestring. I will explain how this shoestring construction leads to mathematical formulas for data processing techniques known as Kirchhoff migration and inversion. (To me, the distinction in those designations is a matter of whether the objective is a reflector map only [migration], or a map plus some information about medium parameter variations across reflectors [inversion].) These formulas, in turn, point the way towards computer processing of reflection data to produce images of reflectors (in the earth). Asymptotic analysis of the processing algorithm leads to an interpretation of the formula as a local Fourier synthesis. The source/receiver configuration and frequency domain bandwidth generate an aperture in wave-vector domain (k-domain). Imaging occurs when the normal direction to a reflector lies in a direction in the k-domain aperture.
J.U. Brackbill <jub@lanl.gov> Los Alamos National Laboratory, NM
Abstract
Peskin's [10] computational models of heart valves are an innovation in the treatment of flows in regions with immersed boundaries that allows one to model complex, time-dependent geometries using simple, rectangular grids. The algorithm is similar to the particle-in-cell method in its use of finite-size, Lagrangian markers to resolve contact discontinuities that move with the flow [9]. The immersed boundary technique has provided straightforward approaches to the solution of many difficult and cumbersome problems. Among them are the modeling of phase change using front tracking [6], suspension flows [5], plasma reactors [2], hydrodynamic rams [4], and granular flows [3]. In modeling granular flows, for example, the immersed boundary technique allows one to model deformable grains as they roll and slide past one another with an O(N) algorithm, where N is the number of computational points.
Although the immersed boundary technique is enormously useful, it could be improved, especially in the accuracy of the treatment of interfaces. Its low order is most evident when there are strong contrasts in the properties of the media in contact at the interface. A review of the ideas that have been explored, specifically for the immersed boundary technique and other related techniques, reveals that the interface problem bears similarity to that posed by the modeling of shocks, but without the benefit of a Riemann solution.
References
- 1
- D. Juric and G. Tryggvason, ``Computations of Boiling Flows,'' Int. J. Multiphase Flow, in press.
- 2
- G. Lapenta and J. U. Brackbill, ``Immersed boundary method for plasma simulation in complex geometries'', IEEE Trans. Plasma Sci., 24, 105, 1996.
- 3
- S. G. Bardenhagen and J. U. Brackbill, ``Dynamic stress bridging in granular material'', J. Appl. Phys, 83, 5732, 1998.
- 4
- B. Kashiwa, Los Alamos National Laboratory, 1999, private communication.
- 5
- D. Sulsky and J. U. Brackbill, ``A numerical method for suspension flows'', J. of Comput. Phys., 96, 339, 1991.
- 6
- D. Juric and G. Tryggvason, ``A front-tracking method for dendritic solidification'', J. Comput. Phys., 123, 127, 1996.
- 7
- D. Juric, ``On the Interface Instability During Rapid Evaporation in Microgravity,'' in Heat Transfer in Microgravity Systems, , ASME Heat Transfer Division (1997).
- 8
- D. Torres and J. U. Brackbill, "An Improved Continuum Surface Force Method", J. Comput. Phys., in preparation.
- 9
- F. H. Harlow, Los Alamos Scientific Laboratory Report No. LA-2301, 1959 (unpublished).
- 10
- C. S. Peskin, ``Flow patterns around heart valves: A numerical method'', J. Comput. Phys., 10, 252, 1972.
Reuben Hersh <rhersh@math.unm.edu> University of New Mexico, NM
Abstract
Motivated by the first-order accurate Euler method and the second-order Crank-Nicholson scheme, we find stable rational approximations of arbitrarily high order to all well-posed linear initial-value problems. Co-author: T. Kato.
James M. Hyman <jh@lanl.gov> Los Alamos National Laboratory, NM
Abstract
We have been investigating new approaches to reduce the errors introduced by splitting methods for the numerical solution of partial differential equations (PDEs) arising in coupled radiation-hydrodynamics, reaction diffusion systems and nonlinear dispersive wave equations. The numerical solution of these problems often can be simplified by splitting the differential operator into several pieces, where each piece is much easier to solve than the original full system. These splitting algorithms are accurate and efficient, except when the spilt equations are strongly coupled and act in nearly equal, but opposite directions. When this happens, it can lead to instabilities, or even worse, the split problem can be stable but inaccurate, causing the simulation to mislead an investigator studying the dynamics of the mathematical model.
The new approach is based on using an approximate solution operator (ASO) to reformulate the splitting method so the error is proportional to the product of the original splitting error and the difference between the true solution and an approximate solution. Thus, the more accurate the ASO, the smaller the numerical error will be. This information (such as a back of the envelope calculation of the asymptotic shape of a boundary layer or the WKB analysis of an oscillating system) can be explicitly input by the programmer or dynamically generated during the solution process. The better the mathematical analysis, the more efficient the PDE code will be.
Furthermore, because the solution at every stage of an ASO-splitting method is an accurate approximation of the solution of the original unsplit problem, many of the inaccuracies and instabilities the plague splitting methods with time dependent coefficients or time dependent boundary conditions are eliminated.
Herb B. Keller <hbk@newvortex.ama.caltech.edu> Applied Mathematics 217-50, California Institute of Technology, Pasadena, CA 91125
Abstract
A complete existence theory for initial value problems for differential algebraic systems in the form
has been developed. The basic requirement is that the total Jacobian matrix:
have maximal rank over an open set containing the initial data. We show that regular systems of ordinary differential equations can be obtained which have unique solutions over the open set. Depending upon the rank of
, the dimensionality
of these regular systems may be less than that of 
.
Extensions of this work allow the rank of 
to change and then
nonuniqueness occurs as the solutions may exhibit folds or bifurcations. An
interesting example illustrating this and even more weird behavior has been
constructed.
Herb Kranzer <kranzer@panther.adelphi.edu> Adelphi University, Garden City, NY
Abstract
We present some examples of systems of conservation laws for which the Riemann problem has no classical weak solution, and we show how a solution can be obtained by using singular shocks. A singular shock is a moving discontinuity with an unbounded, delta-function-like internal structure. We discuss the nature of singular shocks and explore conditions for their admissibility as weak solutions.
Heinz Kreiss <kreiss@joshua.math.ucla.edu> University of California Los Angeles, CA
Abstract
In practice, we often solve partial differential equations numerically without any knowledge whether the underlying analytic problem has a solution. The question therefore is: Can we use the numerical results to decide or, at least, make plausible that the analytic problem, in fact, has a solution which is close to the numerical solution. We will see that, in many cases, we can interpolate the numerical solution and show that the interpolant solves the problem for slightly perturbed data. This is as close as numerical methods can get us. To decide whether the original problem has a solution is a matter of the local stability. We shall discuss different techniques to decide this question and apply it to recent turbulence calculations.
Peter Lax <lax@cims.nyu.edu> Courant Institute of Mathematical Sciences New York University, NY
Abstract
Much effort was spent on trying to devise TVD schemes for multidimensional conservation laws. The search is, of course, doomed to failure, since TV does not D in more than one dimension. The only thing in more than one dimension that does D is energy. According to an important theorem of Friedrichs, for symmetric hyperbolic systems with smooth coefficients schemes with positive weights coserve energy.The author and Xu-Dong Liu have devised such a scheme for the Euler equations of compressible flow that is of second order accuracy in smooth regions; a switch is used to turn off the accurate scheme in regions of discontinuities. We make use of the Roe matrix.
The basic code is very compact, less than 200 lines. Test calculations on basic benchmark problems (diffraction of a strong shock wave by a wedge, shock propagating in a partially obstructed tunnel, two dimensional Riemann problems) show our method to be competitive for accuracy and resolution with the best methods available.
Randy LeVeque <rjl@amath.washington.edu> University of Washington, Seattle, WA
Abstract
A class of wave propagation algorithms for three-dimensional conservation laws and other hyperbolic systems has been developed with Jan Olav Langseth and implemented in the CLAWPACK software. These unsplit finite volume methods are based on solving one-dimensional Riemann problems at the cell interfaces. Waves emanating from the Riemann problem are further split by solving Riemann problems in the transverse directions to model cross-derivative terms. These methods are related to the multi-dimensional Lax-Wendroff method, but proper upwinding can yield stability for Courant numbers up to one. The methods are also written in a form where limiters can be applied to suppress spurious oscillations and give sharp shock capturing.
Richard Liska <liska@siduri.fjfi.cvut.cz> Czech Technical University, Prague, Czech Republic
Abstract
Composite schemes are given by several steps of dispersive second order Lax-Wendroff (LW) difference scheme followed by one step of diffusive first order Lax-Friedrichs (LF) step which serves as a consistent filter removing well-known Lax-Wendroff dispersive oscillations behind the shocks. We have found this simple idea while experimenting with filters on 1D shallow water problems. It turned out that the composite schemes work well also in higher dimensions and for other systems of conservation laws. New versions of non-split two step LF and LW type schemes have been developed in 2D and 3D on orthogonal grids. In 2D the schemes are optimally stable while in 3D we had troubles with stability of the LW type scheme. In 2D the schemes can be generalized to arbitrary logically rectangular grids. Recently the schemes have been formulated also on unstructured meshes with triangular and staggered n-lateral grids. The usefulness of the composite schemes will be demonstrated on several examples from Euler equations in 2D and 3D, shallow water equations in 2D in planar geometry and on the surface of a rotating sphere. This is a joint work with Burt Wendroff. Mikhail Shashkov, Karel Kozel and Michal Janda have contributed to the part using triangular grids.
Len Margolin <len@lanl.gov> Los Alamos National Laboratory, NM
Piotr Smolarkiewicz <smolar@ucar.edu> National Center for Atmospheric Research, Boulder, CO
Abstract
We explore the ability of nonoscillatory advection schemes to represent the effects of the unresolved scales of motion in numerical simulation of turbulent flows. We demonstrate that a nonoscillatory fluid solver can accurately reproduce the dynamics of an atmospheric convective boundary layer. When an explicit turbulence model is implemented, the solver does not add any significant numerical diffusion. Of greater interest, when no explicit turbulence model is implemented, the solver itself appears to include an effective subgrid scale model. Other researchers have reported similar success, simulating turbulent flows in a variety of regimes while using only nonoscillatory advection schemes to model subgrid effects. At this point there is no theory to justify this success, but we offer some speculations.
Bill Morton <bill.morton@comlab.ox.ac.uk> Oxford University Computing Laboratory, Oxford, UK
Abstract
In the design of difference schemes one clearly needs to preserve key properties of the underlying PDE's; and the massive impact that Lax-Wendroff methods have had on our field springs from the emphasis that they placed on using the conservation law form, i.e. conserving the primary variables. Then in 1966 Arakawa showed how important quadratic functions of the primary variables such as energy and enstrophy could be preserved, before Jesperson demonstrated that this followed naturally from using finite element approximations in a Galerkin formulation.
Little attention seems to have been paid, however, to the conservation of derivative functions of the primary variables, such as vorticity. The generic model problem is the wave equation system - as distinct from the single second order equation. One of the difficulties is that one has to prescribe the discretisation of the vorticity that is to be conserved. It turns out that on a uniform square mesh, with the most compact representation of discrete vorticity, there is a unique second order, nine-point, vorticity-preserving scheme and it is the Rotated-Richtmyer two-step Lax-Wendroff method.
Moreover, this scheme can be generalised and applied on a completely unstructured mesh in two or three dimensions. It can also be applied to a variable coefficient problem; so, for instance, for a wave equation with variable density one obtains a discrete version of Kelvin's Theorem.
There are two dual interpretations of the Rotated-Richtmyer scheme on a uniform mesh. But only one of these generalises as described above; and this is a finite volume scheme using time-averaged fluxes through the mesh boundaries. For a Lax-Wendroff scheme these are, of course, obtained from Taylor series expansions. There are alternative means of doing this, however. One possibility is to use characteristics or, for the wave equation, the bicharacteristic cone. This links the vorticity-preserving ideas into a rich variety of approximation methods based on a finite volume formulation. We will indicate some of these in the lecture.
Stan Osher <sjo@joshua.math.ucla.edu> University of California Los Angeles, CA
Abstract
We shall discuss the construction of high resolution methods for compressible multiphase flow. We begin with the LW scheme obtained via interpolation, then add TVD, ENO, an WENO ideas to increase accuracy and remove spurious oscillations. Next the interface is treated using an elementary version of the level set method together with th ghost fluid method which also removes spurious oscillations, this time at material interfaces. We also discuss a recent result which implies that some advected quantities can be conveniently aproximated in nonconsevartive form and still converge to a weak solution, thus generalizing the classical Lax-Wendroff theorem. Most of the recent work discussed above is joint with Ron Fedkiw, Barry Merriman, and Tariq Aslam.
Seymour Parter <parter@cs.wisc.edu> University of Wisconsin, Madison, WI
Abstract
Spectral Collocation Methods for Elliptic Problems are very accurate.In the case of smooth solutions they provide much better approximations than either traditional Finite Difference or traditional Finite Element Methods. Unfortunately the accompaning linear systems are not at all nice.The matrices are full,or at least have no regular patterns of zeros.Further,these matrices are badly ill conditioned.For over twenty years the preffered approach to the solution of these problems have been based on preconditioning by finite difference or finite element operators based on the collocation points.In practice,these methods have proven to be very effective.But, until recently there has been little in the way of mathematical justification. In this talk we discuss the method,what is known, and the challenges for the future.
James J. Quirk <quirk@lanl.gov> Computer Research and Applications Group (CIC-3), Los Alamos National Laboratory, NM
Abstract
Computational fluid dynamics (CFD) is a complex amalgam of numerical analysis, flow physics, computer science, not to mention common sense. And as a consequence, important practical observations are often slow to percolate through the community, resulting in unnecessary scientific confusion and wasted programming effort. This talk will provide an historical perspective on one such communication breakdown: the emergence of pathological failings - such as the so-called carbuncle phenomenon - that plague popular shock-capturing schemes. The aim will be to look beyond the classical debate of whether scheme A is better than scheme B (or vice versa) so as to identify concrete mechanisms which could be used to inject much needed rigour into the darker recesses of CFD. The hope being that such an injection could, in the fullness of time, elevate key activities - like validation & verification - from art forms into bona fide sciences.
Bill Rider <wjr@lanl.gov> Los Alamos National Laboratory, NM
Len Margolin <len@lanl.gov> Los Alamos National Laboratory, NM
Abstract
The advent of nonlinear numerical methods for hyperbolic PDEs has heralded an era where nonoscillatory results are commonplace. Here, we discuss a class of nonlinear methods that are sign-preserving. These methods can be developed from several seemingly independent starting points. The original method in this class is MPDATA. We will show how this method is related to both the standard viscous and to a nonlinearly dispersive-dissipative regularization. Our discussion will touch upon the relation of this form of regularization to several physical mechanisms for dissipation. Additionally, our method is strongly related to classical schemes using artificial viscosity. We will discuss how one can transform our scheme into a monotone (TVD) method, and the consequent production of (locally) entropy-satisfying solutions. Finally, we demonstrate the effectiveness of this method on a broad class of problems.
Mikhail Shashkov <misha@t7.lanl.gov> MS-B284, Group T-7, Theoretical Division, Los Alamos National Laboratory, Los Alamos, New Mexico, USA, 87545
Abstract
In past 6 years we have developed new high-quality, mimetic finite-difference methods based on discrete analog of vector and tensor analysis (DVTA). The basis of DVTA is the design of discrete operators that preserve certain essential properties of, and relationships between, the corresponding analytic operators. The DVTA is the basis for new techniques for large-scale numerical simulations approximating the solution of partial differential equations. The new methods provide a significant extension of the well known and useful finite volume methods and are designed to more faithfully represent important properties of physical processes and the continuum mathematical models of such processes. Algorithms based on these techniques are used for modeling high-speed flows, porous media flows, diffusion processes, and electromagnetic problems. In this presentation we will describe DVTA and demonstrate how it can be used to construct high-quality finite-difference methods for Maxwell's equations.
Joel Smoller <smoller@umich.edu> University of Michigan, Ann Arbor, MI
Abstract
We consider the simplest solution of the Einstein equations that incorporates a shock-wave into a standard Friedmann-Robertson-Walker meric whose equation of state accounts for the Hubble constant and the microwave background radiation temperature. This produces a new solution of the Einstein equations from which we are able to show that the distance from the shock-wave to the center of the explosion at present time is comparable to the Hubble distance. We are motivated by the idea that the expansion of the universe as measured by the Hubble constant might be accounted for by an event more similiar to a classical explosion than by the well-accepted scenario of the Big-Bang.
Bruce Stewart <bstewart@users.buoy.com>
Abstract
Exploratory research sometimes benefits from the use of custom software. We describe the principles and practicalities of software for the geometric study of dynamical systems in two- and three-dimensional phase space. Standard numerical algorithms are replaced by visual algorithms created interactively using graphic interface tools that embody topological principles.
This approach turns out to be more flexible and powerful for both research and teaching, and contributed to two unanticipated discoveries: an index property of chaotic attractors (Japan J. Indust. Appl. Math., 8, pp. 487-504, 1991) with predictive power (Phys. Rev. E, 49, pp. 1019-1027, 1994); and bifurcation with indeterminate outcome (Proc. R. Soc. Lond. A, 432, pp. 113-123, 1991), which completes the triad of behaviors involving sensitive dependence.
Blair Swartz <bks@lanl.gov> Los Alamos National Laboratory, NM
Abstract
An iterative algorithm is proposed to solve the following form of the "ham sandwich" problem: given n bodies in n dimensions, along with a mass fraction specified for each body, find a hyperplane simultaneously separating each body into two pieces so that each body's specified fraction lies on the same side of the hyperplane. A stage of the iteration begins with a guess of the hyperplane's unit normal vector. Next, and for each body: move the hyperplane across that body in the direction of its normal until it cuts off the fraction specified for that body; and find the center of mass of the section of the body contained in this particular hyperplane. Stop iterating if the n (parallel) hyperplanes just constructed are close enough to each other. Otherwise, these n centers of mass determine a new hyperplane; and its unit normal initiates the next stage of the iteration. Quadratic convergence takes place under the assumption that the difference between the mass fractions on the same side of two hyperplanes, each of which contains the center of mass of the body's intersection with one of them, is second-order in the angle between their normals. Numerical examples in two and three dimensions are provided. The algorithm and the principles underlying it can be helpful in constructing local, usually second-order accurate, discontinuous, piecewise linear approximations of a smooth curved interface using the mass fractions it cuts off from small neighboring cells---e.g., from irregular mesh cells in computational fluid dynamics. Caution: for some shapes or configurations of cells there exist n-tuples of mass fractions that cannot be simultaneously sliced from the cells.
Eitan Tadmor <tadmor@math.ucla.edu> University of California Los Angeles, CA
Abstract
Stability implies convergence, but for nonlinear problems we can do with less. To make our point, we focus on non-oscillatory central schemes as prototype for Godunov-type projection methods. A variety of numerical experiments demonstrate that the proposed central schemes offer simple, robust, Riemann-solver-free alternatives to the the more expensive upwind schemes. Our new convergence results apply to general projection methods, including those which are based on modern nonlinear projections employed by both central as well as high-resolution schemes.
It is shown that if the approximate solution possess a weak regularity, then it converges to the unique entropy solution. The convergence result is obtained by interpolation between this weak regularity and well-posedness measured in an appropriate negative norm. Here, weak regularity is quantified in terms of a priori BV-bound, a weaker entropy-production bound, ... Most notably, since our theory does not require a stronger stability condition, one is able to apply these results to high-resolution approximations, and in this context we establish first convergence results of fully-discrete high-resolution approximations. In particular, we answer the long open question concerning the convergence of the (fully-discrete) second-order MUSCL scheme and related high-resolution central approximations. These convergence results are complemented by the corresponding error estimates.
Vidar Thomee <thomee@math.chalmers.se> Chalmers University of Technology, Goteborg, Sweden
Abstract
We consider spatially semidiscrete and fully discrete finite element approximations of the solution of an initial-boundary value problem for a parabolic equation. We survey known results concerning stability and error estimates with respect to the maximum-norm, including some recent developments. The presentation relates to properties of analytic semigroups and their approximation.