Search

We present the Super-Localized Orthogonal Decomposition (SLOD) method for the numerical homogenization of linear elasticity problems with multiscale microstructures modeled by a heterogeneous coefficient field without any periodicity or scale separation assumptions. Compared to the established Localized Orthogonal Decomposition (LOD) and its linear localization approach, SLOD achieves significantly improved sparsity properties through a nonlinear superlocalization technique, leading to computationally efficient solutions with significantly less oversampling - without compromising accuracy. We generalize the method to vector-valued problems and provide a supporting numerical analysis. We also present a scalable implementation of SLOD using the deal.II finite element library, demonstrating its feasibility for high-performance simulations. Numerical experiments illustrate the efficiency and accuracy of SLOD in addressing key computational challenges in multiscale elasticity.

Traditionally, the geometric multigrid method is used with nested levels. However, the construction of a suitable hierarchy for very fine and unstructured grids is, in general, highly non-trivial. In this scenario, the non-nested multigrid method could be exploited in order to handle the burden of hierarchy generation, allowing some flexibility on the choice of the levels. We present a parallel, matrix-free, implementation of the non-nested multigrid method for continuous Lagrange finite elements, where each level may consist of independently partitioned triangulations. Our algorithm has been added to the multigrid framework of the C++ finite-element library deal.II. Several 2D and 3D numerical experiments are presented, ranging from Poisson problems to linear elasticity. We test the robustness and performance of the proposed implementation with different polynomial degrees and geometries., Comment: 24 pages, 23 figures, 14 tables

This work presents a matrix-free finite element solver for finite-strain elasticity adopting an $hp$-multigrid preconditioner. Compared to classical algorithms relying on a global sparse matrix, matrix-free solution strategies significantly reduce memory traffic by repeated evaluation of the finite element integrals. Following this approach in the context of finite-strain elasticity, the precise statement of the final weak form is crucial for performance, and it is not clear a priori whether to choose problem formulations in the material or spatial domain. With a focus on hyperelastic solids in biomechanics, the arithmetic costs to evaluate the material law at each quadrature point might favor an evaluation strategy where some quantities are precomputed in each Newton iteration and reused in the Krylov solver for the linearized problem. Hence, we discuss storage strategies to balance the compute load against memory access in compressible and incompressible neo-Hookean models and an anisotropic tissue model. Additionally, numerical stability becomes increasingly important using lower/mixed-precision ingredients and approximate preconditioners to better utilize modern hardware architectures. Application of the presented method to a patient-specific geometry of an iliac bifurcation shows significant speed-ups, especially for higher polynomial degrees, when compared to alternative approaches with matrix-based geometric or black-box algebraic multigrid preconditioners.

Quality assessment algorithms measure the quality of a captured biometric sample. Since the sample quality strongly affects the recognition performance of a biometric system, it is essential to only process samples of sufficient quality and discard samples of low-quality. Even though quality assessment algorithms are not intended to yield very different quality scores across demographic groups, quality score discrepancies are possible, resulting in different discard ratios. To ensure that quality assessment algorithms do not take demographic characteristics into account when assessing sample quality and consequently to ensure that the quality algorithms perform equally for all individuals, it is crucial to develop a fairness measure. In this work we propose and compare multiple fairness measures for evaluating quality components across demographic groups. Proposed measures, could be used as potential candidates for an upcoming standard in this important field.

We present a space-time multigrid method based on tensor-product space-time finite element discretizations. The method is facilitated by the matrix-free capabilities of the {\ttfamily deal.II} library. It addresses both high-order continuous and discontinuous variational time discretizations with spatial finite element discretizations. The effectiveness of multigrid methods in large-scale stationary problems is well established. However, their application in the space-time context poses significant challenges, mainly due to the construction of suitable smoothers. To address these challenges, we develop a space-time cell-wise additive Schwarz smoother and demonstrate its effectiveness on the heat and acoustic wave equations. The matrix-free framework of the {\ttfamily deal.II} library supports various multigrid strategies, including $h$-, $p$-, and $hp$-refinement across spatial and temporal dimensions. Extensive empirical evidence, provided through scaling and convergence tests on high-performance computing platforms, demonstrate high performance on perturbed meshes and problems with heterogeneous and discontinuous coefficients. Throughputs of over a billion degrees of freedom per second are achieved on problems with more than a trillion global degrees of freedom. The results prove that the space-time multigrid method can effectively solve complex problems in high-fidelity simulations and show great potential for use in coupled problems., Comment: 33 pages, 10 figures, 4 tables

