SIAM Digital Library
 
 
 

SIAM J. on Scientific Computing

All Journal content published prior to 1997 is part of LOCUS.

Search Issue | RSS Feeds RSS
Previous Issue

2006

Volume 28, Issue 6, pp. 2001-2389


A Fully Mass and Volume Conserving Implementation of a Characteristic Method for Transport Problems

Todd Arbogast and Chieh‐Sen Huang

SIAM J. Sci. Comput. 28, pp. 2001-2022 (22 pages) | Cited 6 times

Online Publication Date: December 05, 2006

Full Text: | Download PDF

Show Abstract
The characteristics‐mixed method considers the transport not of a single point or fluid particle, but rather the mass in an entire region of fluid. This mass is transported along the characteristic curves of the hyperbolic part of the transport equation, and the scheme thereby produces very little numerical dispersion, conserves mass locally, and can use long time steps. However, since the shape of a characteristic trace‐back region must be approximated in numerical implementation, its volume may be incorrect, resulting in inaccurate concentration densities and, further, inaccurate reaction dynamics. We present a simple modification to the characteristics‐mixed method that conserves both mass and volume of the transported fluid regions. Our algorithm also handles boundary conditions through a space‐time change of variables in the trace‐back routines, which allows the boundary to be treated as if it were interior to the domain. Nearly point sources, such as wells, present special difficulties, since characteristic trace‐back curves converge in their vicinity. We also present techniques that allow one to conservatively implement wells. The techniques are illustrated in five numerical examples.

Adaptive Finite Element Methods for Elliptic PDEs Based on Conforming Centroidal Voronoi–Delaunay Triangulations

Lili Ju, Max Gunzburger, and Weidong Zhao

SIAM J. Sci. Comput. 28, pp. 2023-2053 (31 pages) | Cited 8 times

Online Publication Date: December 05, 2006

Full Text: | Download PDF

Show Abstract
A new triangular mesh adaptivity algorithm for elliptic PDEs that combines a posteriori error estimation with centroidal Voronoi–Delaunay tessellations of domains in two dimensions is proposed and tested. The ability of the first ingredient to detect local regions of large error and the ability of the second ingredient to generate superior unstructured grids result in a mesh adaptivity algorithm that has several very desirable features, including the following. Errors are very well equidistributed over the triangles; at all levels of refinement, the triangles remain very well shaped, even if the grid size at any particular refinement level, when viewed globally, varies by several orders of magnitude; and the convergence rates achieved are the best obtainable using piecewise linear finite elements. This methodology can be easily extended to higher‐order finite element approximations or mixed finite element formulations although only the linear approximation is considered in this paper.

An Embedded Boundary Method for the Wave Equation with Discontinuous Coefficients

Heinz‐Otto Kreiss and N. Anders Petersson

SIAM J. Sci. Comput. 28, pp. 2054-2074 (21 pages) | Cited 4 times

Online Publication Date: December 05, 2006

Full Text: | Download PDF

Show Abstract
A second order accurate embedded boundary method for the two‐dimensional wave equation with discontinuous wave propagation speed is described. The wave equation is discretized on a Cartesian grid with constant grid size and the interface (across which the wave speed is discontinuous) is allowed to intersect the mesh in an arbitrary fashion. By using ghost points on either side of the interface, previous embedded boundary techniques for the Neumann and Dirichlet problems are generalized to satisfy the jump conditions across the interface to second order accuracy. The resulting discretization of the jump conditions has the desirable property that each ghost point can be updated independently of all other ghost points, resulting in a fully explicit time‐integration method. We prove that the one‐dimensional restriction of the method is stable without damping for arbitrary locations of the interface relative to the grid. For the two‐dimensional case, the previously developed fourth order $A^T A$‐dissipation is generalized to handle jump conditions. We demonstrate that this operator provides sufficient stabilization to enable long‐time simulations while being weak enough to preserve the accuracy of the solution. Numerical examples are given where the method is used to study electromagnetic scattering of a plane wave by a dielectric cylinder. The numerical solutions are evaluated against the analytical solution due to Mie, and pointwise second order accuracy is confirmed.

