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

2009

Volume 31, Issue 6, pp. 4013-4813


Strengthened Linear Sampling Method with a Reference Ball

Jingzhi Li, Hongyu Liu, and Jun Zou

SIAM J. Sci. Comput. 31, pp. 4013-4040 (28 pages)

Online Publication Date: November 13, 2009

Full Text: | Download PDF

Show Abstract
By adding a reference ball as an extra artificial obstacle component to the underlying scattering system, we propose a simple but robust and effective technique to choose the crucial cut-off value required in the linear sampling method (LSM) for inverse scattering problems. The reference ball technique causes few extra computational costs to the LSM, but brings in a practically important byproduct, i.e., it eliminates the interior eigenvalue problem automatically, a well-known barrier when applying the LSM. Some mathematical justifications of the new technique are provided, and numerical experiments are also presented to demonstrate its feasibility and effectiveness.

Robin Based Semi-Implicit Coupling in Fluid-Structure Interaction: Stability Analysis and Numerics

Matteo Astorino, Franz Chouly, and Miguel A. Fernández

SIAM J. Sci. Comput. 31, pp. 4041-4065 (25 pages) | Cited 2 times

Online Publication Date: November 13, 2009

Full Text: | Download PDF

Show Abstract
In this paper, we propose a semi-implicit coupling scheme for the numerical simulation of fluid-structure interaction systems involving a viscous incompressible fluid. The scheme is stable irrespective of the so-called added-mass effect and allows for conservative time-stepping within the structure. The efficiency of the scheme is based on the explicit splitting of the viscous effects and geometrical/convective nonlinearities through the use of the Chorin–Temam projection scheme within the fluid. Stability comes from the implicit pressure-solid coupling and a specific Robin treatment of the explicit viscous-solid coupling, derived from Nitsche's method.

A Finite Volume Method for Solving Parabolic Equations on Logically Cartesian Curved Surface Meshes

Donna A. Calhoun and Christiane Helzel

SIAM J. Sci. Comput. 31, pp. 4066-4099 (34 pages)

Online Publication Date: November 13, 2009

Full Text: | Download PDF

Show Abstract
We present a second order finite volume scheme for the constant-coefficient diffusion equation on curved parametric surfaces. While our scheme is applicable to general quadrilateral surface meshes based on smooth or piecewise smooth coordinate transformations, our primary motivation for developing the present scheme is to solve diffusion problems on a particular set of circular and spherical meshes introduced in [D. A. Calhoun, C. Helzel, and R. J. LeVeque, SIAM Rev., 50 (2008), pp. 723–752] for the discretization of hyperbolic problems. These grids are generated from mappings of a single Cartesian grid and were designed to have nearly uniform cells sizes and avoid the pole singularity associated with polar or spherical grid mappings. The present method for parabolic equations offers several advantages. It does not require analytic metric terms, shows second order accuracy on our disk and sphere grids, can be easily coupled to existing finite volume solvers for logically Cartesian meshes, and handles general mixed boundary conditions. Our parabolic scheme should appeal to researchers in the fields of geophysical fluid dynamics, computational biology, and any other discipline that requires the solution of parabolic equations on quadrilateral surface meshes. In this article, we present several numerical examples demonstrating the accuracy of the scheme, and then use the scheme to solve advection-reaction-diffusion equations modeling biological pattern formation on surfaces.

Optimal Partitions for Eigenvalues

Blaise Bourdin, Dorin Bucur, and Édouard Oudet

SIAM J. Sci. Comput. 31, pp. 4100-4114 (15 pages)

Online Publication Date: November 13, 2009

Full Text: | Download PDF

Show Abstract
We introduce a new numerical method for approximating partitions of a domain minimizing the sum of Dirichlet–Laplacian eigenvalues of any order. First we prove the equivalence of the original problem and a relaxed formulation based on measures. Using this result, we build a numerical algorithm to approximate optimal configurations. We describe numerical experiments aimed at studying the asymptotic behavior of optimal partitions with large numbers of cells.

Extrapolation Algorithms for Solving Mixed Boundary Integral Equations of the Helmholtz Equation by Mechanical Quadrature Methods

Jin Huang and Zhu Wang

SIAM J. Sci. Comput. 31, pp. 4115-4129 (15 pages)

Online Publication Date: November 20, 2009

Full Text: | Download PDF

Show Abstract
This paper presents mechanical quadrature methods with high accuracy for solving mixed boundary integral equations of the Helmholtz equation. By estimating the range of eigenvalues for the discretization matrix of the integral equations and applying the collectively compact convergent theory, we prove the stability and convergence of numerical solutions, which is a challenging task for this method. Moreover, the asymptotic error expansions show the method is of order $h^3$. Hence, extrapolation algorithms can be introduced to achieve higher approximation accuracy degree $(\mathcal{O}(h^5))$. Meanwhile, an a posteriori asymptotic error estimate is derived, which can be used to construct self-adaptive algorithms. The numerical examples support our theoretical analysis.

