This is joint work with David Horak and Lubomir Riha.
We present our novel software package for solution of the constrained quadratic programming problems (QP) called FLLOP (FETI Light Layer On top of PETSc). It is an extension of PETSc framework. PETSc (Portable, Extensible Toolkit for Scientific Computation) is a suite of data structures and routines for the parallel solution of scientific applications modeled by PDE.
FLLOP was originally developed as an implementation of the FETI domain decomposition (DD) method. However, most recent advances in the design of FLLOP allow more general use. There are essentially two levels of generalization: (i) it allows to apply FETI to variational inequalities (e.g. contact problems), (ii) the algebraic part of FETI computation is generalized to a specific combination of data structures, QP transformations, direct and iterative solvers. These ingredients can make sense also out of the original FETI method. Additionally, any combination of these constraints can be specified: equality, inequality, box.
The typical workflow looks like this: natural specification of the QP by the user, a user-specified series of QP transformations, automatic or manual choice of a sensible solver, solution of the most derived QPs by the chosen solver, a series of backward transformations to get a solution of the original QP.
The efficiency of FLLOP solver implementations on several architectures is demonstrated on the large sparse linear systems of equations solved iteratively using conjugate gradients (CG). The resulting problem can be solved by CG in primal variables or after transformation by CG in dual variables. Grouping of the nodes based on DD can be even used to generate the subspace for deflated CG (DCG) applied to the undecomposed problem.
This talk gives introductory information about the ANSELM supercomputer to all perspective users. It will briefly cover infrastructure, development environment and available tools. It will also show different ways of getting access to the supercomputer. Second half of the talk will give an overview about the most significant project that used or still using our supercomputer.
For evolution equations we present a space-time method based on Discontinuous Galerkin finite elements. Space-time methods have advantages when we have to deal with moving domains and if we need to do local refinement in the space-time domain. For this method we present a multigrid approach based on space-time slaps. This method allows the use of parallel solution algorithms. In particular it is possible to solve parallel in time and space. Numerical examples will be given which show the performance of this space–time multigrid approach.
The need for parallel mesh generation of an unstructured Cartesian grid is outlined, as done within the code ZFS, developed by a group of the Institute of Aerodynamics, RWTH Aachen University. The function of the parallel mesh generation based on a triangle surface geometry is described, whose data, which is structured within an alternating digital tree, may be distributed among the working nodes. Properties of the geometry data structure and the mesh generation algorithms are revealed, which are supposed to include temporal and spatial locality. Due to these data locality properties, an initial distribution of geometry data in conjunction with a caching/paging system of the data is proposed. Test results of an implementation within ZFS of the caching/paging system are presented and discussed.
The modeling of many real phenomena from electromagnetism, elasticity, fluid flow, etc., are based on partial differential equations. Efficient software simulation of these phenomena requires accurate and computationally efficient methods for solving the equations which describe them. The difficulties that arise in finding analytical solutions of these equations lead to the necessity of discretization, one of the most powerful method for this purpose being the finite element method. When the physical quantities of interest are not the primary unknowns of the differential equation describing the problem but some derivative of them (like partial derivatives, gradient, symmetric gradient) their computation requires additional numerical differentiation which leads to a loss of accuracy. One of the strategies to increase the solution accuracy is to compute directly the desired quantities of interest.
The aim of this talk is to present a new intrinsic discretization approach for the direct computation of the fluxes in potential theory problems (electric, magnetic and gravitational potential problem). In the proposed intrinsic approach we develop local conditions for the approximation of the gradient vector field and then construct finite element spaces, i.e., local basis functions, directly from these conditions. In order to incorporate the essential boundary conditions we construct lifting operators as the left inverse of elementwise gradient operators. In the non-conforming case we infer from the proof of the second Strang lemma appropriate compatibility conditions at the interfaces between elements of the mesh so that the non-conforming perturbation of the original bilinear form can be estimated in a consistent way. We prove that these conditions lead to an optimal order of convergence. A main characteristic of the method is the decomposition of the finite element spaces in a direct sum of vertex-, edge- and triangle supported subspaces for which basis functions are obtained. Moreover we prove that the non-conforming spaces are spanned by the gradients of standard hp-finite element basis functions enriched by some edge-oriented non-conforming basis functions for even polynomial degree and by some triangle-supported non-conforming basis functions for odd polynomial degree.
We consider retarded boundary integral formulations of the three-dimensional wave equation in unbounded domains. Our goal is to apply a Galerkin method in space and time in order to solve these problems numerically. In this approach the accurate computation of the system matrix entries is the major bottleneck. In order to simplify the arising quadrature problem we use globally smooth and compactly supported basis functions for the time discretization. This furthermore easily allows the use of a variable time-stepping and a variable order of the approximation in time. In order to obtain a scheme that automatically adapts the time grid to local irregularities of the solution we use suitable a posteriori error estimators. Various numerical experiments show the behavior of the adaptive algorithm.
In this talk we address key issues concerning research software engineering. In the first part, an overview of open source libraries relevant to engineering applications is given, such as finite element, linear solver, and meshing libraries. In the second part, we focus on software engineering aspects, particularly dealing with the problem of inflexible monolithic software design and ways to avoid it. To this end, we discuss modular and library-centric software design as well as publicly available component frameworks as a means to decouple scientific applications into reusable components. We conclude the talk by introducing the ViennaStar library collection, a librarycentric-based suite of open source libraries.
We present analytical solution of the Stokes problem in rotationally symmetric domains using polar coordinates and separation of variables. This is then used to find the asymptotic behaviour of the solution in the vicinity of corners, also for Navier-Stokes equations. We also present alternative solution based on the Fourier transform. We apply the result to construct very precise numerical finite element solution.
This is joint work with T. Brzobohaty, T. Kozubek, O. Vlach.
The variationally consistent discretization of the non-penetration condition that was introduced by B. I. Wohlmuth became an important ingredient of effective algorithms for the solution of multibody contact problems of elasticity, especially when it is necessary to deal with the non-matching grids on contact interface or when the contact interface is large and curved. Here we discuss the problems related to the application of this discretization method in the FETI based algorithms. We first show that the resulting constraint matrices with normalized normals are well conditioned under mild conditions. The estimates are then used to extend the optimality results obtained earlier for simple node-to-node description of non-penetration to the mortar discretization. The results are illustrated by numerical solution of 2D and 3D academic benchmarks and of some realistic the problems.
Demand from end users who need to solve their problems which are in many cases very complex is and always has been driving force for developing new efficient algorithms. This is even more apparent in era of supercomputers. Nowadays high performance computers give their users computational power unimaginable few years ago. Demand for algorithms able to tame and utilize this power has been lately driving force for parallelization of existing and development of new parallel algorithms. In this talk overview of algorithms and codes deployed on Anselm used for solving real industrial problems together with their scalability will be presented.
Modeling the transport processes in a vadose zone plays an important role in predicting the reactions of soil biotopes to anthropogenic activity, e.g. modeling contaminant transport, the effect of the soil water regime on changes in soil structure and composition, etc. Water flow is governed by the Richards equation, while the constitutive laws are typically supplied by the van Genuchten model, which can be understood as a pore size distribution function. Certain materials with dominantly uniform pore sizes (e.g. coarse-grained materials) can exhibit ranges of constitutive function values within several orders of magnitude, possibly beyond the computer real numbers length. Thus a numerical approximation of the Richards equation often requires the solution of systems of equations that cannot be solved on computer arithmetic. An appropriate domain decomposition into subdomains that cover only a limited range of constitutive function values, and that will change adaptively, reflecting the time progress of the model, will enable an effective and reliable solution of this problem. This presentation focuses on improving the performance of a nonlinear solver by locating the areas with abrupt changes in the solution.
A continuation problem for finding successive solutions of discretised plane quasi-static contact problems with friction is proposed and studied. The corresponding first-order system is derived and results of existence and uniqueness of its solutions are obtained. Conditions guaranteeing continuation of a solution curve along a direction from the first-order system are given. Numerical continuation of solution curves is performed on examples with non-smooth folds.
A boundary-element counterpart of the domain decomposition vertex solver is proposed and tested for a 2-dimensional Poisson’s equation. While the standard theory has been developed only for triangular or quadrilateral subdomains, where harmonic base functions are available, the practical mesh-partitioners generate complex polygonal subdomains. We aim at bridging this gap. Being inspired by Hofreither, Langer, and Pechstein (2010) we construct the coarse solver on general polygonal subdomains so that the local coarse stiffness matrix is approximated by a boundary-element discretization of the Steklov-Poincaré operator. The efficiency of our approach is documented by a substructuring into L-shape domains, which is robust with respect to material coefficient jumps.
Domain Decomposition Methods (DDMs) is a research area which offers one of the best techniques of constructing effective parallel preconditioners. The term Domain Decomposition refers to the splitting of an original differential equation problem on a given domain or an approximation thereof into a collection of coupled subproblems on subdomains forming decomposition of the original domain. In this talk we discuss some basic linear DDMs frameworks like multiplicative Schwarz method (MSM), additive Schwarz method (ASM) and some substructuring methods like e.g. Neumann-Neumann type of DDMs.
In the second part of the talk we also present nonlinear DDM namely Nonlinearly Preconditioned Inexact Newton Methods (ASPIN) for solving the nonlinear systems of equations. The method proved itself to be more efficient than classical inexact Newton algorithms for many problems.
Discretizations of differential equation problems usually are constructed on the basis of one global mesh (triangulation) of the domain, in which this differential problem is defined. However, in many applications there is a need to construct discretizations on independent nonmatching meshes in subdomains or to use different discretizations methods in the given subdomains (substructures) of the original domain. One of the methods of constructing discretizations of differential equation problems on nonmatching meshes is the mortar method proposed by Ch. Bernardi, Y. Madey, and T. Patera. In this talk we will present the construction of mortar discretizations for few types of conforming and nonconforming finite elements for second and fourth order problems. We also present some domain decomposition parallel preconditioners for the arising algebraic systems of equations.
This is joint work with Jan Zapletal.
In this talk we present the first results of our new parallel boundary element library. BEM offers a suitable alternative to FEM especially in the cases when the computational domain is unbounded (e.g. in sound scattering or electromagnetism problems) since it only introduces unknowns on the boundary of the computational area. However the system matrices resulting from BEM are dense with O(n^2) entries requiring O(n^2) operations to assemble. To reduce these main setbacks we combine the vectorization of the computationally most demanding parts with their shared memory parallelization by OpenMP and sparsification of the resulting matrices by the fast multipole method. We present the first scalability results, give the comparison of efficiency obtained by various compilers, and discuss the possibilities of future development.
We deal with a pressure-correction method for solving unsteady incompressible flows. In this approach, five subsequent equations are solved within each time step. These correspond to three scalar convection-diffusion problems, one for each component of velocity, a pure Neumann problem for the correction of pressure, and a problem of the L2 projection for pressure update. We present a comparative study of several parallel preconditioners and Krylov subspace methods from the PETSc library and investigate their suitability for solving the arising linear systems after discretizing by the finite element method. The target application are large-scale simulations of flows around wings of insects. This is a joint work with Fehmi Cirak.
This is joint work with Victorita Dolean, Patrice Hauret, Pierre Jolivet, Frédéric Nataf, Clemens Pechstein, Daniel J. Rixen and Robert Scheichl.
Domain decomposition methods are a popular way to solve large linear systems. For problems arising from practical applications it is likely that the equations will have highly heterogeneous coefficients. For example a tire is made both of rubber and steel, which are two materials with very different elastic behaviour laws. Many domain decomposition methods do not perform well in this case, specially if the decomposition into subdomains does not accommodate the coefficient variations. For three popular domain decomposition methods (Additive Schwarz, BDD and FETI) we propose a remedy to this problem based on local spectral decompositions. Numerical investigations will confirm robustness with respect to heterogeneous coefficients and automatic (non regular) partitions into subdomains. We will also present large scale computations (over a billion unknowns) conducted by Jolivet in Freefem++ which show that strong scalability is achieved.
[1] F. Nataf, H. Xiang, and V. Dolean. A two level domain decomposition preconditioner based on local Dirichlet-to-Neumann maps. Comptes Rendus Mathématique, 348(21-22) :1163–1167, 2010.This is joint work with M. Cermak and J.Haslinger
The contribution deals with a static case of discretized elasto-perfectly plastic problems obeying Hencky's law in combination with frictionless contact boundary conditions. The main interest is focused on the analysis of the formulation in terms of displacements, limit load analysis and related numerical methods. This covers the study of: i) the dependence of the solution set on the loading parameter ζ, ii) relation between ζ and the parameter α representing the work of external forces, iii) loading process controlled by ζ and by α, iv) numerical methods for solving two different problems, with prescribed values of ζ and α, respectively. The numerical methods are implemented in combination with the Total-FETI domain decomposition method within the MatSol library in Matlab. Some theoretical results and convergence properties are illustrated on numerical examples.
Numerical dispersion, or what is often referred to as the pollution effect, presents a challenge to an efficient finite element discretization of the Helmholtz, Navier equations, and related equations describing acoustic wave propagation in fluids and structures in the medium frequency regime. As a result, a fine discretization is required leading to the solution of large systems of algebraic equations. Two strategies are discussed here to alleviate the difficulty involved in a discretization and a repetitive solution of these equations.
The first part of the talk focuses on frequency sweep problems that arise in many structural dynamic, acoustic, and structural acoustic applications. In each case, they incur the evaluation of a frequency response function for a large number of frequencies. Since each function evaluation requires the solution of an often large-scale system of equations, frequency sweep problems can easily become prohibitively expensive. Interpolatory model order reduction methods have proven to be a powerful tool for reducing their cost. The performance of interpolatory MOR methods depends on the location and number of the interpolation frequency points, and the number of consecutive frequency derivatives matched at each frequency point. An automatic adaptive strategy based on monitoring a relative residual in the frequency band of interest is discussed. The robustness, accuracy, and computational efficiency of this adaptive strategy are highlighted with the solution of several frequency sweep problems associated with large-scale structural dynamic, acoustic, and structural acoustic finite element models.
The second part of the talk reviews the discontinuous enrichment method. To alleviate the pollution effect and improve the unsatisfactory pre-asymptotic convergence of the classical Galerkin finite element method with piecewise polynomial basis functions, several discretization methods based on plane wave bases have been proposed. The discontinuous enrichment method is one such method and has has been shown to offer superior performance to the classical Galerkin finite element method for a number of constant wavenumber Helmholtz problems and has also outperformed two representative methods that use plane waves - the partition of unity and the ultra-weak variation formulation methods. Applications of DEM ranging from the Helmholtz equation for acoustics in fluids, through the vibration of plates, to the wave equation indicate the potential of such alternative discretization strategies.
Numerical realization of mathematical homogenization of elliptic equations is based on homogenization theorem which is in fact quite unfriendly to boundary element methods. Homogenized coefficients of a material are given by an integral formula which includes derivation of a function - solution of an auxiliary periodic equation. In general even if we use BEM for solving this auxiliary problem, we still need values of the solution in volume for computing the homogenized coefficient matrix. Fortunately it can be shown that for some composite materials we can compute the coefficients directly using only values of the auxiliary function at the border between the parts of the composite material. This leads to effective use of indirect BEM in mathematical homogenization.
This is extension of paper of P. Harasim and J. Valdman on Verification of functional a posteriori error estimates for obstacle problem in 1D. Benchmarks with known analytical solutions for contact problem will be implemented numerically and error of discrete approximations will be estimated in energetic norm based on concepts introduced by S. Repin.
This is joint work with Roman Kuzel
Our method is a kind of exact interpolation scheme proposed by Achi Brandt In the exact interpolation scheme, for x being the fine-level approximation of the solution, the coarse-space V=V(x) is constructed so that x satisfies x ∈ V. We achieve it simply by adding the vector x as a first column of the prolongator corresponding to a general purpose coarse-space. (The columns of the prolongator P form a computationally relevant basis of the coarse-space V=Im(P).) The advantages of this construction become obvious when solving non-linear problems. The cost of enriching the coarse-space V by the current approximation x is a single dense column of the prolongator that has to be updated each iteration. Our method can be used for multilevel acceleration of virtually any iterative method used for solving both linear and non-linear systems. The method is tested, with excellent results, on nuclear reactor criticality computations.
