In an influential article, Marfurt suggested that the best scheme for computational wave propagation would involve an averaging of the consistent and lumped finite element approximations. Many authors have considered how this might be accomplished for first order approximation, but the case of higher orders remained unsolved.  We describe recent work on the dispersive and dissipative properties of a novel scheme for computational wave propagation obtained by averaging the consistent (finite element) mass matrix and lumped (spectral element) mass matrix. The objective is to obtain a hybrid scheme whose dispersive accuracy is superior to both of the schemes. We present the optimal value of the averaging constant for all orders of finite elements and proved that for this value the scheme is two orders more accurate compared with finite and spectral element schemes, and, in addition, the absolute accuracy is of this scheme is better than that of finite and spectral element methods.   Joint work with Hafiz Wajid, COMSATS Institute of Technology, Pakistan.

We address the problem of efficiently solving the algebraic linear systems of equations arising from the discretization of a symmetric, elliptic boundary value problem using hp discontinuous Galerkin finite element methods. We introduce a class of domain decomposition preconditioners based on the Schwarz framework, and prove bounds on the condition number of the resulting iteration operators. Numerical results confirming the theoretical estimates are also presented.

Joint work with Paul Houston, University of Nottingham, UK.

We propose iterative methods for the solution of the linear systems resulting from Interior Penalty (IP)  discontinuous Galerkin (DG)  approximations of elliptic problems. The precoonditioners are derived from a natural decomposition of the DG finite element spaces. We present the convergence analysis of the solvers for both symmetric and non-symmetric IP schemes. Extension to problems with jumps in the coefficients and linear elasticity will also be discussed. We describe in detail the preconditioning techniques for low order (piece-wise linear ) IP methods and we indicate how to proceed in the case of high order (odd degree) approximations. The talk is based on joint works with Ludmil T. Zikatanov from  Penn State University (USA).

We will discuss our on-going investigation of hp-adaptive finite elements. We will focus on a posteriori error estimates based on superconvergent derivative recovery. Besides providing both global error estimates and local error indicators, this family of error estimates also provides information that forms the basis of our hp-adaptive refinement and coarsening strategies. In particular, these a posteriori error estimates address in a cost efficient and natural way the critical issue of deciding between h or p refinement/coarsening. Some numerical examples will be provided.

We propose and analyze two-level overlapping additive Schwarz preconditioners for the local discontinuous Galerkin discretization. We prove a condition number estimate and show numerically that the method is scalable in terms of linear iterations. We also present numerical evidence that a parallel implementation of the method shows good scalability and speedup.

There are two main difficulties in solving higher order elliptic boundary value problems: the discretization schemes are more complicated and the discrete problems are very ill-conditioned.  In this talk we will discuss an approach that can solve higher order problems with an efficiency similar to that for second order problems. It is based on discontinuous Galerkin methods, embedded multigrid algorithms and domain decomposition techniques.  

We are interested in solving coupled systems of partial differential equations on computers with a large number of processors. With some combinations of domain decomposition and multigrid methods, one can easily design algorithms that are highly scalable in terms of the number of linear and nonlinear iterations. However, if the goal is to minimize the total compute time and keep the near ideal scalability at the same time, then the task is more difficult. We discuss some recent experience in solving the shallow water equations on the sphere for the modeling of the global climate.

This is a joint work with C. Yang.

Five classes of up to 13 a posteriori error estimators compete in three second-order model cases, namely the conforming and non-conforming first-order approximation of the Poisson-Problem plus some conforming obstacle problem. Since it is the natural first step, the error is estimated in the energy norm exclusively — hence the competition has limited relevance. The competition allows merely guaranteed error control and excludes the question of the best error guess. Even nonsmooth problems can be included. For a variational inequality, Braess considers Lagrange multipliers and some resulting auxiliary equation to view the a posteriori error control of the error in the obstacle problem as computable terms plus errors and residuals in the auxiliary equation. Hence all the former a posteriori error estimators apply to this nonlinear benchmark example as well and lead to surprisingly accurate guaranteed upper error bounds. This approach allows an extension to more general boundary conditions and a discussion of efficiency for the affine benchmark examples. The Luce-Wohlmuth and the least-square error estimators win the competition in several computational benchmark problems. Novel equilibration of nonconsistency residuals and novel conforming averaging error estimators win the competition for Crouzeix-Raviart nonconforming finite element methods. Our numerical results provide sufficient evidence that guaranteed error control in the energy norm is indeed possible with efficiency indices between one and two. Furthermore, accurate error control is slightly more expensive but pays off in all applications under consideration while adaptive mesh-refinement is sufficiently pleasant as accurate when based on explicit residual-based error estimates. Details of our theoretical and empirical ongoing investigations will be found in the papers quoted below.