Efficient Assembly of $H(\mathrm{div})$ and $H(\mathrm{curl})$ Conforming Finite Elements

Marie E. Rognes, Robert C. Kirby, and Anders Logg

SIAM J. Sci. Comput. 31, pp. 4130-4151 (22 pages)

Online Publication Date: November 20, 2009

Full Text: | Download PDF

Show Abstract
In this paper, we discuss how to efficiently evaluate and assemble general finite element variational forms on $H(\mathrm{div})$ and $H(\mathrm{curl})$. The proposed strategy relies on a decomposition of the element tensor into a precomputable reference tensor and a mesh-dependent geometry tensor. Two key points must then be considered: the appropriate mapping of basis functions from a reference element, and the orientation of geometrical entities. To address these issues, we extend here a previously presented representation theorem for affinely mapped elements to Piola-mapped elements. We also discuss a simple numbering strategy that removes the need to contend with directions of facet normals and tangents. The result is an automated, efficient, and easy-to-use implementation that allows a user to specify finite element variational forms on $H(\mathrm{div})$ and $H(\mathrm{curl})$ in close to mathematical notation.

Multilevel Algorithms for Large-Scale Interior Point Methods

Michele Benzi, Eldad Haber, and Lauren Taralli

SIAM J. Sci. Comput. 31, pp. 4152-4175 (24 pages) | Cited 1 time

Online Publication Date: November 20, 2009

Full Text: | Download PDF

Show Abstract
We develop and compare multilevel algorithms for solving constrained nonlinear variational problems via interior point methods. Several equivalent formulations of the linear systems arising at each iteration of the interior point method are compared from the point of view of conditioning and iterative solution. Furthermore, we show how a multilevel continuation strategy can be used to obtain good initial guesses (“hot starts”) for each nonlinear iteration. Some minimal surface and volume-constrained image registration problems are used to illustrate the various approaches.

A Globalized Newton Method for the Accurate Solution of a Dipole Quantum Control Problem

G. von Winckel, A. Borzì, and S. Volkwein

SIAM J. Sci. Comput. 31, pp. 4176-4203 (28 pages)

Online Publication Date: November 20, 2009

Full Text: | Download PDF

Show Abstract
A theoretical and computational framework is presented to obtain accurate controls for fast quantum state transitions that are needed in a host of applications such as nanoelectronic devices and quantum computing. This method is based on a reduced Hessian Krylov–Newton scheme applied to a norm-preserving discrete model of a dipole quantum control problem. The use of second-order numerical methods for solving the control problem is justified, proving the existence of optimal solutions and analyzing first- and second-order optimality conditions. Criteria for the discretization of the nonconvex optimization problem and for the formulation of the Hessian are given to ensure accurate gradients and a symmetric Hessian. Robustness of the Newton approach is obtained using a globalization strategy with a robust linesearch procedure. Results of numerical experiments demonstrate that the Newton approach presented in this paper is able to provide fast and accurate controls for high-energy state transitions.

GPU-Based Volume Reconstruction from Very Few Arbitrarily Aligned X-Ray Images

Daniel Gross, Ulrich Heil, Ralf Schulze, Elmar Schoemer, and Ulrich Schwanecke

SIAM J. Sci. Comput. 31, pp. 4204-4221 (18 pages)

Online Publication Date: November 25, 2009

Full Text: | Download PDF

Show Abstract
This paper presents a three-dimensional GPU-accelerated algebraic reconstruction method in a few-projection cone-beam setting with arbitrary acquisition geometry. To achieve artifact-reduced reconstructions in the challenging case of unconstrained geometry and extremely limited input data, we use linear methods and an artifact-avoiding projection algorithm to provide high reconstruction quality. We apply the conjugate gradient method in the linear case of Tikhonov regularization and the two-point-step-size gradient method in the nonlinear case of total variation regularization to solve the system of equations. By taking advantage of modern graphics hardware we achieve acceleration of up to two orders of magnitude over classical CPU implementations.

High-Order Narrow Stencil Finite-Difference Approximations of Second-Order Derivatives Involving Variable Coefficients

Ramji Kamakoti and Carlos Pantano

SIAM J. Sci. Comput. 31, pp. 4222-4243 (22 pages) | Cited 1 time

Online Publication Date: November 25, 2009

Full Text: | Download PDF

Show Abstract
High-order-accuracy finite-difference approximations are developed for problems involving arbitrary variable coefficients in the second-order derivatives, e.g., the heat equation or turbulence modeling. The methods investigated are discretely conservative, use narrow stencils, and provide stable approximations for these problems. It is known that high-order finite-difference approximations for these types of equations using the chain rule approach may be inadequate for approximating partial differential equations with certain types of variable coefficients. The new approximations are constructed to alleviate this problem by requiring that the operators are stable when the variable coefficients are positive. Examples in heat-transfer problems with variable coefficients are shown to retain the designed order of accuracy and stability with lower error norms than the usual alternative discretizations. Finally, an application of the new stencils is presented for the large-eddy simulation of compressible turbulent flows.

