Ahead of Publication / Just Accepted

The aim of this paper is to propose an efficient spectral-Galerkin method for the numerical approximation of the time-space fractional diffusion equation in an unbounded domain by using the fractional-order generalized Jacobi functions and the mapped Chebyshev functions, and theoretically prove the high convergence rate. The reason for using the fractional-order generalized Jacobi functions is to approximate the solution with singularity at the initial time. We prove that the proposed approximation scheme has a spectral convergence rate when the solution of a given problem satisfies a particular condition. We also establish the stability of the method.

We study the popular modularity matrix and respective functional used in connection with graph clustering and derive some properties useful when performing vertex aggregation of the associated graph. These properties are employed in the derivation of a multilevel parallel pairwise aggregation algorithm. Comparative performance results of the studied algorithm applied to graph clustering tested against the popular Louvain algorithm are presented. Some illustrative examples show that the resulting aggregates if used in an adaptive algebraic multigrid (AMG) are able to follow strong direction of anisotropy in finite element problems.

In this paper we formulate and analyze adaptive (space-time) least-squares finite element methods for the solution of convection-diffusion equations. The convective derivative 𝒗 ⋅ ∇ ⁡ u {\boldsymbol{v}\cdot\nabla u} is considered as part of the total time derivative d d ⁢ t ⁢ u = ∂ t ⁡ u + 𝒗 ⋅ ∇ ⁡ u {\frac{d}{dt}u=\partial_{t}u+\boldsymbol{v}\cdot\nabla u} , and therefore we can use a rather standard stability and error analysis for related space-time finite element methods. For stationary problems we restrict the ansatz space H 0 1 ⁢ ( Ω ) {H^{1}_{0}(\Omega)} such that the convective derivative is considered as an element of the dual H - 1 ⁢ ( Ω ) {H^{-1}(\Omega)} of the test space H 0 1 ⁢ ( Ω ) {H^{1}_{0}(\Omega)} , which also allows unbounded velocities 𝒗 {\boldsymbol{v}} . While the discrete finite element schemes are always unique solvable, the numerical solutions may suffer from a bad approximation property of the finite element space when considering convection dominated problems, i.e., small diffusion coefficients. Instead of adding suitable stabilization terms, we aim to resolve the solutions by using adaptive (space-time) finite element methods. For this we introduce a least-squares approach where the discrete adjoint defines local a posteriori error indicators to drive an adaptive scheme. Numerical examples illustrate the theoretical considerations.

In this work, we present a multicontinuum homogenization approach for a two-phase model with linear transport, governed by coupled flow and transport equations with high-contrast coefficients. The method introduces multiple macroscopic variables to represent distinct physical characteristics of the solution, such as local averages over subregions. To capture fine-scale heterogeneities, multiscale basis functions are constructed by solving local constrained energy minimization problems. The resulting reduced-order model is defined on a coarse scale and involves the macroscopic variables and homogenized coefficients. It retains the essential features of the original fine-scale system while reducing computational costs. Numerical examples are provided to validate the proposed model. The results demonstrate that the solutions obtained from the proposed model accurately approximate those from the original model.

Computational technologies for the approximate solution of multidimensional boundary value problems often rely on irregular computational meshes and finite-volume approximations. In this framework, the discrete problem represents the corresponding conservation law for control volumes associated with the nodes of the mesh. This approach is most naturally and consistently implemented using Delaunay triangulations together with Voronoi diagrams as control volumes. In this paper, we employ meshes with nodes located both at the vertices of Delaunay triangulations and at the generators of Voronoi partitions. The cells of the merged Voronoi–Delaunay mesh are orthodiagonal quadrilaterals. On such meshes, scalar and vector functions, as well as invariant gradient and divergence operators of vector calculus, can be conveniently approximated. We illustrate the capabilities of this approach by solving a steady–state diffusion–reaction problem in an anisotropic medium.

In this work, we present a semi-implicit structure-preserving nonstaggered central scheme (SP-NCS) that is both well-balanced and positivity-preserving, designed to solve the two-dimensional shallow water Exner equations with time-dependent bottom topography. The SP-NCS can be viewed as a two-dimensional extension of the numerical method presented in [D. Li and J. Dong, A robust hybrid unstaggered central and Godunov-type scheme for Saint-Venant–Exner equations with wet/dry fronts, Comput. & Fluids 235 2022, Article ID 105284], which is uncoupled for the shallow water Exner equations. The SP-NCS is also inspired by the nonstaggered central scheme presented in [J. Dong and X. Qian, Structure-preserving nonstaggered central difference schemes at wet-dry fronts for the shallow water equations, Commun. Appl. Math. Comput. 8 2026, 1, 366–410], but their method is significantly adapted to preserve the still-water steady state when the computational domain contains wetting and drying transitions. Retaining the stationary solution in the two-dimensional case, especially when the domain include wetting and drying transitions, presents a significant challenge. This is because the backward step and the discretization of the source term in the corrector step become more complex when dealing with such fronts. A key innovation is the introduction of a structure-preserving parameter that helps retain the stationary solution even in the presence of wet-dry fronts. This task is nontrivial due to the robust requirements involved in ensuring the stability and accuracy of the solution. In particular, the SP-NCS demonstrates robustness in handling the relatively strong interactions between the two coupled models even with a large time step. This is particularly relevant when the rate of change in the bed topography is much slower than the speed of the water surface waves. We rigorously prove the positivity-preserving and well-balanced properties of the SP-NCS. Finally, we conduct several numerical experiments to validate the theoretical results and demonstrate the effectiveness of the proposed method in preserving the stationary solution and in scenarios involving wet-dry fronts, especially for the relatively strong interactions between the two coupled models.