References:

1. S. Bartels, C. Carstensen, R. Klose, An experimental survey of a posteriori Courant finite element error control for the Poisson equation, Adv. Comput. Math., 15 (2001), pp. 79-106.
2. C. Carstensen, C. Merdon, Estimator competition for Poisson problems, J. Comp. Math., 28 (2010), pp. 309-330.
3. C. Carstensen, C. Merdon, Computational survey on a posteriori error estimators for nonconforming finite element methods, Part I: Poisson problems (in preparation)
4. C. Carstensen, C. Merdon, A posteriori error estimator competition for conforming obstacle problems, Part I: Theoretical findings (in preparation), Part II: Numerical results (in preparation)

---

^* This work was supported by DFG Research Center MATHEON and by the World Class University (WCU) program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education, Science and Technology R31-2008-000-10049-0.

In this work we investigate the numerical solution for two-dimensional Maxwell's equations on graded meshes. The approach is based on the Hodge decomposition. The solution u of Maxwell's equations is approximated by solving standard second order elliptic problems. Quasi-optimal error estimates for both u and curl of u in the L₂ norm are obtained on graded meshes. We prove the uniform convergence of the W-cycle and full multigrid algorithms for the resulting discrete problem.

The first part of this talk provides a basic introduction to the building blocks of domain decomposition solvers. Specific details are given for both the classical overlapping Schwarz (OS) algorithm and a recent iterative substructuring (IS) approach called balancing domain decomposition by constraints (BDDC). A more recent hybrid OS-IS approach is also described. The success of domain decomposition solvers depends critically on the coarse space. Similarities and differences between the coarse spaces for OS and BDDC approaches are discussed, along with how they can be obtained from discrete harmonic extensions. Connections are also made between coarse spaces and multiscale modeling approaches from computational mechanics. As a specific example, details are provided on constructing coarse spaces for incompressible fluid problems. The next part of the talk deals with a variety of implementation details for domain decomposition solvers.  These include mesh partitioning options, local and global solver options, reducing the coarse space dimension, dealing with constraint equations, residual weighting to accelerate the convergence of OS methods, and recycling of Krylov spaces to efficiently solve problems with multiple right hand sides. Some potential bottlenecks and remedies for domain decomposition solvers are also discussed. The final part of the talk concerns some recent theoretical advances, new algorithms, and open questions in the analysis of domain decomposition solvers. The focus will be primarily on the work of the speaker and his colleagues on elasticity, fluid mechanics, problems in H(curl), and the analysis of subdomains with irregular boundaries. The speaker gratefully acknowledges contributions of Jan Mandel and Olof Widlund to many topics discussed in this talk.

Radial basis functions provide flexible, meshfree approximations to functions and solutions of differential equations. Naive algorithms suffer from dense linear algebra and severe ill conditioning. Simple multiscale adaptive techniques for the nodes and shape parameters have previously proven very effective in controlling ill conditioning for small node sets. We present a new fast summation method suitable for adaptively generated basis functions with varying shape parameters. When coupled with an easily parallelized restricted additive Schwarz preconditioner, the method can find RBF coefficients in near O(N log N) time for N nodes.

Algebraic multigrid (AMG) is an important method for solving the large sparse linear systems that arise in many PDE-based scientific simulation codes.  A major component of algebraic multigrid methods is the selection of coarse grids and the construction of interpolation.  The idea of compatible relaxation (CR) was introduced by Brandt as a cheap method for measuring the quality of a coarse grid.  In this talk, we will review the theory behind CR, describe our CR-based coarsening algorithm, and discuss aspects of the method that require additional development such as coarsening for systems of PDEs.  We will also discuss CR's ability to predict the convergence behavior of the AMG method and ideas for improving the accuracy of its predictions.  Finally, we will talk about issues of parallelizing these methods to run on massively parallel computers.