A Finite Element Splitting Extrapolation for Second Order Hyperbolic Equations

Xiaoming He and Tao Lü

SIAM J. Sci. Comput. 31, pp. 4244-4265 (22 pages)

Online Publication Date: November 25, 2009

Full Text: | Download PDF

Show Abstract
Splitting extrapolation is an efficient technique for solving large scale scientific and engineering problems in parallel. This article discusses a finite element splitting extrapolation for second order hyperbolic equations with time-dependent coefficients. This method possesses a higher degree of parallelism, less computational complexity, and more flexibility than Richardson extrapolation while achieving the same accuracy. By means of domain decomposition and isoparametric mapping, some grid parameters are chosen according to the problem. The multiparameter asymptotic expansion of the $d$-quadratic finite element error is also established. The splitting extrapolation formulas are developed from this expansion. An approximation with higher accuracy on a globally fine grid can be computed by solving a set of smaller discrete subproblems on different coarser grids in parallel. Some a posteriori error estimates are also provided. Numerical examples show that this method is efficient for solving discontinuous problems and nonlinear hyperbolic equations.

Continuous Block $\theta$-Methods for Ordinary and Delay Differential Equations

Hongjiong Tian, Kaiting Shan, and Jiaoxun Kuang

SIAM J. Sci. Comput. 31, pp. 4266-4280 (15 pages)

Online Publication Date: November 25, 2009

Full Text: | Download PDF

Show Abstract
Continuous numerical methods have many applications in the numerical solution of discontinuous ordinary differential equations (ODEs), delay differential equations, neutral differential equations, integro-differential equations, etc. This paper deals with a continuous extension for the discrete approximate solution of ODEs generated by a class of block $\theta$-methods. Existence and uniqueness for the continuous extension are discussed. Convergence and absolute stability of the continuous block $\theta$-methods for ODEs are studied. As an application, we adopt the continuous block $\theta$-methods to solve delay differential equations and prove that the continuous block $\theta$-methods are $GP$-stable if and only if they are $A_{\omega}$-stable for ODEs. Several numerical experiments are given to illustrate the performance of the continuous block $\theta$-methods.

Sparse Tensor Discretization of Elliptic sPDEs

Marcel Bieri, Roman Andreev, and Christoph Schwab

SIAM J. Sci. Comput. 31, pp. 4281-4304 (24 pages) | Cited 4 times

Online Publication Date: December 04, 2009

Full Text: | Download PDF

Show Abstract
We propose and analyze sparse deterministic-stochastic tensor Galerkin finite element methods (sparse sGFEMs) for the numerical solution of elliptic partial differential equations (PDEs) with random coefficients in a physical domain $D\subset\mathbb{R}^d$. In tensor product sGFEMs, the variational solution to the boundary value problem is approximated in tensor product finite element spaces $V^\Gamma\otimes V^D$, where $V^\Gamma$ and $V^D$ denote suitable finite dimensional subspaces of the stochastic and deterministic function spaces, respectively. These approaches lead to sGFEM algorithms of complexity $O(N_\Gamma N_D)$, where $N_\Gamma=\dim V^\Gamma$ and $N_D=\dim V^D$. In this work, we use hierarchic sequences $V^\Gamma_1\subset V^\Gamma_2\subset\ldots$ and $V^D_1\subset V^D_2\subset\ldots$ of finite dimensional spaces to approximate the law of the random solution. The hierarchies of approximation spaces allow us to define sparse tensor product spaces $V^\Gamma_\ell\hat{\otimes}V^D_\ell$, $\ell=1,2,\dots$, yielding algorithms of $O(N_\Gamma\log N_D+N_D\log N_\Gamma)$ work and memory. We estimate the convergence rate of sGFEM for algebraic decay of the input random field Karhunen–Loève coefficients. We give an algorithm for an input adapted a-priori selection of deterministic and stochastic discretization spaces. The convergence rate in terms of the total number of degrees of freedom of the proposed method is superior to Monte Carlo approximations. Numerical examples illustrate the theoretical results and demonstrate superiority of the sparse tensor product discretization proposed here versus the full tensor product approach.

Numerical Simulation of Diffusive and Aggregation Phenomena in Nonlinear Continuity Equations by Evolving Diffeomorphisms

J. A. Carrillo and J. S. Moll

SIAM J. Sci. Comput. 31, pp. 4305-4329 (25 pages) | Cited 1 time

Online Publication Date: December 04, 2009

Full Text: | Download PDF

