Scalable domain decomposition solvers for stochastic PDEs in high performance computing

A Domain Decomposition (DD) method for solving time-harmonic aeroacoustic problems is presented. The computational domain is decomposed into subdomains, and the aerodynamic-aeroacoustic boundary-value problem is solved independently for each subdomain. On the artificially introduced subdomain boundaries impedance-type transmission conditions are imposed, which assure the uniqueness of the solution. The Finite Element method is used to discretize the subdomain problems. A direct solver (LU factorization) is used to solve the subdomain systems of equations and an iterative scheme updates the transmission conditions to recover the global solution. This DD method is applied to solve the classical problem of a thin airfoil in a transverse gust. The problem is formulated in terms of the unsteady pressure, and exact non-reflecting conditions are imposed on the outer computational boundary. The accuracy of the numerical scheme is validated by comparison with existing solutions for both the near-field unsteady pressure and the farfiled radiated sound. The convergence as well as the computational time and memory requirements of the present method are investigated. It is shown that by combining the subdomain direct solvers with global iterations, this DD method significantly reduces both the computational time and the memory requirements. The present algorithm is particularly suited for parallel computing.

A domain decomposition method with Lagrange multipliers for the Stokes problem is developed and analysed. A common approach to solve the Stokes problem, termed the Uzawa algorithm, is to decouple the velocity and the pressure. This approach yields the Schur complement system for the pressure Lagrange multiplier which is solved with an iterative solver. Each outer iteration of the Uzawa procedure involves the inversion of a Laplacian in each spatial direction. The objective of this paper is to effectively solve this inner system (the vector Laplacian system) by applying the finite‐element tearing and interconnecting (FETI) method. Previously calculated search directions for the FETI solver are reused in subsequent outer Uzawa iterations. The advantage of the approach proposed in this paper is that pressure is continuous across the entire computational domain. Numerical tests are performed by solving the driven cavity problem. An analysis of the number of outer Uzawa iterations and inner FETI iterations is reported. Results show that the total number of inner iterations is almost numerically scalable since it grows asymptotically with the mesh size and the number of subdomains. Copyright © 2003 John Wiley & Sons, Ltd.

SummaryThis article deals with the computation of guaranteed lower bounds of the error in the framework of finite element and domain decomposition methods. In addition to a fully parallel computation, the proposed lower bounds separate the algebraic error (due to the use of a domain decomposition iterative solver) from the discretization error (due to the finite element), which enables the steering of the iterative solver by the discretization error. These lower bounds are also used to improve the goal‐oriented error estimation in a substructured context. Assessments on 2D static linear mechanic problems illustrate the relevance of the separation of sources of error and the lower bounds' independence from the substructuring. We also steer the iterative solver by an objective of precision on a quantity of interest. This strategy consists in a sequence of solvings and takes advantage of adaptive remeshing and recycling of search directions. Copyright © 2016 John Wiley & Sons, Ltd.

This article presents an analysis of deterministic models of computer virus epidemic propagation (SIR models) and their classification. The main areas of research into these models are highlighted. An analysis of existing stochastic models based on the SIR model and their diversity are presented. A method for constructing a stochastic SIR model based on the classical SIR model, represented by a system of Ito stochastic differential equations with a Wiener process, is proposed. The possibility of constructing a stochastic and deterministic model of an information system protected against computer virus infection with probability 1 is demonstrated: a stochastic model, in which infection by viruses occurs continuously, and a deterministic model, in which the virus is present in the information system. A mathematical stochastic model of an information system protected from computer virus outbreaks is constructed as a system of stochastic differential equations whose first integrals are invariants preserved with probability 1. A certain functional relationship between model variables, maintaining a constant value, is considered as the system's security indicator. Introducing a compensator (program control with probability 1 (PCP1)) into the model allows the specified security indicator, described by the model variables, to be maintained with probability 1. Introducing a compensator (program control with probability 1 (PCP1)) into the model allows the specified security indicator, described by the model variables, to be maintained with probability 1. Similarly, based on the proposed algorithm, a deterministic model of an information system protected from computer virus infection is constructed. A control similar to programmed control with probability 1 (PCP1) is introduced into the constructed model, which allows the invariants to be maintained. A distinctive feature of the proposed models is that they preserve invariants associated with the properties that ensure the security of the information system. The behavior of the constructed models is studied using numerical simulation in the MathCad environment. Based on the research results, conclusions were drawn on the possibility of using the proposed method in constructing stochastic models based on other models of epidemic spread, as well as for models of protecting an information system from the spread of a computer virus epidemic.