In lattice QCD, a standard discretization of the Dirac operator is given by the Wilson-Dirac operator, representing a nearest neighbor coupling on a 4d torus with 3x4 variables per grid point. The operator is non-symmetric but (usually) positive definite. Its small eigenmodes are non-smooth due to the stochastic nature of the coupling coefficients. Solving systems with the Wilson-Dirac operator on state-of-the-art lattices, typically in the range of 32-64 grid points in each of the four dimensions, is one of the prominent supercomputer applications today. In this talk we will describe our experience with the domain decomposition principle as one approach to solve the Wilson-Dirac equation in parallel. We will report results on scaling with respect to the size of the overlap, on deflation techniques that turned out to be very useful for implementations on QPACE, the no 1 top green 500 special purpose computer developed tiogethe with IBM  by the SFB-TR 55 in Regensburg and Wuppertal, and on first results on adaptive approaches for obtaining an appropriate coarse system. 

In contrast to the positive definite Helmholtz equation, the deceivingly similar looking indefinite Helmholtz equation is difficult to solve using classical iterative methods. Applying directly a Krylov method to the discretized equations without preconditioning leads in general to stagnation and very large iteration counts. Using classical incomplete LU preconditioners can even make the situation worse. Classical domain decomposition and multigrid methods also fail to converge when applied to such systems. The purpose of this presentation is to investigate in each case where the problems lie, and to explain why classical iterative methods have such difficulties to solve indefinite Helmholtz problems. I will also present remedies that have been proposed over the last decade, for incomplete LU type preconditioners, domain decomposition and also multigrid methods.  

Recently a new numerical method for the two-dimensional Maxwell's equation based on the Hodge decomposition for divergence-free vector fields has been introduced by Brenner, Cui, Nan and Sung. The advantage of this new approach is that an approximation of the vector field is obtained by solving several standard second order scalar elliptic boundary value problems instead of using more complicated methods. For the linear Courant finite elements standard energy residual type a posteriori error techniques can be applied to obtain guaranteed upper bounds for the L²-error. For smooth solutions a duality argument shows reliability of an L² residual type a posteriori error estimator for the H(curl)-error. A dual weighted residual error estimator is derived for singular solutions. Numerical experiments verify reliability and show empirically efficiency of the proposed error estimators. It is shown that adaptive mesh-refinement numerically leads to optimal convergence rates for general domains that are non-convex and may include holes.

Eigenvalue problems for partial differential equations (PDEs) arise in a large number of current technological applications, e.g., in the computation of the acoustic field inside vehicles (such as cars, trains or airplanes). Another current key application is the noise compensation in highly efficient motors and turbines. For the analysis of standard adaptive finite element methods an exact solution of the discretized algebraic eigenvalue problem is required, and the error and complexity of the algebraic eigenvalue problems are ignored. In the context of eigenvalue problems these costs often dominate the overall costs and because of that, the error estimates for the solution of the algebraic eigenvalue problem with an iterative method have to be included in the adaptation process. The goal of our work is to derive adaptive methods of optimal complexity for the solution of PDE-eigenvalue problems including problems with parameter variations in the context of homotopy methods. In order to obtain low (or even optimal) complexity methods, we derive and analyse methods that adapt with respect to the computational grid, the accuracy of the iterative solver for the algebraic eigenvalue problem, and also with respect to the parameter variation. Such adaptive methods require the investigation of a priori and a posteriori error estimates in all three directions of adaptation. As a model problem we study eigenvalue problems that arise in convection-diffusion problems. We developed robust a posteriori error estimators for the discretization as well as for the iterative solver errors, first for self-adjoint second order eigenvalue problems (undamped problem, diffusion problem), and then bring in the non-selfadjoint part (damping, convection) via a homotopy, where the step-size control for the homotopy is included in the adaptation process.

Two parallel domain decomposition methods for solving the coupled Stokes-Darcy model are proposed and analyzed. One is for a steady state Stokes-Darcy model with BJ condition and the other one is for a time-dependent Stokes-Darcy model with BJS condition. Robin boundary conditions are utilized to decouple the Stokes and Darcy parts of the system. Two numerical examples are used to illustrate the features of the two methods and confirm the theoretical results respectively.