Show Abstract
We propose a numerical algorithm for solving nonlinear continuity equations written in Lagrangian coordinates. This transformation is intimately related to variational approaches for the well-posedness of gradient flows of energy functionals with respect to the quadratic transportation distance in optimal transport theory. These schemes allow the numerical approximation of both diffusive and aggregation regimes of different models. Positivity, energy decreasing, and mesh adaptation are built-in in the numerical scheme, and thus we are capable of capturing blow-up densities and of dealing with vacuum regions and merging of mass patches in a natural way.

The Implicit Closest Point Method for the Numerical Solution of Partial Differential Equations on Surfaces

Colin B. Macdonald and Steven J. Ruuth

SIAM J. Sci. Comput. 31, pp. 4330-4350 (21 pages)

Online Publication Date: December 16, 2009

Full Text: | Download PDF

Show Abstract
Many applications in the natural and applied sciences require the solutions of partial differential equations (PDEs) on surfaces or more general manifolds. The closest point method is a simple and accurate embedding method for numerically approximating PDEs on rather general smooth surfaces. However, the original formulation is designed to use explicit time stepping. This may lead to a strict time-step restriction for some important PDEs such as those involving the Laplace–Beltrami operator or higher-order derivative operators. To achieve improved stability and efficiency, we introduce a new implicit closest point method for surface PDEs. The method allows for large, stable time steps while retaining the principal benefits of the original method. In particular, it maintains the order of accuracy of the discretization of the underlying embedding PDE, it works on sharply defined bands without degrading the accuracy of the method, and it applies to general smooth surfaces. It also is very simple and may be applied to a rather general class of surface PDEs. Convergence studies for the in-surface heat equation and a fourth-order biharmonic problem are given to illustrate the accuracy of the method. We demonstrate the flexibility and generality of the method by treating flows involving diffusion, reaction-diffusion, and fourth-order spatial derivatives on a variety of interesting surfaces including surfaces of mixed codimension.

An SDP-Based Divide-and-Conquer Algorithm for Large-Scale Noisy Anchor-Free Graph Realization

Ngai-Hang Z. Leung and Kim-Chuan Toh

SIAM J. Sci. Comput. 31, pp. 4351-4372 (22 pages)

Online Publication Date: December 16, 2009

Full Text: | Download PDF

Show Abstract
We propose the DISCO algorithm for graph realization in $\mathbb{R}^d$, given sparse and noisy short-range intervertex distances as inputs. Our divide-and-conquer algorithm works as follows. When a group has a sufficiently small number of vertices, the basis step is to form a graph realization by solving a semidefinite program. The recursive step is to break a large group of vertices into two smaller groups with overlapping vertices. These two groups are solved recursively, and the subconfigurations are stitched together, using the overlapping atoms, to form a configuration for the larger group. At intermediate stages, the configurations are improved by gradient descent refinement. The algorithm is applied to the problem of determining protein moleculer structure. Tests are performed on molecules taken from the Protein Data Bank database. For each molecule, given 20–30% of the inter-atom distances less than 6Å that are corrupted by a high level of noise, DISCO is able to reliably and efficiently reconstruct the conformation of large molecules. In particular, given 30% of distances with 20% multiplicative noise, a 13000-atom conformation problem is solved within an hour with a root mean square deviation of 1.6Å.

An Adaptive Wavelet Method for the Chemical Master Equation

Tobias Jahnke

SIAM J. Sci. Comput. 31, pp. 4373-4394 (22 pages) | Cited 1 time

Online Publication Date: January 08, 2010

Full Text: | Download PDF

Show Abstract
An adaptive wavelet method for the chemical master equation is constructed. The method is based on the representation of the solution in a sparse Haar wavelet basis, the time integration by Rothe's method, and an iterative procedure which in each time-step selects those degrees of freedom which are essential for propagating the solution. The accuracy and efficiency of the approach is discussed, and the performance of the adaptive wavelet method is demonstrated by numerical examples.

Spectral Representation and Reduced Order Modeling of the Dynamics of Stochastic Reaction Networks via Adaptive Data Partitioning

Khachik Sargsyan, Bert Debusschere, Habib Najm, and Olivier Le Maître

SIAM J. Sci. Comput. 31, pp. 4395-4421 (27 pages) | Cited 1 time

Online Publication Date: January 08, 2010

Full Text: | Download PDF

Show Abstract
Dynamical analysis tools are well established for deterministic models. However, for many biochemical phenomena in cells the molecule count is low, leading to stochastic behavior that causes deterministic macroscale reaction models to fail. The main mathematical framework representing these phenomena is based on jump Markov processes that model the underlying stochastic reaction network. Conventional dynamical analysis tools do not readily generalize to the stochastic setting due to nondifferentiability and absence of explicit state evolution equations. We developed a reduced order methodology for dynamical analysis that relies on the Karhunen–Loève decomposition and polynomial chaos expansions. The methodology relies on adaptive data partitioning to obtain an accurate representation of the stochastic process, especially in the case of multimodal behavior. As a result, a mixture model is obtained that represents the reduced order dynamics of the system. The Schlögl model is used as a prototype bistable process that exhibits time scale separation and leads to multimodality in the reduced order model.