We improve the performance of multigrid solvers on many-core architectures with cache hierarchies by reorganizing operations in the smoothing step to minimize memory transfers. We focus on patch smoothers, which offer robust convergence rates with respect to the finite element degree for various equations, in the setting of multiplicative subspace correction for numerical efficiency. By combining the computation of local residuals with local solvers, we increase the locality of the problem and thus reduce data transfers. The thread-parallel implementation of this algorithm is based on coloring, which contradicts cache efficiency. We improve data locality by rearranging the loop into batches so that more data can be reused. The organization of consecutive batches prioritizes data locality.

Unfitted finite element methods, like CutFEM, have traditionally been implemented in a matrix-based fashion, where a sparse matrix is assembled and later applied to vectors while solving the resulting linear system. With the goal of increasing performance and enabling algorithms with polynomial spaces of higher degrees, this contribution chooses a more abstract approach by matrix-free evaluation of the operator action on vectors instead. The proposed method loops over cells and locally evaluates the cell, face, and interface integrals, including the contributions from cut cells and the different means of stabilization. The main challenge is the efficient numerical evaluation of terms in the weak form with unstructured quadrature points arising from the unfitted discretization in cells cut by the interface. We present design choices and performance optimizations for tensor-product elements and demonstrate the performance by means of benchmarks and application examples. We demonstrate a speedup of more than one order of magnitude for the operator evaluation of a discontinuous Galerkin discretization with polynomial degree three compared to a sparse matrix-vector product and develop performance models to quantify the performance properties over a wide range of polynomial degrees.

Computational modeling of the melt pool dynamics in laser-based powder bed fusion metal additive manufacturing (PBF-LB/M) promises to shed light on fundamental mechanisms of defect generation. These processes are accompanied by rapid evaporation so that the evaporation-induced recoil pressure and cooling arise as major driving forces for fluid dynamics and temperature evolution. The magnitude of these interface fluxes depends exponentially on the melt pool surface temperature, which, therefore, has to be predicted with high accuracy. The present work utilizes a diffuse interface finite element model based on a continuum surface flux (CSF) description of interface fluxes to study dimensionally reduced thermal two-phase problems representative for PBF-LB/M in a finite element framework. It is demonstrated that the extreme temperature gradients combined with the high ratios of material properties between metal and ambient gas lead to significant errors in the interface temperatures and fluxes when classical CSF approaches, along with typical interface thicknesses and discretizations, are applied. It is expected that this finding is also relevant for other types of diffuse interface PBF-LB/M melt pool models. A novel parameter-scaled CSF approach is proposed, which is constructed to yield a smoother temperature field in the diffuse interface region, significantly increasing the solution accuracy. The interface thickness required to predict the temperature field with a given level of accuracy is less restrictive by at least one order of magnitude for the proposed parameter-scaled approach compared to classical CSF, drastically reducing computational costs. Finally, we showcase the general applicability of the parameter-scaled CSF to a 3D simulation of stationary laser melting of PBF-LB/M considering the fully coupled thermo-hydrodynamic multi-phase problem, including phase change.

