CS 267 Applications of Parallel Computers Solving Linear Systems Arising from PDEsOutlineRecap of “Sources of Parallelism” LectureContinuation of RecapSolving PDEsPoisson’s equation arises in many modelsRelation of Poisson’s to Gravity, ElectrostaticsPoisson’s equation in 1D: 2u/x2 = f(x)2D Poisson’s equation: 2u/x2 + 2u/y2 = f(x,y)Algorithms for 2D Poisson with N UnknownsShort explanations of AlgorithmsComments on practical meshesJacobi’s Method on Regular 2D MeshConvergence of Nearest Neighbor MethodsSlide 21Parallelizing Jacobi’s MethodLocality Optimization in JacobiRedundant Ghost Nodes in JacobiSuccessive Overrelaxation (SOR)Gauss-SeidelSlide 27PowerPoint PresentationSlide 29Conjugate Gradient (CG) for solving A*x = bSummary of Jacobi, SOR and CGSolving the Poisson equation with the FFTSerial FFTUsing the 1D FFT for filteringUsing the 2D FFT for image compressionRelated TransformsSerial Algorithm for the FFTDivide and Conquer FFTDivide-and-Conquer FFTAn Iterative AlgorithmParallel 1D FFTBlock Layout of 1D FFTCyclic Layout of 1D FFTParallel ComplexityFFT With “Transpose”Is the Communication Step a Transpose?Complexity of the FFT with TransposeComment on the 1D Parallel FFTHigher Dimension FFTsFFTW – Fastest Fourier Transform in the West01/14/19 CS267 Lecture 15 1CS 267Applications of Parallel ComputersSolving Linear Systems Arising from PDEsKathy Yelickwww.cs.berkeley.edu/~yelick/cs26701/14/19 CS267 Lecture 15 2Outline•Review Poisson equation•Overview of Methods for Poisson Equation•Jacobi’s method•Red-Black SOR method•Conjugate Gradient•FFT•Multigrid (later lecture this semester)•Reduce to sparse-matrix-vector multiply•Need them to understand Multigrid01/14/19 CS267 Lecture 15 3Recap of “Sources of Parallelism” Lecture•Discrete event systems:•Examples: “Game of Life,” logic level circuit simulation. •Particle systems:•Examples: billiard balls, semiconductor device simulation, galaxies.•Lumped variables depending on continuous parameters:•ODEs, e.g., circuit simulation (Spice), structural mechanics, chemical kinetics.•Continuous variables depending on continuous parameters:•PDEs, e.g., heat, elasticity, electrostatics.•A given phenomenon can be modeled at multiple levels.•Many simulations combine more than one of these techniques.01/14/19 CS267 Lecture 15 4Continuation of RecapExamples of PDE systems include•Hyperbolic problems (waves):•Quantum mechanics: Wave-function(position,time)•Elliptic (steady state) problems:•Electrostatic or Gravitational Potential: Potential(position)•Parabolic (time-dependent) problems:•Heat flow: Temperature(position, time)•Diffusion: Concentration(position, time)Many problems combine features of above•Fluid flow: Velocity,Pressure,Density(position,time)•Elasticity: Stress,Strain(position,time)01/14/19 CS267 Lecture 15 6Solving PDEsSolution strategies:•Hyperbolic problems (waves):•Use explicit time-stepping •Solution at each point depends on neighbors at previous time•Elliptic (steady state) problems:•Everything depends on everything else•This means locality is harder to find than in hyperbolic problems•Parabolic (time-dependent) problems:•Involves an elliptic solve at each time-step•Focus on elliptic problems•Canonical example is the Poisson equation2u/x2 + 2u/y2 + 2u/z2 = f(x,y,z)01/14/19 CS267 Lecture 15 7Poisson’s equation arises in many models•Heat flow: Temperature(position, time)•Diffusion: Concentration(position, time)•Electrostatic or Gravitational Potential: Potential(position)•Fluid flow: Velocity,Pressure,Density(position,time)•Quantum mechanics: Wave-function(position,time)•Elasticity: Stress,Strain(position,time)•Variations of pure Poisson involve other coefficients3D: 2u/x2 + 2u/y2 + 2u/z2 = f(x,y,z)2D: 2u/x2 + 2u/y2 = f(x,y)1D: 2u/x2 = f(x)f represents the sources and boundary conditions01/14/19 CS267 Lecture 15 8Relation of Poisson’s to Gravity, Electrostatics•Force on particle at (x,y,z) due to particle at 0 is -(x,y,z)/r^3, where r = sqrt(x +y +z )•Force is also gradient of potential V(x,y,z) = -1/r Force = -(d/dx V, d/dy V, d/dz V) = -grad V•V(x,y,z) satisfies Poisson’s equation (try it!)2 2 201/14/19 CS267 Lecture 15 9Poisson’s equation in 1D: 2u/x2 = f(x)2 -1 -1 2 -1 -1 2 -1 -1 2 -1 -1 2T =2-1 -1Graph and “stencil”Discretize 2u/x2 = f(x) on regular mesh ui = u(i*h) to get [ u i+1 – 2*u i + u i-1 ] / h2 = f(x)Write as solving Tu = -h2 * f for u where01/14/19 CS267 Lecture 15 102D Poisson’s equation: 2u/x2 + 2u/y2 = f(x,y)•Similar to the 1D case, but the matrix T is now•Grid points numbered left to right, top row to bottom row•3D is analogous•Similar “adjacency matrix” for arbitrary graph4 -1 -1-1 4 -1 -1 -1 4 -1 -1 4 -1 -1 -1 -1 4 -1 -1 -1 -1 4 -1 -1 4 -1 -1 -1 4 -1 -1 -1 4T =4-1-1-1-1Graph and “stencil”01/14/19 CS267 Lecture 15 11Algorithms for 2D Poisson with N UnknownsAlgorithm Serial PRAM Memory #Procs•Dense LU N3N N2N2•Band LU N2N N3/2N•Jacobi N2 N N N•Explicit Inv. N2 log N N2N2•Conj.Grad. N 3/2N 1/2 *log N N N•RB SOR N 3/2N 1/2N N•Sparse LU N 3/2N 1/2N*log N N•FFT N*log N log N N N•Multigrid N log2 N N N•Lower bound N log N NPRAM is an idealized parallel model with zero cost communication01/14/19 CS267 Lecture 15 12Short explanations of Algorithms•Sorted in two orders (roughly):•from slowest to fastest on sequential machines•from most general (any matrix) to most specialized (matrices “like” Poisson)•Dense LU: Gaussian elimination; works on any N-by-N matrix•Band LU: exploit fact that T is nonzero only on sqrt(N) diagonals nearest main diagonal, so faster•Jacobi: essentially does matrix-vector multiply by T in inner loop of iterative algorithm. (May be used in Multigrid)•Explicit Inverse: assume we want to solve many systems with T, we can pre-compute and store inv(T) “amortized”, and just multiply by it.
View Full Document