Finite-Element Preconditioning of G-NI Spectral Methods

Claudio Canuto, Paola Gervasio, and Alfio Quarteroni

SIAM J. Sci. Comput. 31, pp. 4422-4451 (30 pages)

Online Publication Date: January 15, 2010

Full Text: | Download PDF

Show Abstract
Several old and new finite-element preconditioners for nodal-based spectral discretizations of $-\Delta u=f$ in the domain $\Omega=(-1,1)^d$ ($d=2$ or 3), with Dirichlet or Neumann boundary conditions, are considered and compared in terms of both condition number and computational efficiency. The computational domain covers the case of classical single-domain spectral approximations (see [C. Canuto et al., Spectral Methods. Fundamentals in Single Domains, Springer, Heidelberg, 2006]), as well as that of more general spectral-element methods in which the preconditioners are expressed in terms of local (upon every element) algebraic solvers. The primal spectral approximation is based on the Galerkin approach with numerical integration (G-NI) at the Legendre–Gauss–Lobatto (LGL) nodes in the domain. The preconditioning matrices rely on either $\mathbb{P}_1$, $\mathbb{Q}_1$, or $\mathbb{Q}_{1,NI}$ (i.e., with numerical integration) finite elements on meshes whose vertices coincide with the LGL nodes used for the spectral approximation. The analysis highlights certain preconditioners, which yield the solution at an overall cost proportional to $N^{d+1}$, where $N$ denotes the polynomial degree in each direction.

Extrapolated Implicit-Explicit Time Stepping

Emil M. Constantinescu and Adrian Sandu

SIAM J. Sci. Comput. 31, pp. 4452-4477 (26 pages)

Online Publication Date: January 15, 2010

Full Text: | Download PDF

Show Abstract
This paper constructs extrapolated implicit-explicit time stepping methods that allow one to efficiently solve problems with both stiff and nonstiff components. The proposed methods are based on Euler steps and can provide very high order discretizations of ODEs, index-1 DAEs, and PDEs in the method-of-lines framework. Implicit-explicit schemes based on extrapolation are simple to construct, easy to implement, and straightforward to parallelize. This work establishes the existence of perturbed asymptotic expansions of global errors, explains the convergence orders of these methods, and studies their linear stability properties. Numerical results with stiff ODE, DAE, and PDE test problems confirm the theoretical findings and illustrate the potential of these methods to solve multiphysics multiscale problems.

Computational Solution of Blow-Up Problems for Semilinear Parabolic PDEs on Unbounded Domains

Hermann Brunner, Xiaonan Wu, and Jiwei Zhang

SIAM J. Sci. Comput. 31, pp. 4478-4496 (19 pages)

Online Publication Date: January 20, 2010

Full Text: | Download PDF

Show Abstract
This paper is concerned with the numerical solution of semilinear parabolic PDEs on unbounded spatial domains whose solutions blow up in finite time. The focus of the presentation is on the derivation of the nonlinear absorbing boundary conditions for one-dimensional and two-dimensional computational domains and on a simple but efficient adaptive time-stepping scheme. The theoretical results are illustrated by a broad range of numerical examples, including problems with multiple blow-up points.

Well-Centered Triangulation

Evan VanderZee, Anil N. Hirani, Damrong Guoy, and Edgar A. Ramos

SIAM J. Sci. Comput. 31, pp. 4497-4523 (27 pages)

Online Publication Date: January 20, 2010

Full Text: | Download PDF

Show Abstract
Meshes composed of well-centered simplices have nice orthogonal dual meshes (the dual Voronoi diagram). This is useful for certain numerical algorithms that prefer such primal-dual mesh pairs. We prove that well-centered meshes also have optimality properties and relationships to Delaunay and minmax angle triangulations. We present an iterative algorithm that seeks to transform a given triangulation in two or three dimensions into a well-centered one by minimizing a cost function and moving the interior vertices while keeping the mesh connectivity and boundary vertices fixed. The cost function is a direct result of a new characterization of well-centeredness in arbitrary dimensions that we present. Ours is the first optimization-based heuristic for well-centeredness and the first one that applies in both two and three dimensions. We show the results of applying our algorithm to small and large two-dimensional meshes, some with a complex boundary, and obtain a well-centered tetrahedralization of the cube. We also show numerical evidence that our algorithm preserves gradation and that it improves the maximum and minimum angles of acute triangulations created by the best known previous method.

Fast Molecular Solvation Energetics and Forces Computation

Chandrajit Bajaj and Wenqi Zhao

SIAM J. Sci. Comput. 31, pp. 4524-4552 (29 pages) | Cited 2 times

