Assessing the algorithmic scaling behavior of PDE solvers Matthias Bolten University of Wuppertal 2015/09/17
Motivation Stationary problems Time-dependent problems Wrap-up Outline Motivation Stationary problems Model problem Discretization Direct solvers Iterative methods Time-dependent problems Model problem Explicit methods Implicit methods Wrap-up Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 1/57
Motivation Stationary problems Time-dependent problems Wrap-up Outline Motivation Stationary problems Model problem Discretization Direct solvers Iterative methods Time-dependent problems Model problem Explicit methods Implicit methods Wrap-up Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 2/57
Motivation Stationary problems Time-dependent problems Wrap-up Motivation ◮ Many problems in scientific computing involve solution of a PDE ◮ A lot of optimized methods exist, providing good efficiency in terms of FLOPS ◮ Several scalable methods exist, providing good efficiency in terms of—well—efficiency ◮ A lot of effort is put on optimizing scientific applications ◮ The underlying method is rarely questioned ◮ But: Many methods’ algorithmic scalability does not scale ◮ Bad algorithmic scaling leads to bad weak scaling behavior Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 3/57
Motivation Stationary problems Time-dependent problems Wrap-up Motivation: Example Assume there is an application solving an elliptic partial differential equation in 3d. Now, consider two cases: 1. Solution of the PDE using Gaussian elimination at 80% peak 2. Solution of the PDE using full multigrid at 0.5% peak Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 4/57
Motivation Stationary problems Time-dependent problems Wrap-up Motivation: Example Assume there is an application solving an elliptic partial differential equation in 3d. Now, consider two cases: 1. Solution of the PDE using Gaussian elimination at 80% peak 2. Solution of the PDE using full multigrid at 0.5% peak Which is the better approach? Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 4/57
Motivation Stationary problems Time-dependent problems Wrap-up Motivation: Example Assume there is an application solving an elliptic partial differential equation in 3d. Now, consider two cases: 1. Solution of the PDE using Gaussian elimination at 80% peak 2. Solution of the PDE using full multigrid at 0.5% peak Which is the better approach? Answer: Method 1, as it provides more impressive performance numbers Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 4/57
Motivation Model problem Stationary problems Discretization Time-dependent problems Direct solvers Wrap-up Iterative methods Outline Motivation Stationary problems Model problem Discretization Direct solvers Iterative methods Time-dependent problems Model problem Explicit methods Implicit methods Wrap-up Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 5/57
Motivation Model problem Stationary problems Discretization Time-dependent problems Direct solvers Wrap-up Iterative methods Model problem: Poisson equation Let Ω = [0 , 1] 2 . Consider the Poisson equation given by − ∆ u ( x ) = f ( x ) , for x ∈ Ω and (1) u ( x ) = 0 for x ∈ ∂ Ω . This is the prototypical elliptic partial differential equation (PDE). A vast amount of other elliptic PDEs exist with ◮ higher dimensionality, ◮ more complex geometries, ◮ varying coefficients, ◮ additional (convective) terms, ◮ . . . . Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 6/57
Motivation Model problem Stationary problems Discretization Time-dependent problems Direct solvers Wrap-up Iterative methods Discretization: Finite differences I ◮ The simplest discretization is obtained by approximating the differential quotient d f ( x 0 + h ) − f ( x 0 ) ≈ f ( x 0 + h ) − f ( x 0 ) dxf ( x ) | x = x 0 = lim h h h → 0 ◮ Discretization on regular grid yields stencil in each grid point ◮ For the Laplace operator ∆ we obtain the stencil 1 1 1 − 4 1 h 2 1 Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 7/57
Motivation Model problem Stationary problems Discretization Time-dependent problems Direct solvers Wrap-up Iterative methods Discretization: Finite differences II ◮ Discretization on grid with N = ( n + 1) · ( n + 1) grid points using 5-point scheme yields linear system Lu = f , where − 1 h 2 (4 u i,j − u i − 1 ,j − u i +1 ,j − u i,j − 1 − u i,j +1 ) = f i,j , for i, j = 1 , . . . , n k , where h = 1 /n and u i,j = 0 for i ∈ { 0 , n + 1 } or j ∈ { 0 , n + 1 } . ◮ Eigenvalues and eigenvectors given by λ l,m = 4 − 2 cos( lπh ) − 2 cos( mπh ) , (2) ( ϕ l,m ) i,j = sin( lπih ) sin( mπjh ) , (3) Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 8/57
Motivation Model problem Stationary problems Discretization Time-dependent problems Direct solvers Wrap-up Iterative methods Discretization: Finite elements I Instead of (5) one can consider its weak formulation, i.e., find u vanishing on ∂ Ω s.t. � � − (∆ u ) v dx = f vdx Ω Ω for all test functions v . Using Green’s identity we obtain � � − ∇ u · ∇ v dx = f vdx. Ω Ω The functions u and v have to be chosen from proper function spaces. Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 9/57
Motivation Model problem Stationary problems Discretization Time-dependent problems Direct solvers Wrap-up Iterative methods Discretization: Finite elements II (Conforming) finite element methods use finite subspaces of these function spaces for both u and v . To obtain this subspace one proceeds as follows: The domain is decomposed into polyhedrons, often triangles are chosen, resulting in a triangulation like the following: Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 10/57
Motivation Model problem Stationary problems Discretization Time-dependent problems Direct solvers Wrap-up Iterative methods Discretization: Finite elements III The function space under consideration consists of all continuous functions vanishing on the boundary that restricted to the triangles are polynomials of a certain degree p , e.g., p = 1 . As basis we choose functions vanishing at all nodes but one, e.g.: Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 11/57
Motivation Model problem Stationary problems Discretization Time-dependent problems Direct solvers Wrap-up Iterative methods Discretization: Finite elements IV For a regular grid we can choose the following triangulation: Choosing piecewise linear basis functions yields the same stencil as for finite differences: 1 1 . 1 − 4 1 h 2 1 Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 12/57
Motivation Model problem Stationary problems Discretization Time-dependent problems Direct solvers Wrap-up Iterative methods Discretization: General considerations ◮ In either case, problem leads to linear system Ax = b, (4) where A ∈ C n × n , x, b ∈ C n ◮ Even exact solution x is only approximation to continuous problem, not a sampling of the continuous solution ◮ Error of exact solution is known as discretization error ◮ Error estimates exist ◮ Solution of (4) up to machine precision not required ◮ Additional algebraic error affordable Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 13/57
Motivation Model problem Stationary problems Discretization Time-dependent problems Direct solvers Wrap-up Iterative methods Direct solvers ◮ Usually based on factorization of system matrix A ◮ Well-known methods: ◮ A = LU : LU factorization (Gaussian elimination) ◮ A = LDL T : LDL factorization ◮ A = LL ∗ : Cholesky factorization (A hermitian, positive (semi-)definite) ◮ Direct solvers for sparse matrices exist, including heuristics to minimize fill-in ◮ Parallel implementations of the latter available ◮ Also include special solvers, e.g. using FFTs Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 14/57
Motivation Model problem Stationary problems Discretization Time-dependent problems Direct solvers Wrap-up Iterative methods Direct solvers ◮ Usually based on factorization of system matrix A ◮ Well-known methods: ◮ A = LU : LU factorization (Gaussian elimination) ◮ A = LDL T : LDL factorization ◮ A = LL ∗ : Cholesky factorization (A hermitian, positive (semi-)definite) ◮ Direct solvers for sparse matrices exist, including heuristics to minimize fill-in ◮ Parallel implementations of the latter available ◮ Also include special solvers, e.g. using FFTs ◮ Provide exact solution of (4) (up to machine precision, subject to round-off) Matthias Bolten, Assessing the algorithmic scaling behavior of PDE solvers 14/57
Recommend
More recommend