We present accurate and mathematically consistent formulations of a diffuse-interface model for two-phase flow problems involving rapid evaporation. The model addresses challenges including discontinuities in the density field by several orders of magnitude, leading to high velocity and pressure jumps across the liquid-vapor interface, along with dynamically changing interface topologies. To this end, we integrate an incompressible Navier-Stokes solver combined with a conservative level-set formulation and a regularized, i.e., diffuse, representation of discontinuities into a matrix-free adaptive finite element framework. The achievements are three-fold: First, we propose mathematically consistent definitions for the level-set transport velocity in the diffuse interface region by extrapolating the velocity from the liquid or gas phase. They exhibit superior prediction accuracy for the evaporated mass and the resulting interface dynamics compared to a local velocity evaluation, especially for strongly curved interfaces. Second, we show that accurate prediction of the evaporation-induced pressure jump requires a consistent, namely a reciprocal, density interpolation across the interface, which satisfies local mass conservation. Third, the combination of diffuse interface models for evaporation with standard Stokes-type constitutive relations for viscous flows leads to significant pressure artifacts in the diffuse interface region. To mitigate these, we propose to introduce a correction term for such constitutive model types. Through selected analytical and numerical examples, the aforementioned properties are validated. The presented model promises new insights in simulation-based prediction of melt-vapor interactions in thermal multiphase flows such as in laser-based powder bed fusion of metals., Comment: 58 pages, 23 figures

This article proposes a novel high-performance computing approach for the prediction of the temperature field in powder bed fusion (PBF) additive manufacturing processes. In contrast to many existing approaches to part-scale simulations, the underlying computational model consistently resolves physical scan tracks without additional heat source scaling, agglomeration strategies or any other heuristic modeling assumptions. A growing, adaptively refined mesh accurately captures all details of the laser beam motion. Critically, the fine spatial resolution required for resolved scan tracks in combination with the high scan velocities underlying these processes mandates the use of comparatively small time steps to resolve the underlying physics. Explicit time integration schemes are well-suited for this setting, while unconditionally stable implicit time integration schemes are employed for the interlayer cool down phase governed by significantly larger time scales. These two schemes are combined and implemented in an efficient fast operator evaluation framework providing significant performance gains and optimization opportunities. The capabilities of the novel framework are demonstrated through realistic AM examples on the centimeter scale including the first scan-resolved simulation of the entire NIST AM Benchmark cantilever specimen, with a computation time of less than one day. Apart from physical insights gained through these simulation examples, also numerical aspects are thoroughly studied on basis of weak and strong parallel scaling tests. As potential applications, the proposed thermal PBF simulation framework can serve as a basis for microstructure and thermo-mechanical predictions on the part-scale, but also to assess the influence of scan pattern and part geometry on melt pool shape and temperature, which are important indicators for well-known process instabilities.

This work compares two Nitsche-type approaches to treat non-conforming triangulations for a high-order discontinuous Galerkin (DG) solver for the acoustic conservation equations. The first approach (point-to-point interpolation) uses inexact integration with quadrature points prescribed by a primary element. The second approach uses exact integration (mortaring) by choosing quadratures depending on the intersection between non-conforming elements. In literature, some excellent properties regarding performance and ease of implementation are reported for point-to-point interpolation. However, we show that this approach can not safely be used for DG discretizations of the acoustic conservation equations since, in our setting, it yields spurious oscillations that lead to instabilities. This work presents a test case in that we can observe the instabilities and shows that exact integration is required to maintain a stable method. Additionally, we provide a detailed analysis of the method with exact integration. We show optimal spatial convergence rates globally and in each mesh region separately. The method is constructed such that it can natively treat overlaps between elements. Finally, we highlight the benefits of non-conforming discretizations in acoustic computations by a numerical test case with different fluids.

We present an implementation of a fully stage-parallel preconditioner for Radau IIA type fully implicit Runge--Kutta methods, which approximates the inverse of $A_Q$ from the Butcher tableau by the lower triangular matrix resulting from an LU decomposition and diagonalizes the system with as many blocks as stages. For the transformed system, we employ a block preconditioner where each block is distributed and solved by a subgroup of processes in parallel. For combination of partial results, we either use a communication pattern resembling Cannon's algorithm or shared memory. A performance model and a large set of performance studies (including strong scaling runs with up to 150k processes on 3k compute nodes) conducted for a time-dependent heat problem, using matrix-free finite element methods, indicate that the stage-parallel implementation can reach higher throughputs when the block solvers operate at lower parallel efficiencies, which occurs near the scaling limit. Achievable speedup increases linearly with number of stages and are bounded by the number of stages. Furthermore, we show that the presented stage-parallel concepts are also applicable to the case that $A_Q$ is directly diagonalized, which requires complex arithmetic or the solution of two-by-two blocks and sequentializes parts of the algorithm. Alternatively to distributing stages and assigning them to distinct processes, we discuss the possibility of batching operations from different stages together.