Online Publication Date: January 20, 2010

Full Text: | Download PDF

Show Abstract
The total free energy of a molecule includes the classical molecular mechanical energy (which is understood as the free energy in vacuum) and the solvation energy, which is caused by the change of the environment of the molecule (solute) from vacuum to solvent. The solvation energy is important to the study of the intermolecular interactions. In this paper we develop a fast surface-based generalized Born method to compute the electrostatic solvation energy along with the energy derivatives for the solvation forces. The most time-consuming computation is the evaluation of the surface integrals over an algebraic spline molecular surface (ASMS), and the fast computation is achieved by the use of the nonequispaced fast Fourier transform (NFFT) algorithm. The main results of this paper involve (a) an efficient sampling of quadrature points over the molecular surface by using nonlinear patches, (b) fast linear time estimation of energy and intermolecular forces, (c) error analysis, and (d) efficient implementation combining fast pairwise summation and the continuum integration using nonlinear patches.

Simultaneously Sparse Solutions to Linear Inverse Problems with Multiple System Matrices and a Single Observation Vector

Adam C. Zelinski, Vivek K. Goyal, and Elfar Adalsteinsson

SIAM J. Sci. Comput. 31, pp. 4533-4579 (47 pages)

Online Publication Date: January 20, 2010

Full Text: | Download PDF

Show Abstract
A problem that arises in slice-selective magnetic resonance imaging (MRI) radio-frequency (RF) excitation pulse design is abstracted as a novel linear inverse problem with a simultaneous sparsity constraint. Multiple unknown signal vectors are to be determined, where each passes through a different system matrix and the results are added to yield a single observation vector. Given the matrices and lone observation, the objective is to find a simultaneously sparse set of unknown vectors that approximately solves the system. We refer to this as the multiple-system single-output (MSSO) simultaneous sparse approximation problem. This manuscript contrasts the MSSO problem with other simultaneous sparsity problems and conducts an initial exploration of algorithms with which to solve it. Greedy algorithms and techniques based on convex relaxation are derived and compared empirically. Experiments involve sparsity pattern recovery in noiseless and noisy settings and MRI RF pulse design.

A Micro-Macro Decomposition-Based Asymptotic-Preserving Scheme for the Multispecies Boltzmann Equation

Shi Jin and Yingzhe Shi

SIAM J. Sci. Comput. 31, pp. 4580-4606 (27 pages)

Online Publication Date: January 27, 2010

Full Text: | Download PDF

Show Abstract
In this paper we extend the micro-macro decomposition-based asymptotic-preserving scheme developed by Bennoune, Lemou, and Mieussens [J. Comput. Phys., 227 (2008), pp. 3781–3803] for the single species Boltzmann equation to the multispecies problems. An asymptotic-preserving scheme for the kinetic equation is very efficient in the fluid regime where the Knudsen number is small and the collision term becomes stiff. It allows a coarse (independent of the Knudsen number) mesh size and a large time step in the fluid regime. The difficulty associated with multispecies problems is that there are no local conservation laws for each species, resulting in extra stiff nonlinear source terms that need to be discretized properly in order to (1) avoid Newton-type solvers for nonlinear algebraic systems, and (2) to be asymptotic-preserving. We show that these extra nonlinear source terms can be solved using only linear system solvers, and the scheme preserves the correct Euler and Navier–Stokes limits. Numerical examples are used to demonstrate the efficiency and applicability of the schemes for both Euler and Navier–Stokes regimes.

Adaptive Discontinuous Galerkin Methods for Eigenvalue Problems Arising in Incompressible Fluid Flows

K. Andrew Cliffe, Edward J. C. Hall, and Paul Houston

SIAM J. Sci. Comput. 31, pp. 4607-4632 (26 pages)

Online Publication Date: January 27, 2010

Full Text: | Download PDF

Show Abstract
In this article we consider the a posteriori error estimation and adaptive mesh refinement of discontinuous Galerkin finite element approximations of the hydrodynamic stability problem associated with the incompressible Navier–Stokes equations. Particular attention is given to the reliable error estimation of the eigenvalue problem in channel and pipe geometries. Here, computable a posteriori error bounds are derived based on employing the generalization of the standard dual-weighted-residual approach, originally developed for the estimation of target functionals of the solution, to eigenvalue/stability problems. The underlying analysis consists of constructing both a dual eigenvalue problem and a dual problem for the original base solution. In this way, errors stemming from both the numerical approximation of the original nonlinear flow problem and the underlying linear eigenvalue problem are correctly controlled. Numerical experiments highlighting the practical performance of the proposed a posteriori error indicator on adaptively refined computational meshes are presented.

Toward a Mathematical Analysis for Drift-Flux Multiphase Flow Models in Networks

Mapundi K. Banda, Michael Herty, and Jean-Medard T. Ngnotchouye

SIAM J. Sci. Comput. 31, pp. 4633-4653 (21 pages)