We present a method for the numerical approximation of an optimal control problem constrained by a full-space transmission problem. The problem can be controlled by the interior source term, the jump over the interface, the jump of the flux over the interface, or any combination thereof. We complete the first-order optimality condition by the Johnson–Nédélec coupling, and a Lagrange multiplier which corresponds to the Bielak–MacCamy coupling. Our final formulation fulfills the LBB conditions on the continuous as well as discrete level, without restrictions. We develop a reliable a posteriori error estimator, propose an adaptive algorithm, and show its rate optimality. Numerical experiments confirm our theoretical findings.

This article presents a MATLAB software package for solving a three-dimensional (3D) second-order elliptic problem with mixed boundary conditions using the primal hybrid finite element method (FEM). First, we introduce a novel fast 3D uniform finite element mesh refinement technique implemented in MATLAB and establish the Nodes-to-Edge and Faces-to-Tetrahedron connectivity through an efficient and systematic approach. We then describe an efficient MATLAB assembly procedure for the 3D lowest-order primal hybrid finite element matrices. Furthermore, we develop a vectorized Schur complement solver, where the computational improvement over MATLAB’s default direct solver ( mldivide ) arises from reducing the original block system to a substantially smaller Schur complement system. The run-time performance of the software is demonstrated through numerical experiments.

We develop multipoint stress mixed finite element methods for linear elasticity with weakly enforced stress symmetry on distorted quadrilateral grids, which can be reduced to positive definite cell-centered systems. The methods utilize the lowest-order Brezzi–Douglas–Marini finite element spaces for the stress and employ vertex quadrature rules to localize the interaction of degrees of freedom. This approach allows for local stress elimination around each vertex. We introduce two methods. The first method uses a piecewise constant rotation, resulting in a cell-centered system for the displacement and the rotation. The second method employs a continuous piecewise bilinear rotation, enabling further elimination of the rotation and resulting in a cell-centered system for the displacement only. The methods utilize a non-symmetric vertex quadrature rule for the stress bilinear form and both non-symmetric and symmetric vertex quadrature rules for the asymmetry bilinear forms. Stability and error analysis are performed for both methods. First-order convergence is established for all variables in the L 2 {L^{2}} -norm. Numerical results are presented that verify the theoretical results.

Large linear systems are ubiquitous in modern computational science and engineering. Their efficient solution is frequently challenging, especially in conjunction with problems described by parametric PDEs, where many such systems have to be solved. The primary approach for solving large linear systems is the use of Krylov subspace iterative methods with well-designed preconditioners. Recently, graph neural networks (GNNs) have been used in conjunction with solving parametric PDEs, and have been shown to be promising tool for designing preconditioners. Nevertheless, the GNN-based preconditioners reported in the literature have worse effect on system’s spectrum than analogous preconditioners from classical linear algebra. Here we employ well-established preconditioners from linear algebra as starting point for training GNN to obtain preconditioners that yield a more substantial reduction in the condition number of the systems when compared to classical preconditioners. Numerical experiments demonstrate the efficiency of our approach in comparison to classical and known neural network-based methods for parametric PDEs. In addition, a heuristic justification for the loss function employed in this study is provided, and we demonstrate that preconditioners obtained by learning with this loss function reduce the condition number in a way that is more desirable for the conjugate gradient iterative method.

We introduce a framework for repurposing error estimators for source problems to compute an estimator for the gap between eigenspaces and their discretizations. Of interest are eigenspaces of finite clusters of eigenvalues of unbounded nonselfadjoint linear operators with compact resolvent. Eigenspaces and eigenvalues of rational functions of such operators are studied as a first step. Under an assumption of convergence of resolvent approximations in the operator norm and an assumption on global reliability of source problem error estimators, we show that the gap in eigenspace approximations can be bounded by a globally reliable and computable error estimator. Also included are applications of the theoretical framework to first-order system least squares (FOSLS) discretizations and discontinuous Petrov–Galerkin (DPG) discretizations, both yielding new estimators for the error gap. Numerical experiments with a selfadjoint model problem and with a leaky nonselfadjoint waveguide eigenproblem show that adaptive algorithms using the new estimators give refinement patterns that target the cluster as a whole instead of individual eigenfunctions.

Roos and Stynes [Some open questions in the numerical analysis of singularly perturbed differential equations, Comput. Methods Appl. Math. 15 2015, 4, 531–550] posed several open problems, one of which is obtaining uniform convergence of higher-order finite element methods on Bakhvalov-type meshes. In this article, we propose a finite element analysis of any order on an eXp-Bakhvalov mesh to solve a two-dimensional singularly perturbed boundary value problem, whose solution exhibits exponential layers. A careful selection of the interpolation operator, considering the characteristics of the layers, allows the finite element method to achieve optimal-order convergence with respect to the singular perturbation parameter. Numerical results are presented to support the theoretical findings.

The paper is concerned with a posteriori estimates for approximations of boundary value problems generated by the spectral fractional Laplace operator. The derivation is based upon the Stinga–Torrea extension, which generalizes the Caffarelli–Silvestre extension and transfers the corresponding nonlocal problem in a bounded domain to a local problem of higher dimensionality. A posteriori estimates are first derived for this local problem. Two–sided error bounds for the original problem follow from them. The estimates are fully computable and contain no conditions and constants depending on a method or mesh used to compute an approximation. They are valid for any energy admissible approximation of the extended problem.