The optimal design of structures and systems described by partial differential equations (PDEs) often gives rise to large-scale optimization problems, in particular if the underlying system of PDEs represents a multi-scale, multi-physics problem. Therefore, reduced order modeling techniques such as balanced truncation model reduction, proper orthogonal decomposition, or reduced basis methods are used to significantly decrease the computational complexity while maintaining the desired accuracy of the approximation. In particular, we are interested in such shape optimization problems where the design issue is restricted to a relatively small area of the computational domain. In this case, it appears to be natural to rely on a full order model only in that specific part of the domain and to use a reduced order model elsewhere. An appropriate methodology to realize this idea consists in a combination of domain decomposition techniques and balanced truncation model reduction. We will consider such an approach for shape optimization problems associated with the time-dependent Stokes system and derive explicit error bounds for the modeling error. As an application in life sciences, we will be concerned with the optimal design of surface acoustic wave driven microfluidic biochips that are used in clinical diagnostics, pharmacology, and forensics for high-throughput screening and hybridization in genomics and protein profiling in proteomics.

In this poster, an efficient rearrangement algorithm is proposed to find the optimal shape and topology for eigenvalue problems in an inhomogeneous media. The method is based on Rayleigh quotient formulation of eigenvalue and a monotone iteration process to achieve the optimality. The common numerical approach for these problems is to start with an initial guess for the shape and then gradually evolve it, until it morphs into the optimal shape. One of the difficulties is that the topology of the optimal shape is unknown. Developing numerical techniques which can automatically handle topology changes becomes essential for shape and topology optimization problems. The level set approach based on both shape derivatives and topological derivatives has been well known for its ability to handle topology changes. However, CFL constrain significantly slows down the algorithm when the mesh is further refined. Due to the efficient rearrangement, the new method not only has the ability of topological changes but also is exempt from CFL condition. We provides numerous numerical examples to demonstrate the robustness and efficiency of our approach.

Eulerian formulations of problems with interfaces avoid the subtleties of tracking and remeshing, but do they complicate solution of the discrete equations, relative to domain decomposition methods that respect the interface? We consider two different interface problems – one involving cracking and one involving phase separation. Crack problems can be formulated by extended finite element methods (XFEM), in which discontinuous fields are represented via special degrees of freedom. These DOFs are not properly handled in a typical AMG coarsening process, which leads to slow convergence. We propose a Schwarz approach that retains AMG advantages on the standard DOFs and avoids coarsening the enriched DOFs. This strategy allows reasonably mesh-independent convergence rates, though the convergence degradation of the (lower dimensional) set of crack DOFs remains to be addressed. Phase separation problems can be formulated by the Cahn-Hilliard approach, in which the level set of a continuous Eulerian field demarcates the phases. Here, scalable preconditioners follow naturally, once the subtlety of the temporal discretization is sorted out. The first project is joint with R. Tuminaro and H. Waisman and the second with X.-C. Cai and C. Yang.

Second-order definite Maxwell problems arise in many practical applications, such as the modeling of electromagnetic diffusion in ALE-MHD simulations. Typically, these problems are discretized with Nedelec finite elements resulting in a large sparse linear system which is challenging for linear solvers due to the large nullspace of curl-operator. In this poster we describe our work on the Auxiliary-space Maxwell Solver (AMS) which is a provably efficient scalable code for solving definite Maxwell problems based on the recent Hiptmair-Xu (HX) decomposition of the lowest-order Nedelec space. We demonstrate the scalability of the method and its robustness with respect to jumps in material coefficients. We also report some results from recent work on the algebraic extension of the AMS algorithm and theory to linear systems obtained by explicit element reduction.

The numerical simulation of the coarsening of binary alloys based on the Cahn-Larch`e equations requires fast, reliable and robust solvers for Cahn-Hilliard equations with logarithmic potential. After semi-implicit time discretization (cf.  Blowey and Elliott 92), the resulting spatial problem can be reformulated as a non-smooth pde-constrained optimal control problem with cost functional induced by the Laplacian. The associated Karush-Kuhn-Tucker conditions take the form of a nonsmooth saddle point problem degenerating to a variational inclusion in the deep quench limit. Our considerations are based on recent work of Gr¨aser & Kornhuber 09 and the upcoming dissertation of Gr¨aser 10. The starting point is the elimination of the primal variable leading to a nonlinear Schur complement which turns out to be the Fr´ech`et derivative of a convex functional. Now so-called nonsmooth Schur-Newton methods can be derived as gradient-related descent methods applied to this functional. In the discrete case we can show global convergence for an exact and an inexact version independent of any regularization parameters. Local quadratic convergence or finite termination can be shown for piecewise smooth nonlinearities or in the deep quench limit respectively. The algorithm can be reinterpreted as a preconditioned Uzawa method and generalizes the well-known primal-dual active set strategy by Kunisch, Ito, and Hinterm¨uller 03. A (discrete) Allen-Cahn-type problem and a linear saddle point problem have to be solved (approximately) in each iteration step. In numerical computations we observe mesh-independent local convergence for initial iterates provided by nested iteration. In the deep quench limit, the numerical complexity of the (approximate) solution of the arising linear saddle point problem dominates the detection of the actual active set.