This work investigates a variant of the conjugate gradient (CG) method and embeds it into the context of high-order finite-element schemes with fast matrix-free operator evaluation and cheap preconditioners like the matrix diagonal. Relying on a data-dependency analysis and appropriate enumeration of degrees of freedom, we interleave the vector updates and inner products in a CG iteration with the matrix-vector product with only minor organizational overhead. As a result, around 90% of the vector entries of the three active vectors of the CG method are transferred from slow RAM memory exactly once per iteration, with all additional access hitting fast cache memory. Node-level performance analyses and scaling studies on up to 147k cores show that the CG method with the proposed performance optimizations is around two times faster than a standard CG solver as well as optimized pipelined CG and s-step CG methods for large sizes that exceed processor caches, and provides similar performance near the strong scaling limit., Comment: 19 pages, 14 figures

This work studies three multigrid variants for matrix-free finite-element computations on locally refined meshes: geometric local smoothing, geometric global coarsening, and polynomial global coarsening. We have integrated the algorithms into the same framework-the open-source finite-element library deal.II-, which allows us to make fair comparisons regarding their implementation complexity, computational efficiency, and parallel scalability as well as to compare the measurements with theoretically derived performance models. Serial simulations and parallel weak and strong scaling on up to 147,456 CPU cores on 3,072 compute nodes are presented. The results obtained indicate that global coarsening algorithms show a better parallel behavior for comparable smoothers due to the better load balance particularly on the expensive fine levels. In the serial case, the costs of applying hanging-node constraints might be significant, leading to advantages of local smoothing, even though the number of solver iterations needed is slightly higher., Comment: 34 pages, 17 figures

Approximately $75 \%$ of the raw material and $50 \%$ of the products in the chemical industry are granular materials. The Discrete Element Method (DEM) provides detailed insights of phenomena at particle scale and it is therefore often used for modeling granular materials. However, because DEM tracks the motion and contact of individual particles separately, its computational cost increases non-linearly $O(n_p\log(n_p))$ -- $O(n_p^2)$ depending on the algorithm) with the number of particles ($n_p$). In this article, we introduce a new open-source parallel DEM software with load balancing: Lethe-DEM. Lethe-DEM, a module of Lethe, consists of solvers for two-dimensional and three-dimensional DEM simulations. Load-balancing allows Lethe-DEM to significantly increase the parallel efficiency by $\approx 25 - 70 \%$ depending on the granular simulation. We explain the fundamental modules of Lethe-DEM, its software architecture, and the governing equations. Furthermore, we verify Lethe-DEM with several tests including analytical solutions and comparison with other software. Comparisons with experiments in a flat-bottomed silo, wedge-shaped silo, and rotating drum validate Lethe-DEM. We investigate the strong and weak scaling of Lethe-DEM with $1 \leq n_c \leq 192$ and $32 \leq n_c \leq 320$ processes, respectively, with and without load-balancing. The strong-scaling analysis is performed on the wedge-shaped silo and rotating drum simulations, while for the weak-scaling analysis, we use a dam break simulation. The best scalability of Lethe-DEM is obtained in the range of $5000 \leq n_p/n_c \leq 15000$. Finally, we demonstrate that large scale simulations can be carried out with Lethe-DEM using the simulation of a three-dimensional cylindrical silo with $n_p=4.3 \times 10^6$ on 320 cores.

