Search

The development of modern computing architectures with ever-increasing amounts of parallelism has allowed for the solution of previously intractable problems across a variety of scientific disciplines. Despite these advances, multiscale computing problems continue to pose an incredible challenge to modern architectures because they require resolving scales that often vary by orders of magnitude in both space and time. Such complications have led us to consider alternative discretizations for partial differential equations (PDEs) which use expansions involving integral operators to approximate spatial derivatives (Christlieb et al. in J Comput Phys 379:214–236, 2019; Christlieb et al. J Sci Comput 82:52(3):1–29, 2020; Christlieb et al. J Comput Phys 415:1–25, 2020). These constructions use explicit information within the integral terms, but treat boundary data implicitly, which contributes to the overall speed of the method. This approach is provably unconditionally stable for linear problems and stability has been demonstrated experimentally for nonlinear problems. Additionally, it is matrix-free in the sense that it is not necessary to invert linear systems and iteration is not required for nonlinear terms. Moreover, the scheme employs a fast summation algorithm that yields a method with a computational complexity of $${\mathcal {O}}(N)$$ , where N is the number of mesh points along a coordinate direction. While much work has been done to explore the theory behind these methods, their practicality in large scale computing environments is a largely unexplored topic. In this work, we explore the performance of these methods by developing a domain decomposition algorithm suitable for distributed memory systems along with shared memory algorithms. As a first pass, we derive an artificial Courant–Friedrichs–Lewy condition that enforces a nearest-neighbor (N-N) communication pattern and briefly discuss possible generalizations. We also analyze several approaches for implementing the parallel algorithms by optimizing predominant loop structures and maximizing data reuse. Using a hybrid design that employs MPI and Kokkos (Edwards and Trott in J Parallel Distrib Comput 74:3202–3216, 2014) for the distributed and shared memory components of the algorithms, respectively, we show that our methods are efficient and can sustain an update rate $$> 1\times 10^8$$ DOF/node/s. We provide results that demonstrate the scalability and versatility of our algorithms using several different PDE test problems, including a nonlinear example, which employs an adaptive time-stepping rule.

In \cite{christlieb2019kernel}, the authors developed a class of high-order numerical schemes for the Hamilton-Jacobi (H-J) equations, which are unconditionally stable, yet take the form of an explicit scheme. This paper extends such schemes, so that they are more effective at capturing sharp gradients, especially on nonuniform meshes. In particular, we modify the weighted essentially non-oscillatory (WENO) methodology in the previously developed schemes by incorporating an exponential basis and adapting the previously developed nonlinear filters used to control oscillations. The main advantages of the proposed schemes are their effectiveness and simplicity, since they can be easily implemented on higher-dimensional nonuniform meshes. We perform numerical experiments on a collection of examples, including H-J equations with linear, nonlinear, convex and non-convex Hamiltonians. To demonstrate the flexibility of the proposed schemes, we also include test problems defined on non-trivial geometry.