Online Publication Date: January 27, 2010

Full Text: | Download PDF

Show Abstract
Dynamics of multiphase flows through networks are considered. The dynamics of flow through the connected arcs are governed by an isothermal no-slip drift-flux model. Such problems arise in the context of multicomponent flows or in gas transport in pipe networks in which a phase change takes place due to geometrical or physical forces. Coupling conditions for the vertices (joints) in a network have been proposed. We present conditions at and introduce a mathematical representation of the vertex flow for the no-slip drift-flux case of multiphase flows. Mathematical analysis of coupling conditions at the vertices as well as numerical simulations and comparative studies with theoretical predictions are undertaken.

An Efficient Iterative Approach for Large-Scale Separable Nonlinear Inverse Problems

Julianne Chung and James G. Nagy

SIAM J. Sci. Comput. 31, pp. 4654-4674 (21 pages)

Online Publication Date: January 27, 2010

Full Text: | Download PDF

Show Abstract
We present an efficient iterative approach to solving separable nonlinear least squares problems that arise in large-scale inverse problems. A variable projection Gauss–Newton method is used to solve the nonlinear least squares problem, and Tikhonov regularization is incorporated using an iterative hybrid scheme. Regularization parameters are chosen automatically using a weighted generalized cross validation method, thus providing a nonlinear solver that requires very little input from the user. Applications from image deblurring and digital tomosynthesis illustrate the effectiveness of the resulting numerical scheme.

Continuous and Discrete Composite Adjoints for the Hessian of the Lagrangian in Shooting Algorithms for Dynamic Optimization

Ralf Hannemann and Wolfgang Marquardt

SIAM J. Sci. Comput. 31, pp. 4675-4695 (21 pages)

Online Publication Date: January 27, 2010

Full Text: | Download PDF

Show Abstract
One approach to solve optimal control problems by direct methods is the so-called sequential approach or single shooting. Only the control variables are discretized resulting in a nonlinear program (NLP) which can be solved with SQP or interior point methods. This paper presents a new methodology to efficiently provide the Hessian of the Lagrangian of that resulting NLP. The algorithm is based on the second-order adjoint method and introduces the novel concept of composite adjoints to reduce the computational effort of a Hessian evaluation. Though this contribution is for the sake of simplicity restricted to single shooting, the same methodology can also be easily applied to multiple shooting.

A Fast Time Stepping Method for Evaluating Fractional Integrals

Jing-Rebecca Li

SIAM J. Sci. Comput. 31, pp. 4696-4714 (19 pages)

Online Publication Date: February 03, 2010

Full Text: | Download PDF

Show Abstract
We evaluate the fractional integral $I^\alpha[f](t)=\frac{1}{\Gamma(\alpha)}\int_0^t(t-\tau)^{\alpha-1}\,f(\tau)\,d\tau$, $0<\alpha<1$, at time steps $t=\Delta t,2\Delta t,\dots,N\Delta t$ by making use of the integral representation of the convolution kernel $t^{\alpha-1}=\frac{1}{\Gamma(1-\alpha)}\int_0^{\infty}e^{-\xi\,t}\,\xi^{-\alpha}\,d\xi$. We construct an efficient $Q$-point quadrature of this integral representation and use it as a part of a fast time stepping method. The new method has algorithmic complexity $O(NQ)$ and storage requirement $O(Q)$. The number of quadrature nodes $Q$ is independent of $N$ and grows like $O\bigl(\bigl(-\log\epsilon-\log\Delta t\bigr)^2\bigr)$, where $\epsilon$ is the quadrature error tolerance and $\Delta t$ is the size of the time step. The (possible) singularity of $f$ near $\tau=0$ is taken into account. This new method is particularly well-suited for long time simulations.

Skew-Radial Basis Function Expansions for Empirical Modeling

Arta A. Jamshidi and Michael J. Kirby

SIAM J. Sci. Comput. 31, pp. 4715-4743 (29 pages)

Online Publication Date: February 03, 2010

Full Text: | Download PDF

Show Abstract
We propose a skew-radial basis function (sRBF) expansion for the empirical model fitting problem. sRBFs employ a standard radial term, which is now made asymmetric by modulating, or skewing it with another function. The additional parameters in the skewing function permit the composite radial basis function to more flexibly adapt its shape to the data. We present several examples that illustrate the utility of sRBF representations for both the overdetermined data fitting problem and the data interpolation problem. We derive conditions under which skew perturbations of positive definite interpolation matrices remain positive definite. We observe that the sRBFs are particularly effective for producing uniform approximations and fitting jump discontinuities. We present an application to the time-series prediction of the maximum wind intensity of a hurricane and outline future work in image processing. The resulting sRBF models have reduced order, improved accuracy, and interpolation matrices with lower condition numbers.