Analysis of a Fourth‐Order Scheme for a Three‐Dimensional Convection‐Diffusion Model Problem

Ashvin Gopaul and Muddun Bhuruth

SIAM J. Sci. Comput. 28, pp. 2075-2094 (20 pages) | Cited 2 times

Online Publication Date: December 05, 2006

Full Text: | Download PDF

Show Abstract
We derive closed form expressions for the eigenvalues and discrete solution arising from a 19‐point compact discretization of a three‐dimensional convection‐diffusion problem. It is shown that the coefficient matrix is positive definite when the cell‐Reynolds number is greater than some critical value. By analyzing the terms composing the discrete solution, we prove that an oscillation‐free discrete solution is guaranteed whenever the cell‐Reynolds number exceeds a value which is grid‐size dependent. An interesting result is that as the mesh size is refined, this value approaches the Golden Mean.

An Augmented Lagrangian‐Based Approach to the Oseen Problem

Michele Benzi and Maxim A. Olshanskii

SIAM J. Sci. Comput. 28, pp. 2095-2113 (19 pages) | Cited 32 times

Online Publication Date: December 05, 2006

Full Text: | Download PDF

Show Abstract
We describe an effective solver for the discrete Oseen problem based on an augmented Lagrangian formulation of the corresponding saddle point system. The proposed method is a block triangular preconditioner used with a Krylov subspace iteration like BiCGStab. The crucial ingredient is a novel multigrid approach for the (1,1) block, which extends a technique introduced by Schoberl for elasticity problems to nonsymmetric problems. Our analysis indicates that this approach results in fast convergence, independent of the mesh size and largely insensitive to the viscosity. We present experimental evidence for both isoP2‐P0 and isoP2‐P1 finite elements in support of our conclusions. We also show results of a comparison with two state‐of‐the‐art preconditioners, showing the competitiveness of our approach.

An Adaptive WEM Algorithm for Solving Elliptic Boundary Value Problems in Fairly General Domains

S. Berrone and T. Kozubek

SIAM J. Sci. Comput. 28, pp. 2114-2138 (25 pages) | Cited 2 times

Online Publication Date: December 05, 2006

Full Text: | Download PDF

Show Abstract
In this paper, we introduce a simple adaptive wavelet element algorithm similar to the Cohen–Dahmen–DeVore algorithm [A. Cohen, W. Dahmen, and R. DeVore, Math. Comp., 70 (2001), pp. 27–75]. The main difference is that we do not assume knowledge of the many constants appearing therein. The algorithm is easy to implement and applicable to a large class of problems in fairly general domains. The efficiency is illustrated by several two‐dimensional numerical examples and compared with an adaptive finite element method.

A Predictor‐Corrector Algorithm for Reaction‐Diffusion Equations Associated with Neural Activity on Branched Structures

M. J. Rempe and D. L. Chopp

SIAM J. Sci. Comput. 28, pp. 2139-2161 (23 pages) | Cited 6 times

Online Publication Date: December 05, 2006

Full Text: | Download PDF

Show Abstract
In this paper, we present a numerical method for solving reaction‐diffusion equations on one‐dimensional branched structures. The method effectively decouples the many branches at the nodal junctions so that the equations can be solved as a system of smaller problems that are purely tridiagonal. This strategy enables a number of improvements. In particular, difficulties associated with implicit methods for closed loops in the network are eliminated, and spatial adaptivity can be easily implemented. We apply this method to the electrical activity of a neuron. A neuron can be effectively represented as a one‐dimensional branched structure, and the model equations, based on the Hodgkin–Huxley cable equations, are a set of reaction equations coupled to a single diffusion process. Spatially adaptive methods are particularly effective for neuronal simulations as the region of activity is often very localized in space. Since the algorithm scales with activity instead of network size, it will make a substantial reduction in computational cost over existing methods.