We present convergence and optimality results for a standard AFEM

SOLVE -> ESTIMATE -> MARK -> REFINE.

The results range from plain convergence for inf-sup problems to contraction properties and quasi-optimal rates of an AFEM with D"orfler marking for elliptic problems. Beyond that we show how these results can be used to design convergent adaptive methods for a linear parabolic pde and a nonlinear stationary Stokes problem.

The (mixed finite element) discretization of the linearized Navier-Stokes equations leads to a linear system of equations of saddle point type. The iterative solution of this linear system requires the construction of suitable preconditioners, especially in the case of high Reynolds numbers. In the past, a stabilizing approach has been suggested which does not change the exact solution but influences the accuracy of the discrete solution as well as the effectiveness of iterative solvers. This stabilization technique can be performed on the continuous side before the discretization, where it is known as ``grad-div'' stabilization, as well as on the discrete side where it is known as an ``augmented Lagrangian'' technique (and does not change the discrete solution).

We study the applicability of H-LU factorizations to solve the arising subproblems in the different variants of stabilized saddle point systems.

We compare two domain decomposition methods for a linearized contact problem. The first method we consider has been used in an engineering community; we provide theoretical and numerical evidence that this method is not scalable with respect to the number of subdomains (processors). We propose a scalable alternative and analyze its properties, both theoretically and numerically. We also solve a model problem using a combination of a primal-dual active set method, viewed as a semismooth Newton method, and the scalable alternative.

