# Efficient Discretization of the Laplacian: Application to Moving Boundary Problems
## 1 Introduction  
The numerical solution of PDEs is one of the most important problems in computational mathematics. In real-life cases, these equations serve as a mathematical model of different phenomena. To simulate them, the corresponding equations have to be solved. For this, usually, we have to perform a numerical solution procedure. From a modeling perspective, the primary challenge in these simulations lies in defining and incorporating realistic boundary conditions. At the same time, the basis of any numerical solution is an appropriate discretization. Typically, this is applied to both spatial and temporal variables, provided that we deal with time-dependent problems. Our focus will be spatial discretization. The computational methods include mainly the finite element [1], finite difference (FD) [2], and spectral or Fourier-type methods [3]. Alternatively, in the case of simple linear differential operators, a family of mesh-free methods (e.g., the method of fundamental solutions [4; 5]) can also be applied. However, for complex spatial domains, which is common in realistic scenarios, both Fourier-type and traditional finite difference discretizations often fall short. Accordingly, several attempts were made to extend finite difference methods for non-uniform grids, see [6; 7]. At the same time, in the engineering practice, for the numerical solution of PDEs, usually finite element-based discretizations are used. More precisely, one should call them Galerkin-type discretizations, which include the non-conforming or the discontinuous Galerkin-type discretizations. However, a significant drawback of this approach is its time-consuming nature. Constructing a finite element discretization involves generating acomprehensive data structure that includes a hierarchical organization of vertices, edges, and faces, along with their interrelations. Moreover, within such a framework, one should assemble all data into a large algebraic equation. In the case of (local) linear differential operators, we obtain a large linear algebraic system with a sparse matrix. Here, the full discretization (i.e., the construction of the above matrix) can take more time compared to the solution of the resulting system. This can extremely slow down the simulation procedure, especially when the computational domain changes dynamically over time, such as in moving boundary problems. Addressing these issues remains a significant challenge in computational mathematics. A main source of the enhanced computing time is that the above finite element assembling procedure can hardly be vectorized. Accordingly, to speed up the assembling, several procedures were proposed, which are also called sequencing. In [8], the authors introduced a vectorized algorithm for this purpose, which was then developed further in [9] in another way to reduce computing time. Another conventional way to speed up finite element assembling is to perform parallel processing and computations on GPUs. Regarding this, several attempts were also made, see, e.g., [10; 11]. Whether all of these results led indeed to a significant speed-up or not, they all need some extra effort. In each case, one has to build up an extra structure to support the algorithm. In concrete terms, usually, the authors employ a kind of graph coloring algorithm to partition elements into non-overlapping sets. In this way, the elemental matrices can be assembled into a global finite element matrix in parallel without conflicts. In such procedures, a suitable sparse format of the global matrix is needed [11]. The aim of the present work is to propose a vectorized algorithm for assembling the discretization matrix of the Laplace operator in general domains. The efficient (accurate and fast) discretization of the Laplacian operator has great importance in a number of physical models.