Constructing Embedded Lattice Rules for Multivariate Integration

Ronald Cools, Frances Y. Kuo, and Dirk Nuyens

SIAM J. Sci. Comput. 28, pp. 2162-2188 (27 pages) | Cited 9 times

Online Publication Date: December 05, 2006

Full Text: | Download PDF

Show Abstract
Lattice rules are a family of equal‐weight cubature formulae for approximating high‐dimensional integrals. By now it is well established that good generating vectors for lattice rules having $n$ points can be constructed component‐by‐component for integrands belonging to certain weighted function spaces, and that they can achieve the optimal rate of convergence. Although the lattice rules constructed this way are extensible in dimension, they are not extensible in $n$; thus when $n$ is changed the generating vector needs to be constructed anew. In this paper we introduce a new algorithm for constructing good generating vectors for embedded lattice rules which can be used for a range of $n$ while still being extensible in dimension. By using an adaptation of the fast component‐by‐component construction algorithm (which makes use of fast Fourier transforms), we are able to obtain good generating vectors for thousands of dimensions and millions of points, under both product weight and order‐dependent weight settings, at the cost of $\calO(d n (\log (n))^2)$ operations. With a sufficiently large number of points and good overall quality, these embedded lattice rules can be used for practical purposes in the same way as a low‐discrepancy sequence. We show for a range of weight settings in an unanchored Sobolev space that our embedded lattice rules achieve the same (optimal) rate of convergence $\calO(n^{-1+\delta})$, $\delta>0$, as those constructed for a fixed number of points, and that the implied constant gains only a factor of 1.30 to 1.55.

Modeling Water Transport across Elastic Boundaries Using an Explicit Jump Method

Anita T. Layton

SIAM J. Sci. Comput. 28, pp. 2189-2207 (19 pages) | Cited 6 times

Online Publication Date: December 11, 2006

Full Text: | Download PDF

Show Abstract
A mathematical model is presented to simulate water and solute transport in a highly viscous fluid with a water‐permeable, elastic immersed membrane. In this model, fluid motion is described by Stokes flow, whereas water fluxes across the membrane are driven by transmural pressure and solute concentration differences. The elastic forces, arising from the membrane being distorted from its relaxed configuration, and the transmembrane water fluxes introduce into model solutions discontinuities across the membrane. Such discontinuities are faithfully captured using a second‐order explicit jump method [A. Mayo, SIAM J. Numer. Anal., 21 (1984), pp. 285–299], in which jumps in the solution and its derivatives are incorporated into a finite‐difference scheme. Numerical results suggest that the method exhibits desirable volume accuracy and mass conservation.

Analysis and Comparison of Geometric and Algebraic Multigrid for Convection‐Diffusion Equations

Chin‐Tien Wu and Howard C. Elman

SIAM J. Sci. Comput. 28, pp. 2208-2228 (21 pages) | Cited 2 times

Online Publication Date: December 11, 2006

Full Text: | Download PDF