For the self-developed three-dimensional whole-core High-fidelity NEutron Transport calculation program HNET, although the numerical acceleration algorithms can improve the computation performance of high-fidelity neutron transport in terms of algorithms and models, it still faces some critical issues such as long computing time and enormous memory requirement. The rapid development of high-performance clusters provides a foundation for the application of massively parallel computing. Most current MOC programs are based on a single-type variable to achieve efficient parallelism, and spatial domain decomposition methods are the most common parallel schemes. However, its parallelism is limited and cannot fully utilize the current state-of-the-art computer resources. To solve this problem, hybrid parallel strategies are implemented in HNET to further expand the parallel degree, improve the speed of computation, and reduce the memory requirement. A hybrid MPI/OpenMP method based on domain decomposition and characteristic rays is proposed for the method of characteristics (MOC). For domain decomposition, the simulation domain is divided into spatial subdomains, with each subdomain handled by different processes. On this basis, the characteristic ray parallelism is implemented taking advantage of the inherent parallelism of the characteristic rays. Meanwhile, the optimization of the hybrid parallel strategy further improves the computation speed by eliminating atomic operations and using private pointer arrays and other techniques. In addition, in the framework of generalized equivalence theory (GET) based two-level coarse mesh finite difference method (CMFD), there is also a certain time consumption in solving CMFD linear system. Hence, for CMFD, a hybrid MPI/OpenMP method based on domain decomposition and secondary domain decomposition can be used. Using the secondary domain decomposition method, each subdomain is divided into sub-subdomains, with each sub-subdomain handled by threads, which allows CMFD to utilize the resources of characteristic ray parallelism in MOC, and also increases the speed of CMFD computation. Numerical results show that for both steady-state and transient calculations, the hybrid MPI/OpenMP in HNET can further expand the parallelism and accelerate the computation. It can take full use of parallel resources and achieve large-scale parallelism.

A numerical implementation of the Domain Decomposition (DD) method for parallel computation of time-harmonic waves is presented. The Finite Element method on unstructured mesh is employed to discretize the boundary-value problem. A nonover-lapping mesh partitioning is used to build an overlapping domain decomposition. The algorithm also defines the new subdomain boundaries, other than the original boundaries of the global computational domain, on which transmission conditions are introduced. The independent subdomain problems can be solved concurrently, one problem per processor, on a parallel computer. The iterative scheme used to update the subdomain transmission conditions is computationally efficient and only swaps the nodal values in the overlapping region between adjacent subdomains. For exterior problems, the communication required to update the conditions on the outer boundary is substantially reduced by employing a modified Dirichlet-to-Neumann (DtN) map with truncated support. Numerical results are presented for a radiation model problem. It is shown that the DtN nonreflecting conditions eliminate spurious reflections and that the iterative scheme leads to a continuous global solution. The DD algorithm also has a very good convergence rate and significantly reduces memory requirement in comparison with a direct solver.

Iterative solution techniques for finite element applications

Integrating batched sparse iterative solvers for the collision operator in fusion plasma simulations on GPUs

Unconditional stability of the Crank-Nicolson Finite Difference Time Domain (CN-FDTD) method permits us to use time steps over the Courant-Friedrich-Lewy (CFL) limit of conventional FDTD method. However, in this work it was realized that, when the time step is set above CFL limit the coefficient matrix arising from Crank-Nicolson method is no longer diagonally dominant and iterative solvers require longer solution time in each FDTD iteration. Frequency dependent CN-FDTD (FD-CN-FDTD) scheme for Debye media is formulated and numerical tests are performed with two widely used sparse iterative solvers, Bi-Conjugate Gradient Stabilised (BiCGStab) and Generalised Minimal Residual (GMRES), for comparison. BiCGStab outperforms GMRES in every aspect of the study.

