Numerical Analysis Archives for Fall 2017 to Spring 2018
Adaptive approximation of the Monge-Kantorovich problem
When: Tue, August 2, 2016 - 3:30pmWhere: Kirwan Hall 1311
Speaker: Prof. Soeren Bartels (Department of Mathematics, University of Freiburg, Germany) - http://aam.uni-freiburg.de/abtlg-en/ls/lsbartels/FrontPage?set_language=en
Abstract: Optimal transportation problems define high-dimensional linear
programs. An efficient approach to their numerical solution is
based on reformulations as nonlinear partial differential
equations. If transportation cost is proportional to distance
this leads to the Monge--Kantorovich problem which is a
constrained minimization problem on Lipschitz continuous
functions. We discuss the iterative solution via splitting
methods and devise an adaptive mesh refinement strategy based on
an a~posteriori error estimate for the primal-dual gap. This is
joint work with Patrick Schoen (University of Freiburg).
Solutions of Non-Linear Differential Equations with Feature Detection Using Fast Walsh Transforms
When: Tue, September 6, 2016 - 3:30pmWhere: AV Williams 3258
Speaker: Peter Gnoffo (NASA Langley Research Center) -
Abstract: Walsh functions form an orthonormal basis set consisting of square waves. Square waves make the system well suited for detecting and representing functions with discontinuities. Given a uniform distribution of 2 cells on a one-dimensional element, it is proved that the inner product of the Walsh Root function for group p with every polynomial of degree < (p - 1) across the element is identically zero. It is also proved that the magnitude and location of a discontinuous jump, as represented by a Heaviside function, are explicitly identified by its Fast Walsh Transform (FWT) coefficients. These two proofs enable an algorithm that quickly provides a Weighted Least Squares fit to distributions across the element that include a discontinuity. It is shown that flux reconstruction relative to the FWT fit in partial differential equations provides improved accuracy and eliminates the need for flux limiting in the vicinity of a discontinuity. The detection of a discontinuity further enables analytic relations to locally describe its evolution and provide increased accuracy. Examples are provided for time-accurate advection, Burgers equation, and quasi-one-dimensional nozzle flow.
ShapeFit and ShapeKick for Robust, Scalable Structure from Motion
When: Wed, September 14, 2016 - 3:30pmWhere: AV Williams 3258
Speaker: Paul Hand (Department of Computational and Applied Mathematics, Rice University) - http://www.caam.rice.edu/~hand/
Abstract: We consider the problem of recovering a set of locations given observations of the direction between pairs of these locations. This recovery task arises from the Structure from Motion problem, in which a three-dimensional structure is sought from a collection of two-dimensional images. In this context, the locations of cameras and structure points are to be found from Epipolar geometry and point correspondences among images. These correspondences are often incorrect because of lighting, shadows, and the effects of perspective. Hence, the resulting observations of relative directions contain significant outliers. We introduce a new method for outlier-tolerant location recovery from pairwise directions. This method, called ShapeFit, is a convex Second Order Cone Program that can be efficiently solved. Empirically, ShapeFit can succeed on synthetic data with over 50% corruption. Rigorously, we prove that ShapeFit can recover a set of locations exactly when a fraction of the measurements are adversarially corrupted and when the data model is random. On real data, an ADMM implementation of ShapeFit yields performance comparable to the state-of-the-art with an order of magnitude speed-up. Our proposed numerical framework is flexible in that it accommodates other approaches to location recovery and can be used to speed up other methods. These properties are demonstrated by extensively testing against state-of-the-art methods for location recovery on 13 large, irregular collections of images of real scenes.
Scalable methods for machine learning and sparse signal recovery
When: Wed, October 5, 2016 - 2:00pmWhere: 4122 CSIC
Speaker: Tom Goldstein (University of Maryland, Department of Computer Science) - https://www.cs.umd.edu/~tomg/
Abstract: However, these resources present a host of new algorithmic challenges. Practical algorithms for large-scale data analysis must scale well across many machines, have low communication requirements, and have low (nearly linear) runtime complexity to handle extremely large problems.
In this talk, we discuss alternating direction methods as a practical and general tool for solving a wide range of model-fitting problems in a distributed framework. We then focus on new "transpose reduction" strategies that allow extremely large regression problems to be solved quickly on a single node. We will study the performance of these algorithms for fitting linear classifiers and sparse regression models on tera-scale datasets using thousands of cores.
Algorithms for Flow Ensembles
When: Tue, October 18, 2016 - 3:30pmWhere: 3258 AV Williams
Speaker: Bill Layton (University of Pittsburgh) - http://www.math.pitt.edu/~wjl/
Abstract: It is long known that assessing the reliability and improving the accuracy of fluid flow simulations requires computing ensembles. However, the Navier-Stokes equations are sensitive to mesh resolution: simulations can be plausible but wrong for under resolved meshes. These two issues introduce the inevitable competition in memory, time and other resources of resolution vs. ensemble computations. This talk will present some first steps in resolving this competition. Interestingly, new algorithms mean new turbulence models are possible. This work is joint with Nan Jiang, Missouri Univ. Science and Technology.
Monotone numerical methods for nonlinear parabolic and integro-parabolic problems
When: Tue, October 25, 2016 - 3:30pmWhere: 3258 AV Williams
Speaker: Igor Boglaev (Massey University, University of New Zealand) - http://www.massey.ac.nz/massey/learning/colleges/college-of-sciences/about/fundamental-sciences/staff-list.cfm?stref=754330
Abstract: The talk is concerned with monotone numerical methods for nonlinear
parabolic and integro-parabolic problems. The basic idea of the iterative
methods for the computation of numerical solutions is the monotone ap-
proach which involves the notion of upper and lower solutions and the con-
struction of monotone sequences from a suitable linear discrete system.
The monotone property of the iterations gives improved upper and lower
bounds of the solution in each iteration. Error estimates between the com-
puted approximations and the solutions of the nonlinear discrete problems
are obtained for each monotone iterative method. The monotone conver-
gence property is used to prove the convergence of the nonlinear discrete
problems to the corresponding differential problems as mesh sizes decrease
to zero.
Applications are given to several models arising from physical, chemical
and biological systems. Numerical experiments are given to some of these
models, including a discussion on a rate of convergence of the monotone
sequences.
Mixed-Integer PDE-Constrained Optimization
When: Tue, November 1, 2016 - 3:30pmWhere: 3258 AV Williams
Title: Mixed-Integer PDE-Constrained Optimization
Speaker: Sven Leyffer (Argonne National Laboratory), https://wiki.mcs.anl.gov/leyffer/index.php/Sven_Leyffer
Abstract: Many complex applications can be formulated as optimization problems constrained bypartial differential equations (PDEs) with integer decision variables. Examples include the remediation of contaminated sites and the maximization of oil recovery; the design of next generation solar cells; the layout design of wind-farms; the
design and control of gas networks; disaster recovery; and topology optimization.
We will present emerging applications of mixed-integer PDE-constrained optimization,
review existing approaches to solve these problems, and highlight their computational and
mathematical challenges. We introduce a new set of benchmark set for this challenging
class of problems, and present some early numerical experience using both mixed-integer
nonlinear solvers and heuristic techniques.
Finding Low-Rank Solutions via the Burer-Monteiro Approach, Efficiently and Provably
When: Tue, November 8, 2016 - 3:30pmWhere: 3258 AV Williams
Title: Finding Low-Rank Solutions via the Burer-Monteiro Approach, Efficiently and Provably
Speaker: Anastasios Kyrillidis, University of Texas, http://akyrillidis.github.io/about/
Abstract: A low rank matrix can be described as the outer product of two tall matrices, where the total number of variables is much smaller. One could exploit this observation in optimization: e.g., consider the minimization of a convex function f over low-rank matrices, where the low-rank set is modeled via such factorization. This is not the first time such a heuristic has been used in practice. The reasons of this choice could be three-fold: (i) it might model better an underlying task (e.g., f may have arisen from a relaxation of a rank constraint in the first place), (ii) it might lead to computational gains, since smaller rank means fewer variables to maintain and optimize, (iii) it leads to statistical gains, as it might prevent over-fitting in machine learning or inference problems.
Though, such parameterization comes at a cost: the objective is now a bilinear function over the variables, and thus non-convex with respect to the factors. Such cases are known to be much harder to analyze. In this talk, we will discuss and address some of the issues raised by working directly on such parameterization: How does the geometry of the problem change after this transformation? What can we say about global/local minima? Does this parameterization introduce spurious local minima? Does initialization over the factors play a key role, and how we can initialize in practice? Can we claim any convergence guarantees under mild assumptions on the objective f? And if yes, at what rate?
Modern Multiscale Methods for First-Principles Collisionless Plasma Simulation
When: Tue, November 15, 2016 - 3:30pmWhere: 3258 AV Williams
Speaker: Luis Chacon (Los Alamos National Laboratory) - http://www.lanl.gov/expertise/profiles/view/luis-chacon
Abstract: Collisionless plasmas are described by the Vlasov-Maxwell equations. This set of equations is high-dimensional (spanning three spatial and three velocity dimensions), highly nonlinear, and remarkably multi-scale, supporting disparate time and length scales. These features make its efficient numerical integration extremely challenging. The high-dimensionality of these equations have made particle methods quite attractive. In the plasma context, the method is termed particle-in-cell (PIC). PIC is naturally adaptive in velocity space, and resolves the curse of dimensionality by allowing the estimation of moment integrals be independent of the dimensionality of the underlying phase space. However, PIC can be noisy, and is generally problematic for long-term integrations of
the Vlasov-Maxwell equations due to both accuracy limitations (e.g., lack of energy conservation results in secular energy growth which subtracts fidelity from the simulation) and efficiency ones (particle methods are typically explicit, and feature both temporal and spatial stability constraints that force the resolution of both the fastest frequencies and the smallest length-scales supported by the model). These limitations make the long-term,
system-scale PIC simulation of physical systems intractable, even with the most powerful supercomputers.
Recently, fully implicit, nonlinear algorithms have been proposed for both electrostatic [1,2] and electromagnetic [3,4,5] PIC descriptions that enable for the first time truly multiscale kinetic simulations of collisionless plasmas with particle methods. These algorithms 1) feature exact conservation properties, thus avoiding secular growth of
conserved quantities, and 2) eliminate the numerical stability constraints of explicit PIC, thus allowing the use of large time steps and mesh sizes (compatible with the physics of interest). The approaches, based on modern nonlinear iterative methods, minimize the solver memory footprint by nonlinearly enslaving the kinetic component (particles) to the field equations. Thus, only field equations appear explicitly in the nonlinear residual, resulting in only modest memory requirements for the nonlinear solver. Only a single copy of the particles is needed, as with explicit implementations.
However, despite drastically decreasing the number of degrees of freedom (by allowing larger mesh cells) and the number of time steps required for a given simulation (by allowing large time steps), the resulting algebraic system remains extremely ill-conditioned and requires effective preconditioning for efficiency. A powerful advantage of nonlinear kinetic enslavement is that one can explore moment-based preconditioners. In this talk,
we will introduce the fully implicit, conservative multidimensional PIC method, and our approach to fluid preconditioning. We demonstrate the promise of the approach with various challenging numerical examples.
References:
[1] G. Chen, L. Chacón, D. C. Barnes, “An energy- and charge-conserving, implicit, electrostatic particle-in-cell algorithm,” J. Comput. Phys., 230, 7018–7036 (2011)
[2] G. Chen, L. Chacón, C. Leibs, D. A. Knoll, W. Taitano, “Fluid preconditioning for Newton-Krylov-based, fully implicit, electrostatic particle-in-cell simulations”, J. Comput. Phys., 258, p. 555–567 (2014)
[3] G. Chen, L. Chacón, “An energy- and charge-conserving, nonlinearly implicit, electromagnetic 1D-3V Vlasov–Darwin particle-in-cell algorithm,” Comput. Phys. Commun., 185 (10), 2391-2402 (2014)
[4] G. Chen, L. Chacón, “An energy- and charge-conserving, nonlinearly implicit, multidimensional electromagnetic Vlasov–Darwin particle-in-cell algorithm,” Comput. Phys. Commun., 197, 73-87 (2015)
[5] L. Chacón, G. Chen, “A curvilinear, fully implicit, conservative electromagnetic PIC algorithm in multiple dimensions,” J. Comput. Phys., 316, 578–597 (2016)
TBA
When: Tue, November 29, 2016 - 3:30amWhere: 3258 AV Williams
Speaker: Michele Benzi (Emory University) - http://www.mathcs.emory.edu/~benzi/
Variable coefficients and numerical methods for electromagnetic waves - Joint Numerical Analysis - CSCAMM Seminar
When: Wed, January 25, 2017 - 2:00pmWhere: CSIC 4122
Speaker: Dr. Lise-Marie Imbert-Gerard (Courant Institute, NYU) - http://www.cims.nyu.edu/~imbertgerard/index.html
Abstract: In the first part of the talk, we will discuss a numerical method for wave propagation in inhomogeneous media. The Trefftz method relies on basis functions that are solution of the homogeneous equation. In the case of variable coefficients, basis functions are designed to solve an approximation of the homogeneous equation. The design process yields high order interpolation properties for solutions of the homogeneous equation. This introduces a consistency error, requiring a specific analysis.
In the second part of the talk, we will discuss a numerical method for elliptic partial differential equations on manifolds. In this framework the geometry of the manifold introduces variable coefficients. Fast, high order, pseudo-spectral algorithms were developed for inverting the Laplace-Beltrami operator and computing the Hodge decomposition of a tangential vector field on closed surfaces of genus one in a three dimensional space. Robust, well-conditioned solvers for the Maxwell equations will rely on these algorithms.
Grid-adaptation for large eddy simulations of statistically stationary turbulent flows
When: Tue, February 7, 2017 - 3:30pmWhere: 3258 AV Williams
Speaker: Johan Larsson (Department of Mechanical Engineering, University of Maryland) - http://terpconnect.umd.edu/~jola
Abstract: The grid-spacing directly controls both the numerical and the modeling errors in
coarse-grained simulations of multi-scale problems, one example of which is
the large eddy simulation (LES) technique for computations of turbulence.
The modeling errors can not be estimated exactly, and therefore much of the established machinery for grid-adaptation (e.g., truncation errors, how to estimate residuals, etc) becomes less meaningful in this context. The talk will briefly introduce the LES technique and why we are interested in grid-adaptation in this context, and will then describe the development of a directional error estimator capable of driving an anisotropic grid-adaptation process. The talk will then discuss the process of proving grid-convergence in LES and propose a combined verification/adaptation approach that addresses both problems.
The proposed approach requires the ability to generate grids where the grid-spacing is modified globally by a non-integer factor; one objective of the talk is to seek collaborations in the AMSC community on this specific grid-generation problem.
Multilevel Monte Carlo Analysis for Optimal Control of Elliptic PDEs with Random Coefficients
When: Tue, February 21, 2017 - 3:30pmWhere: 3258 AV Williams
Speaker: Elisabeth Ullmann (Department of Mathematics, TU Munich) - https://www-m2.ma.tum.de/bin/view/Allgemeines/Ullmann
Abstract: This work is motivated by the need to study the impact of data uncertainties and material imperfections on the solution to optimal control problems constrained by partial differential equations. We consider a pathwise optimal control problem constrained by a diffusion equation with random coefficient together with box constraints for the control. For each realization of the diffusion coefficient we solve an optimal control problem using the variational discretization [M. Hinze, Comput. Optim. Appl., 30 (2005), pp. 45-61]. Our framework allows for lognormal coefficients whose realizations are not uniformly bounded away from zero and infinity.
We establish finite element error bounds for the pathwise optimal controls. This analysis is nontrivial due to the limited spatial regularity and the lack of uniform ellipticity and boundedness of the diffusion operator. We apply the error bounds to prove convergence of a multilevel Monte Carlo estimator for the expected value of the pathwise optimal controls. In addition we analyze the computational complexity of the multilevel estimator. We perform numerical experiments in 2D space to confirm the convergence result and the complexity bound.
Some fast solvers for models of fluid flow, linear elasticity and poroelasticity
When: Tue, March 14, 2017 - 3:30pmWhere: 3258 AV Williams
Speaker: Mingchao Cai (Morgan State University) -
Abstract: Fluid model, linear elasticity model and poroelasticity model have wide applications in geosciences and biomechanics. For example, blood-vessel wall interactions are modeled by using both fluid flow model and elasticity model; brain edema and cancellous bones are usually modeled by using poroelastic models. In this presentation, I will discuss the Cahouet-Chabard preconditioner using exact and inexact Multigrid solvers for incompressible fluid flow model, the two-level overlapping Schwarz methods for linear elasticity model, and my recent progress on fast solvers for poroelasticity model.
Accelerating Sparse Factorization Methods with Algorithmic and Hardware Advances
When: Tue, March 28, 2017 - 3:30pmWhere: 3258 AV Williams
Speaker: Xiaoye (Sherry) Li (Lawrence Berkeley Laboratory) - http://crd-legacy.lbl.gov/~xiaoye/
Abstract: Many extreme-scale simulation codes encompass multiphysics components in multiple spatial and length scales. The resulting discretized sparse linear systems can be highly indefinite, nonsymmetric and extremely ill-conditioned. For such problems, factorization based algorithms are often the most robust algorithmic choices among many alternatives. We present our recent research on novel parallel factorization algorithms that are efficient for solving such problems. From algorithm side, we incorporate data-sparse low-rank structures, such as hierarchical matrix algebra, to achieve lower arithmetic and communication complexity as well as robust preconditioner. From parallelization side, we exploit sparse data structures represented by DAGs and trees to schedule coarse-grained tasks and SIMD/SIMT for fine-grained parallelism. We will illustrate both theoretical and practical aspects of the methods, and demonstrate their performance on manycore architectures including GPU clusters and the latest Intel Xeon Phi KNL platforms, using a variety of real world problems.
Taylor Approximation and Variance Reduction for PDE-Constrained Optimal Control Problems Under Uncertainty
When: Tue, April 4, 2017 - 3:30pmWhere: 3258 AV Williams
Speaker: Peng Chen (ICES, University of Texas at Austin) - http://users.ices.utexas.edu/~peng/
Abstract: In this talk, we present an efficient method based on Taylor approximation for PDE-constrained optimal control problems under high-dimensional uncertainty. The computational complexity of the method does not depend on the nominal but only on the intrinsic dimension of the uncertain parameter, thus the curse of dimensionality is broken for intrinsically low-dimensional problems. Further correction for the Taylor approximation is proposed as a variance reduction method, which leads to an unbiased evaluation of the statistical moments in the objective function. We apply the method for a turbulence model with infinite-dimensional random viscosity.
Some fast solvers for models of fluid flow, linear elasticity and poroelasticity
When: Tue, April 18, 2017 - 3:30pmWhere: 3258 AV Williams
Speaker: Mingchao Cai (Morgan State University) -
Abstract: Fluid model, linear elasticity model and poroelasticity model have wide applications in geosciences and biomechanics. For example, blood-vessel wall interactions are modeled by using both fluid flow model and elasticity model; brain edema and cancellous bones are usually modeled by using poroelastic models. In this presentation, I will discuss the Cahouet-Chabard preconditioner using exact and inexact Multigrid solvers for incompressible fluid flow model, the two-level overlapping Schwarz methods for linear elasticity model, and my recent progress on fast solvers for poroelasticity model.
Variance reduction techniques for computing the trace of the inverse of a matrix
When: Tue, May 2, 2017 - 3:30pmWhere: 3258 AV Williams
Speaker: Andreas Stathopolous (College of William and Mary) - http://www.cs.wm.edu/~andreas/
Abstract: The computation of the trace of the inverse of a large matrix, A, appears
in many statistics, machine learning, and quantum mechanics applications.
Our driving application comes from Lattice Quantum Chromodynamics (LQCD).
When the size of A does not allow the use of direct methods, most
applications rely on the Hutchinson method, a Monte Carlo variant that
requires the solution of a linear system with A per step.
The variance of the estimator determines the accuracy and therefore the
computational cost of the method. For Hutchinson, the variance equals the
squared Frobenious norm of the matrix inverse, excluding its diagonal.
In this talk, we present some variance reduction techniques we developed
over the last few years that speed up the Hutchinson method. The goal is
to approximate the off-diagonal elements of the inverse of A based either
on structural or on spectral information.
For LQCD, where the discretization space is a 4D regular lattice torus,
our Hierarchical Probing method produces a sequence of Hadamard vectors
that, if used in the Hutchinson method, hierarchically annihilate elements
of the inverse of A whose vertices in the lattice are increasingly farther
from each other. Based on the decay of Green's function, this approach
has yielded significant variance reduction.
Our second approach is to deflate A from its smallest singular triplets.
The hope is that the variance of the deflated A is (much) smaller. Contrary
to low rank matrix approximations, the above deflation may actually
increase variance. We provide an analysis based on standard random unitary
matrices, and derive criteria on when to expect improvement.
Finally, we touch upon the daunting computational task of computing
the smallest singular triplets of A, and the recent progress our group
has made in this direction which are included in the software PRIMME.
Adaptive methods for high-dimensional elliptic PDEs
When: Tue, May 9, 2017 - 3:30pmWhere: 3258 AV Williams
Speaker: Markus Bachmyer (University of Bonn) -
Abstract: This talk gives an overview of recent results on adaptive solvers for PDEs posed on very high-dimensional domains. These methods combine low-rank tensor decompositions with sparse basis expansions of the lower-dimensional tensor components. The iterative schemes are guaranteed to converge, yield computable a posteriori error bounds, and are shown to be of near-optimal computational complexity in terms of the approximability of the solution. As central ingredient for these complexity estimates, we obtain bounds on the tensor ranks that can arise in the iteration.
Based on joint works with Wolfgang Dahmen and with Reinhold Schneider.