WRF-Fire consists of a fire-spread model, implemented by the level set method, coupled with the Weather Research and Forecasting model (WRF). In every time step, the fire model inputs the surface wind and outputs the heat flux from the fire. The level set method allows submesh representation of the burning region and flexible implementation of various ignition modes. This presentation will address the numerical methods in the fire module, solving the Hamilton-Jacobi level set equation, modeling real wildfire experiments, and visualization.   Visualizations by Bedrich Sousedik, Erik Anderson, and Joel Daniels.   Modeling of real fires with Adam Kochanski.   Jan Mandel, Jonathan D. Beezley, Janice L. Coen, and Minjeong Kim, Data Assimilation for Wildland Fires: Ensemble Kalman filters in coupled atmosphere-surface models, IEEE Control Systems Magazine 29, Issue 3, June 2009, 47-65   Jan Mandel, Jonathan D. Beezley, and Volodymyr Y. Kondratenko, Fast Fourier Transform Ensemble Kalman Filter with Application to a Coupled Atmosphere-Wildland Fire Model. Anna M. Gil-Lafuente, Jose M. Merigo (Eds.) Computational Intelligence in Business and Economics (Proceedings of the MS'10 International Conference, Barcelona, Spain, 15-17 July 2010), World Scientific, pp. 777-784. Also available at arXiv:1001.1588

Magnetohydrodynamics (MHD) is a fluid theory that describes Plasma Physics by treating the plasma as a fluid of charged particles.  Hence, the equations that describe the plasma form a nonlinear system that couples Navier-Stokes with Maxwell's equations. We describe how the FOSLS method can be applied to incompressible resistive MHD to yield a well-posed, H$^1$-equivalent functional minimization. To solve this system of PDEs, a nested-iteration-Newton-FOSLS-AMG-LAR approach is taken. Much of the work is done on relatively coarse grids, including most of the linearizations.  We show that at most one Newton step and a few V-cycles are all that is needed on the finest grid. Estimates of the local error and of relevant problem parameters that are established while ascending through the sequence of nested grids are used to direct local adaptive mesh refinement (LAR), with the goal of obtaining an optimal grid at a minimal computational cost. An algebraic multigrid solver is used to solve the linearization steps. A parallel implementation is described that uses a binning strategy. We show that once the solution is sufficiently resolved, refinement becomes uniform which essentially eliminates load balancing on the finest grids. The ultimate goal is to resolve as much physics as possible with the least amount of computational work. We show that this is achieved in the equivalent of a few dozen work units on the finest grid. (A work unit equals a fine grid residual evaluation). Numerical results are presented for two instabilities in a large aspect-ratio tokamak, the tearing mode and the island coalescence mode.  

An acceleration of a Galerkin boundary integral equation (BIE) method for solving the three-dimensional Laplace equation is investigated in the context of the precorrected fast Fourier transform (PFFT) scheme. The PFFT technique is an algorithm for rapid computation of the dense matrix-vector products arising in an iterative solution of discretized integral equations. In the PFFT method, the problem domain is overlaid with a regular Cartesian grid that serves as an auxiliary platform for computation. With the aid of the fast Fourier transform (FFT) procedure, the necessary influence matrices of the discretized problem are rapidly evaluated on the Fourier grid in a sparse manner resulting in a significant reduction in execution time and computer memory requirements.

---

This research was supported by the Office of Advanced Scientific Computing Research, U.S. Department of Energy, under contract DE-AC05-00OR22725 with UT-Battelle, LLC.

In contrast to most of the existing theory of adaptive finite element methods (AFEM), we design an AFEM for -Δ u = f with right hand side f in H ^-1 instead of L². This has two important consequences. First we formulate our AFEM in the natural space for f, which is nonlocal. Second, we show that decay rates for the data estimator are dominated by those for the solution u in the energy norm. This allows us to conclude that the performance of AFEM is solely dictated by the approximation class of u.

This is joint work with A. Cohen and R. DeVore.

The weakly over-penalized symmetric interior penalty (WOPSIP) method was introduced for second order elliptic problems by Brenner et al. in 2008. It belongs to the family of discontinuous Galerkin methods. We will discuss two-level additive Schwarz preconditioners for the WOPSIP method. The key ingredient of the two-level additive Schwarz preconditioner is the construction of the subdomain solvers and the coarse solver. In our approach, we consider different choices of coarse spaces and intergrid transfer operators. It is shown that the condition number estimates previously obtained for classical finite element methods also hold for the WOPSIP method. Numerical results will be presented, which illustrate the parallel performance of these preconditioners.

This is joint work with Andrew T. Barker, Susanne C. Brenner, and Li-yeng Sung.

We will report on our experiences implementing PDE solvers on Peta-Scale computers, such as the 290 000 core IBM Blue Gene system in the Jülich Supercomputing Center. The talk will have two parts, the first one reporting on our Hierarchical Hybrid Grid  method, a prototype Finite Element Multigrid Solver scaling up to a trillion (10^12)  degrees of freedom on a tetrahedral mesh by using a carefully designed matrix-free implementation.  The second part of the talk will present our work on simulating complex flow phenomena using the Lattice-Boltzmann algorithm. Our software includes parallel algorithms for treating free surfaces with the goal of simulating fully resolved bubbly flows and foams. Additionally, we will report on a parallel fluid-structure-interaction technique with many moving rigid objects. This work is directed towards the modeling of particulate flows that we can represent using fully resolved geometric models of each individual particle embedded in the flow. The talk will end with some remarks on the challenges that algorithm developers will be facing on the path to exascale in the coming decade.

Energy minimization provides a general framework for developing a family of multigrid algorithms. The proposed strategy is applicable to Hermitian, non-Hermitian, definite, and indefinite problems. Each column of the grid transfer operator P is minimized in an energy-based norm while enforcing two types of constraints: a defined sparsity pattern and preservation of specified modes in the range of P. A Krylov-based strategy is used to minimize energy, which is equivalent to solving A P = 0 with the constraints ensuring a nontrivial solution. For the Hermitian positive definite case, a conjugate gradient-based (CG) method is utilized to construct grid transfers, while methods based on generalized minimum residual (GMRES) and CG on the normal equations (CGNR) are explored for the general case.

One of the main advantages of the approach is that it is flexible, allowing for arbitrary coarsenings, unrestricted sparsity patterns, straightforward long distance interpolation, and general use of constraints, either user-defined or auto-generated. We illustrate how this flexibility can be used to adapt an algebraic multigrid scheme to an extended finite element discretization suitable for modeling fracture. Computational results are presented illustrate that this particular energy minimization scheme gives rise to mesh independent convergence rates and is relatively insensitive to the number and location of cracks being modeled.

The quality of a finite element solution hinges in particular on the approximation properties of the finite element space.  In the first part of this talk we will consider the approximation of the gradient of a target function by continuous piecewise polynomials over a simplicial, 'shape-regular' mesh and prove the following result: the global best approximation error is equivalent to an appropriate sum in terms of the local best approximation errors on the elements, which do not overlap.  This means in particular that, for gradient norms, the continuity requirement does not downgrade the local approximation potential on elements and that discontinuous piecewise polynomials do not offer additional approximation power.  In the second part of the talk we will discuss the usefulness of this result in the context of adaptive methods for partial differential equations. Joint work with Francesco Mora (Milan).

 In recent years, variants of the two-level Schwarz algorithm have been developed in collaboration between Clark Dohrmann of Sandia-Albuquerque and a group at the Courant Institute. By a modification of the coarse component of the preconditioner, borrowed in part from older domain decomposition methods of iterative substructuring type, the new methods are easier to implement for general subdomain geometries and can be made insensitive to large variations on the coefficients of the partial differential equation across the interface between the subdomains. After an introduction to the design of these methods, results on applications to almost incompressible elasticity and Reissner-Mindlin plates - solved by using mixed finite elements - and problems posed in H(div) and H(curl) will be discussed. Some of these results will appear in the doctoral dissertations of Jong Ho Lee and Duk-soon Oh, two PhD candidates at the Courant Institute.

Parallel computers with large storage capacities have paved the way for increasing both the fidelity and complexity of large-scale simulations.  Increasing fidelity leads to tighter constraints on time steps for stability of explicit methods.  Increasing complexity tends to also increase the number of physics models and variations in time scales.  Providing both a stable solution process which can accurately capture nonlinear coupling between dynamically relevant phenomena while stepping over fast waves or rapid adjustments leads us toward implicit solution approaches.   This presentation provides an overview of issues arising in large- scale, multiphysics simulations and some of the motivators for looking at implicit approaches.  We discuss several popular implicit nonlinear solver technologies and show examples of uses of them within the context of problems found in supernova, subsurface simulation, fusion, and nonlinear diffusion problems.

The nonlocal continuum dielectric model is an important extension of the classical Poisson dielectric model. This poster will report some recent results we made on the finite element analysis and fast solver development for one commonly-used nonlocal continuum dielectric model. We first prove that the finite element equation of this model has the unique solution but leads to a dense linear system, which is very expansive to be solved. Surprisingly, we then discover and prove that such a dense linear system can be converted to a system of two sparse finite element equations in a form similar to the standard mixed finite element equation. In this way, fast numerical solvers can be developed to solve the nonlocal continuum dielectric model in an optimal order. Some numerical results in free energy calculation will also be presented to demonstrate the great promise of nonlocal dielectric modeling in improving the accuracy of the classic Poisson dielectric model in computing electrostatic potential energies. This project is a joined work with Prof. Ridgeway Scott, Peter Brune, (both from University of Chicago), and Yi Jiang under the support of NSF grant #DMS-0921004.

Several numerical techniques will be presented for solving discretized partial differential equations (PDEs) by special multilevel methods based on one or no grid with nearly optimal computational complexity in a user-friendly fashion.

Hypre is a software library for solving large, sparse linear systems of equations on massively parallel computers. The library was created with the primary goal of providing users with advanced parallel preconditioners. The library features parallel multigrid solvers for both structured and unstructured grid problems. For ease of use, these solvers are accessed from the application code via hypre’s conceptual linear system interfaces, which allow a variety of natural problem descriptions. The motivation for the design of hypre, an overview of its interfaces and some of its performance highlights are presented.

  Multigrid algorithms for discretized nonlinear partial differential equations and systems are nearly as old as multigrid itself. Over the years several approaches and variants of nonlinear multigrid algorithms have been developed. Typically, for relatively easy problems the different approaches exhibit similar performance. However, for difficult problems the behavior varies, and it is not easy to predict which approach may prevail.   In this talk we will consider nonlinear multigrid, focusing on the task of coarse-grid correction, in a general framework of variational coarsening. Such a view reveals clear relations between the various existing approaches and may suggest future variants. This study also sheds light on the choice of inter-grid transfer operators, which are so important for obtaining fast multigrid convergence, and which have received much attention in linear multigrid algorithms but far less so in nonlinear multigrid.     

The multigrid method solves the finite element equations in optimal order, i.e., solving a linear system of $O(N)$ equations in $O(N)$ arithmetic operations. Based on low level solutions, we can use finite element extrapolation to obtain the high-level finite element solution on some coarse-level element boundary, at an higher accuracy $O(h_i4)$. Thus, we can solve higher level $(h_j, junderset{sim}{