  High fidelity scientific simulations modeling physical phenomena typically
require solving large linear systems of equations which result from
discretization of a partial differential equation (PDE) by some numerical
method. This step often takes a vast amount of computational time to complete,
and therefore presents a bottleneck in simulation work. Solving these linear
systems efficiently requires the use of massively parallel hardware with high
computational throughput, as well as the development of algorithms which
respect the memory hierarchy of these hardware architectures to achieve high
memory bandwidth.
  In this paper, we present an algorithm to accelerate Jacobi iteration for
solving structured problems on graphics processing units (GPUs) using a
hierarchical approach in which multiple iterations are performed within on-chip
shared memory every cycle. A domain decomposition style procedure is adopted in
which the problem domain is partitioned into subdomains whose data is copied to
the shared memory of each GPU block. Jacobi iterations are performed internally
within each block's shared memory, avoiding the need to perform expensive
global memory accesses every step. We test our algorithm on the linear systems
arising from discretization of Poisson's equation in 1D and 2D, and observe
speedup in convergence using our shared memory approach compared to a
traditional Jacobi implementation which only uses global memory on the GPU. We
observe a x8 speedup in convergence in the 1D problem and a nearly x6 speedup
in the 2D case from the use of shared memory compared to a conventional GPU
approach.