Providing evidence of finite-time singularities of the incompressible Euler equations in three space dimensions is still an unsolved problem. Likewise, the zeroth law of turbulence has not been proven to date by numerical experiments. We address this issue by high-resolution numerical simulations of the inviscid three-dimensional Taylor-Green vortex problem using a novel high-order discontinuous Galerkin discretization approach. Our main finding is that the kinetic energy evolution does not tend towards exact energy conservation for increasing spatial resolution of the numerical scheme, but instead converges to a solution with nonzero kinetic energy dissipation rate. This implies an energy dissipation anomaly in the absense of viscous dissipation according to Onsager's conjecture, and serves as an indication of finite-time singularities in incompressible inviscid flows. We demonstrate convergence to a dissipative solution for the three-dimensional inviscid Taylor-Green problem with a measured relative $L_2$-error of $0.27 \%$ for the temporal evolution of the kinetic energy and $3.52 \%$ for the kinetic energy dissipation rate.

This work presents the efficient, matrix-free finite-element library hyper.deal for solving partial differential equations in two to six dimensions with high-order discontinuous Galerkin methods. It builds upon the low-dimensional finite-element library deal.II to create complex low-dimensional meshes and to operate on them individually. These meshes are combined via a tensor product on the fly and the library provides new special-purpose highly optimized matrix-free functions exploiting domain decomposition as well as shared memory via MPI-3.0 features. Both node-level performance analyses and strong/weak-scaling studies on up to 147,456 CPU cores confirm the efficiency of the implementation. Results of the library hyper.deal are reported for high-dimensional advection problems and for the solution of the Vlasov--Poisson equation in up to 6D phase space., Comment: 33 pages, 18 figures

The present work develops hybrid multigrid methods for high-order discontinuous Galerkin discretizations of elliptic problems. Fast matrix-free operator evaluation on tensor product elements is used to devise a computationally efficient PDE solver. The multigrid hierarchy exploits all possibilities of geometric, polynomial, and algebraic coarsening, targeting engineering applications on complex geometries. Additionally, a transfer from discontinuous to continuous function spaces is performed within the multigrid hierarchy. This does not only further reduce the problem size of the coarse-grid problem, but also leads to a discretization most suitable for state-of-the-art algebraic multigrid methods applied as coarse-grid solver. The relevant design choices regarding the selection of optimal multigrid coarsening strategies among the various possibilities are discussed with the metric of computational costs as the driving force for algorithmic selections. We find that a transfer to a continuous function space at highest polynomial degree (or on the finest mesh), followed by polynomial and geometric coarsening, shows the best overall performance. The success of this particular multigrid strategy is due to a significant reduction in iteration counts as compared to a transfer from discontinuous to continuous function spaces at lowest polynomial degree (or on the coarsest mesh). The coarsening strategy with transfer to a continuous function space on the finest level leads to a multigrid algorithm that is robust with respect to the penalty parameter of the SIPG method. Detailed numerical investigations are conducted for a series of examples ranging from academic test cases to more complex, practically relevant geometries. Performance comparisons to state-of-the-art methods from the literature demonstrate the versatility and computational efficiency of the proposed multigrid algorithms.

This work proposes a basis for improved throughput of matrix-free evaluation of discontinuous Galerkin symmetric interior penalty discretizations on hexahedral elements. The basis relies on ideas of Hermite polynomials. It is used in a fully discontinuous setting not for higher order continuity but to minimize the effective stencil width, namely to limit the neighbor access of an element to one data point for the function value and one for the derivative. The basis is extended to higher orders with nodal contributions derived from roots of Jacobi polynomials and extended to multiple dimensions with tensor products, which enable the use of sum factorization. The beneficial effect of the reduced data access on modern processors is shown. Furthermore, the viability of the basis in the context of multigrid solvers is analyzed. While a plain point-Jacobi approach is less efficient than with the best nodal polynomials, a basis change via sum-factorization techniques enables the combination of the fast matrix-vector products with effective multigrid constituents. The basis change is essentially for free on modern hardware because these computations can be hidden behind the cost of the data access., Comment: 20 pages, 8 figures