Matrix-oriented FEM formulation for reaction-diffusion PDEs on a large class of 2D domains

For the spatial discretization of elliptic and parabolic partial differential equations (PDEs), we provide a Matrix-Oriented formulation of the classical Finite Element For the spatial discretization of elliptic and parabolic partial differential equations (PDEs), we provide a Matrix-Oriented formulation of the classical Finite Element Method, called MO-FEM, of arbitrary order $k\in\mathbb{N}$. On a quite general class of 2D domains, namely \emph{separable domains}, and even on special surfaces, the discrete problem is then reformulated as a \emph{multiterm Sylvester matrix equation} where the additional terms account for the geometric contribution of the domain shape. By considering the IMEX Euler method for the PDE time discretization, we obtain a sequence of these equations. To solve each multiterm Sylvester equation, we apply the matrix-oriented form of the Preconditioned Conjugate Gradient (MO-PCG) method with a matrix-oriented preconditioner that captures the spectral properties of the Sylvester operator. Solving the Poisson problem and the heat equation on some separable domains by MO-FEM-PCG, we show a gain in computational time and memory occupation wrt the classical vector PCG with same preconditioning or wrt a LU based direct method. As an application, we show the advantages of the MO-FEM-PCG to approximate Turing patterns on some separable domains and cylindrical surfaces for a morphochemical reaction-diffusion PDE system for battery modeling.