Show Abstract
The discrete convection‐diffusion equations obtained from streamline diffusion finite element discretization are solved on both uniform meshes and adaptive meshes. Estimates of error reduction rates for both geometric multigrid (GMG) and algebraic multigrid (AMG) are established on uniform rectangular meshes for a model problem. Our analysis shows that GMG with line Gauss–Seidel smoothing and bilinear interpolation converges if $h\gg \epsilon^{2/3}$, and AMG with the same smoother converges more rapidly than GMG if the interpolation constant $\beta$ in the approximation assumption of AMG satisfies $\beta \ll (\frac{h}{\sqrt{\epsilon}})^{\alpha}, \hs{1mm} \mbox{where $ \alpha = \Big\{ {\begin{smallmatrix} 1, &{h < \sqrt \varepsilon, } \\ 2, &{h \ge \sqrt \varepsilon.} \\ \end{smallmatrix}}$}$ On unstructured triangular meshes, the performance of GMG and AMG, both as solvers and as preconditioners for GMRES, are evaluated. Numerical results show that GMRES with AMG preconditioning is a robust and reliable solver on both type of meshes.

Central WENO Schemes for Hamilton–Jacobi Equations on Triangular Meshes

Doron Levy, Suhas Nayak, Chi‐Wang Shu, and Yong‐Tao Zhang

SIAM J. Sci. Comput. 28, pp. 2229-2247 (19 pages) | Cited 6 times

Online Publication Date: December 11, 2006

Full Text: | Download PDF

Show Abstract
We derive Godunov‐type semidiscrete central schemes for Hamilton–Jacobi equations on triangular meshes. High‐order schemes are then obtained by combining our new numerical fluxes with high‐order WENO reconstructions on triangular meshes. The numerical fluxes are shown to be monotone in certain cases. The accuracy and high‐resolution properties of our scheme are demonstrated in a variety of numerical examples.

Semi‐Implicit Covolume Method in 3D Image Segmentation

S. Corsaro, K. Mikula, A. Sarti, and F. Sgallari

SIAM J. Sci. Comput. 28, pp. 2248-2265 (18 pages) | Cited 8 times

Online Publication Date: December 11, 2006

Full Text: | Download PDF

Show Abstract
We introduce a three‐dimensional (3D) semi‐implicit complementary volume numerical scheme for solving the level set formulation of (Riemannian) mean curvature flow problem. We apply the scheme to segmentation of objects (with interrupted edges) in 3D images. The method is unconditionally stable and efficient regarding computational times. The study of its experimental order of convergence on 3D examples shows its second order accuracy for smooth solutions and first order accuracy for highly singular solutions with vanishing gradients as arising in image segmentation.

A Parallel Multistage ILU Factorization Based on a Hierarchical Graph Decomposition

Pascal Hénon and Yousef Saad

SIAM J. Sci. Comput. 28, pp. 2266-2293 (28 pages) | Cited 7 times

Online Publication Date: December 18, 2006

Full Text: | Download PDF

Show Abstract
PHIDAL (parallel hierarchical interface decomposition algorithm) is a parallel incomplete factorization method which exploits a hierarchical interface decomposition of the adjacency graph of the coefficient matrix. The idea of the decomposition is similar to that of the well‐known wirebasket techniques used in domain decomposition. However, the method is devised for general, irregularly structured, sparse linear systems. This paper describes a few algorithms for obtaining good quality hierarchical graph decompositions and discusses the parallel implementation of the factorization procedure. Numerical experiments are reported to illustrate the scalability of the algorithm and its effectiveness as a general purpose parallel linear system solver.

2–D Parachute Simulation by the Immersed Boundary Method

Yongsam Kim and Charles S. Peskin

SIAM J. Sci. Comput. 28, pp. 2294-2312 (19 pages) | Cited 12 times

Online Publication Date: December 21, 2006

Full Text: | Download PDF

Show Abstract
Parachute aerodynamics involves an interaction between the flexible, elastic, porous parachute canopy and the high speed airflow (relative to the parachute) through which the parachute falls. Computer simulation of parachute dynamics typically simplifies the problem in various ways, e.g., by considering the parachute as a rigid bluff body. Here, we avoid such simplification by using the immersed boundary (IB) method to study the full fluid‐structure interaction. The IB method is generalized to handle porous immersed boundaries, and the generalized method is used to study the influence of porosity on parachute stability.

Discrete Maximum Principle and Adequate Discretizations of Linear Parabolic Problems

István Faragó and Róbert Horváth

SIAM J. Sci. Comput. 28, pp. 2313-2336 (24 pages) | Cited 2 times

Online Publication Date: December 21, 2006

Full Text: | Download PDF

Show Abstract
In this paper, we analyze the connections between the different qualitative properties of numerical solutions of linear parabolic problems with Dirichlet‐type boundary condition. First we formulate the qualitative properties for the differential equations and shed light on their relations. Then we show how the well‐known discretization schemes can be written in the form of a one‐step iterative process. We give necessary and sufficient conditions of the main qualitative properties of these iterations. We apply the results to the finite difference and Galerkin finite element solutions of linear parabolic problems. In our main result we show that the nonnegativity preservation property is equivalent to the maximum‐minimum principle and they imply the maximum norm contractivity. In one, two, and three dimensions, we list sufficient a priori conditions that ensure the required qualitative properties. Finally, we demonstrate the above results on numerical examples.

Iterative Validation of Eigensolvers: A Scheme for Improving the Reliability of Hermitian Eigenvalue Solvers

James R. McCombs and Andreas Stathopoulos

SIAM J. Sci. Comput. 28, pp. 2337-2358 (22 pages) | Cited 2 times

Online Publication Date: December 21, 2006

Full Text: | Download PDF

Show Abstract
Iterative eigenvalue solvers for large, sparse matrices may miss some of the required eigenvalues that are of high algebraic multiplicity or tightly clustered. Block methods, locking, a posteriori validation, or simply increasing the required accuracy are often used to avoid missing or to detect a missed eigenvalue, but each has its own shortcomings in robustness or performance. To resolve these shortcomings, we have developed a postprocessing algorithm, iterative validation of eigensolvers (IVE), that combines the advantages of each technique. IVE detects numerically multiple eigenvalues among the approximate eigenvalues returned by a given solver, adjusts the block size accordingly, then calls the given solver using locking to compute a new approximation in the subspace orthogonal to the current approximate eigenvectors. This process is repeated until no additional missed eigenvalues can be identified. IVE is general and can be applied as a wrapper to any Rayleigh–Ritz‐based, Hermitian eigensolver. Our experiments show that IVE is very effective in computing missed eigenvalues even with eigensolvers that lack locking or block capabilities, although such capabilities may further enhance robustness. By focusing on robustness in a postprocessing stage, IVE allows the user to decouple the notion of robustness from that of performance when choosing the block size or the convergence tolerance.

MultiScale Modeling of Physical Phenomena: Adaptive Control of Models

J. Tinsley Oden, Serge Prudhomme, Albert Romkes, and Paul T. Bauman

SIAM J. Sci. Comput. 28, pp. 2359-2389 (31 pages) | Cited 25 times

Online Publication Date: December 26, 2006

Full Text: | Download PDF

Show Abstract
It is common knowledge that the accuracy with which computer simulations can depict physical events depends strongly on the choice of the mathematical model of the events. Perhaps less appreciated is the notion that the error due to modeling can be defined, estimated, and used adaptively to control modeling error, provided one accepts the existence of a base model that can serve as a datum with respect to which other models can be compared. In this work, it is shown that the idea of comparing models and controlling model error can be used to develop a general approach for multiscale modeling, a subject of growing importance in computational science. A posteriori estimates of modeling error in so‐called quantities of interest are derived and a class of adaptive modeling algorithms is presented. Several applications of the theory and methodology are presented. These include preliminary work on random multiphase composite materials, molecular statics simulations with applications to problems in nanoindentation, and analysis of molecular dynamics models using various techniques for scale bridging.
Close

Your Recent History Recent History Help

Recently Viewed

  1. Multigrid for High-Dimensional Elliptic Partial Differential Equations on Non-equidistant Grids
  2. Additive and Multiplicative Two-Level Spectral Preconditioning for General Linear Systems
  3. Semi‐Implicit Covolume Method in 3D Image Segmentation
  4. Central WENO Schemes for Hamilton–Jacobi Equations on Triangular Meshes
  5. Constructing Embedded Lattice Rules for Multivariate Integration
© SIAM By using SIAM Publications Online you agree to abide by the Terms and Conditions of Use. Banner art adapted from a figure by Hinke M. Osinga and Bernd Krauskopf (University of Bristol, UK).
close