A Novel Multigrid Method for Sn Discretizations of the Mono-Energetic Boltzmann Transport Equation in the Optically Thick and Thin Regimes with Anisotropic Scattering, Part I

B. Lee

SIAM J. Sci. Comput. 31, pp. 4744-4773 (30 pages) | Cited 1 time

Online Publication Date: February 03, 2010

Full Text: | Download PDF

Show Abstract
This paper presents a new multigrid method applied to the most common Sn discretizations (Petrov–Galerkin, diamond-differenced, corner-balanced, and discontinuous Galerkin) of the mono-energetic Boltzmann transport equation in the optically thick and thin regimes, and with strong anisotropic scattering. Unlike methods that use scalar DSA diffusion preconditioners for the source iteration, this multigrid method is applied directly to an integral equation for the scalar flux. Thus, unlike the former methods that apply a multigrid strategy to the scalar DSA diffusion operator, this method applies a multigrid strategy to the integral source iteration operator, which is an operator for 5 independent variables in spatial three dimensions (3-d), 3 in space and 2 in angle, and 4 independent variables in spatial 2-d (2 in space and 2 in angle). The core smoother of this multigrid method involves applications of the integral operator. Since the kernel of this integral operator involves the transport sweeps, applying this integral operator requires a transport sweep (an inversion of an upper triangular matrix) for each of the angles used. As the equation is in 5-space or 4-space, the multigrid approach in this paper coarsens in both angle and space, effecting efficient applications of the coarse integral operators. Although each $V$-cycle of this method is more expensive than a $V$-cycle for the DSA preconditioner, since the DSA equation does not have angular dependence, the overall computational efficiency is about the same for problems where DSA preconditioning is effective. This new method also appears to be more robust over all parameter regimes than DSA approaches. Moreover, this new method is applicable to a variety of Sn spatial discretizations, to problems involving a combination of optically thick and thin regimes, and more importantly, to problems with anisotropic scattering cross-sections, all of which DSA approaches perform poorly or are not applicable at all. This multigrid approach is most effective in neutron scattering applications, where the total cross-section coefficient $\sigma_{t}$ and spatial meshsize $h$ satisfy $\sigma_{t}h\approx1$. For this case, coarsening can be done aggressively. For problems with $\sigma_{t}h\approx10$, this multigrid scheme requires a moderately coarsened multiple-coarsening scheme. An even slow coarsening, an angle semicoarsening, is required for problems with $\sigma_{t}h$ ranges between 100 and 1000, which occur in high-energy photon applications.

Reducing the I/O Volume in Sparse Out-of-core Multifrontal Methods

Emmanuel Agullo, Abdou Guermouche, and Jean-Yves L'Excellent

SIAM J. Sci. Comput. 31, pp. 4774-4794 (21 pages)

Online Publication Date: February 03, 2010

Full Text: | Download PDF

Show Abstract
Sparse direct solvers, and in particular multifrontal methods, are methods of choice to solve the large sparse systems of linear equations arising in certain simulation problems. However, they require a large amount of memory (e.g., in comparison to iterative methods). In this context, out-of-core solvers may be employed: disks are used when the required storage exceeds the available physical memory. In this paper, we show how to process the task dependency graph of multifrontal methods in a way that minimizes the input/output (I/O) requirements. From a theoretical point of view, we show that minimizing the storage requirement can lead to a huge volume of I/O compared to directly minimizing the I/O volume. Then experiments on large real-world problems also show that applying standard algorithms to minimize the storage is not always efficient at reducing the volume of I/O and that significant gains can be obtained with the use of our algorithms to minimize I/O. We finally show that efficient memory management algorithms can be applied to all the variants proposed.

A Second-Order Accurate Conservative Front-Tracking Method in One Dimension

Caroline Gatti-Bono, Phillip Colella, and David Trebotich

SIAM J. Sci. Comput. 31, pp. 4795-4813 (19 pages)

Online Publication Date: February 03, 2010

Full Text: | Download PDF

Show Abstract
This paper presents a conservative front-tracking method for shocks and contact discontinuities that is second-order accurate. It is based on a volume-of-fluid method that treats the moving front with concepts similar to those of an embedded-boundary method. Special care is taken in the centering of the data to ensure the right order of accuracy at the front, and a redistribution step guarantees conservation. A suite of test problems, for tracking both shocks and contact discontinuities, is presented that confirms that the method is second-order accurate.
Close

Your Recent History Recent History Help

Recently Viewed

  1. Iterative Algorithms Based on Decoupling of Deblurring and Denoising for Image Restoration
  2. Constructing Sobol Sequences with Better Two-Dimensional Projections
  3. Asymptotic Sampling Distribution for Polynomial Chaos Representation from Data: A Maximum Entropy an...
  4. Iterated Hard Shrinkage for Minimization Problems with Sparsity Constraints
  5. Continuation of Invariant Subspaces in Large Bifurcation Problems
© 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