  Hierarchical $\mathcal{H}^2$-matrices are asymptotically optimal
representations for the discretizations of non-local operators such as those
arising in integral equations or from kernel functions. Their $O(N)$ complexity
in both memory and operator application makes them particularly suited for
large-scale problems. As a result, there is a need for software that provides
support for distributed operations on these matrices to allow large-scale
problems to be represented. In this paper, we present high-performance,
distributed-memory GPU-accelerated algorithms and implementations for
matrix-vector multiplication and matrix recompression of hierarchical matrices
in the $\mathcal{H}^2$ format.
  The algorithms are a new module of H2Opus, a performance-oriented package
that supports a broad variety of $\mathcal{H}^2$-matrix operations on CPUs and
GPUs. Performance in the distributed GPU setting is achieved by marshaling the
tree data of the hierarchical matrix representation to allow batched kernels to
be executed on the individual GPUs. MPI is used for inter-process
communication. We optimize the communication data volume and hide much of the
communication cost with local compute phases of the algorithms. Results show
near-ideal scalability up to 1024 NVIDIA V100 GPUs on Summit, with performance
exceeding 2.3 Tflop/s/GPU for the matrix-vector multiplication, and 670
Gflops/s/GPU for matrix compression, which involves batched QR and SVD
operations.
  We illustrate the flexibility and efficiency of the library by solving a 2D
variable diffusivity integral fractional diffusion problem with an algebraic
multigrid-preconditioned Krylov solver and demonstrate scalability up to 16M
degrees of freedom problems on 64 GPUs.