Batched linear solvers, which solve many small related but independent problems, are important in several applications. This is increasingly the case for highly parallel processors such as graphics processing units (GPUs), which need a substantial amount of work to keep them operating efficiently and solving smaller problems one-by-one is not an option. Because of the small size of each problem, the task of coming up with a parallel partitioning scheme and mapping the problem to hardware is not trivial. In recent history, significant attention has been given to batched dense linear algebra. However, there is also an interest in utilizing sparse iterative solvers in a batched form, and this presents further challenges. An example use case is found in a gyrokinetic Particle-In-Cell (PIC) code used for modeling magnetically confined fusion plasma devices. The collision operator has been identified as a bottleneck, and a proxy app has been created for facilitating optimizations and porting to GPUs. The current collision kernel linear solver does not run on the GPU-a major bottleneck. As these matrices are well-conditioned, batched iterative sparse solvers are an attractive option. A batched sparse iterative solver capability has recently been developed in the Ginkgo library. In this paper, we describe how the software architecture can be used to develop an efficient solution for the XGC collision proxy app. Comparisons for the solve times on NVIDIA V100 and A100 GPUs and AMD MI100 GPUs with one dual-socket Intel Xeon Skylake CPU node with 40 OpenMP threads are presented for matrices representative of those required in the collision kernel of XGC. The results suggest that GINKGO's batched sparse iterative solvers are well suited for efficient utilization of the GPU for this problem, and the performance portability of Ginkgo in conjunction with Kokkos (used within XGC as the heterogeneous programming model) allows seamless execution for exascale oriented heterogeneous architectures at the various leadership supercomputing facilities.

This paper describes the Domain Decomposition Method (DDM) for high frequency electromagnetic field problems. As an efficient iterative method to solve a complex symmetric system arising from the formulation of the E method, a DDM based on a MINimal RESidual (MINRES)-like iterative procedure is presented. Numerical experiments are demonstrated to report a faster convergence compared to conventional methods.

A new low-memory modification of the density fitting approximation based on a combination of a continuous fast multipole method (CFMM) and a preconditioned conjugate gradient solver is presented. Iterative conjugate gradient solver uses preconditioners formed from blocks of the Coulomb metric matrix that decrease the number of iterations needed for convergence by up to one order of magnitude. The matrix-vector products needed within the iterative algorithm are calculated using CFMM, which evaluates them with the linear scaling memory requirements only. Compared with the standard density fitting implementation, up to 15-fold reduction of the memory requirements is achieved for the most efficient preconditioner at a cost of only 25% increase in computational time. The potential of the method is demonstrated by performing density functional theory calculations for zeolite fragment with 2592 atoms and 121,248 auxiliary basis functions on a single 12-core CPU workstation.

A comparison is made between the probabilistic domain decomposition (DD) method and a certain deterministic DD method for solving linear elliptic boundary-value problems. Since in the deterministic approach the CPU time is affected by intercommunications among the processors, it turns out that the probabilistic method performs better, especially when the number of subdomains (hence, of processors) is increased. This fact is clearly illustrated by some examples. The probabilistic DD algorithm has been implemented in an MPI environment, in order to exploit distributed computer architectures. Scalability and fault-tolerance of the probabilistic DD algorithm are emphasized.

A non-overlapping Schwarz domain decomposition method with high-order finite elements for flow acoustics

A non-symmetric interface Schur complement arises from non- selfcon-jugate second order elliptic problems with domain decomposition methods. The usual numerical methods for solving it are GMRES, OR-THOMIN, and BICGSTAB, but they take a large amount of computer time and memory. The authors find in this paper that the nonsymmetric Schur complement can in fact be changed into a symmetric one by scaling. Then an efficient preconditioner can be provided by which the preconditioned system can be solved iteratively by a modified PCG method. When the problem is imposed on a rectangular region, the condition number is estimated and is nearly one. Numerical experiments are also presented. Non-selfconjugate problems arise in mathematical modeling and numerical simulation of fluid flows and transport in porous media.