Journal of Scientific Computing最新文献

We propose and study conformal integrators for linearly damped stochastic Poisson systems. We analyse the qualitative and quantitative properties of these numerical integrators: preservation of dynamics of certain Casimir and Hamiltonian functions, almost sure bounds of the numerical solutions, and strong and weak rates of convergence under appropriate conditions. These theoretical results are illustrated with several numerical experiments on, for example, the linearly damped free rigid body with random inertia tensor or the linearly damped stochastic Lotka-Volterra system.

We propose and analyze a space-time Local Discontinuous Galerkin method for the approximation of the solution to parabolic problems. The method allows for very general discrete spaces and prismatic space-time meshes. Existence and uniqueness of a discrete solution are shown by means of an inf-sup condition, whose proof does not rely on polynomial inverse estimates. Moreover, for piecewise polynomial spaces satisfying an additional mild condition, we show a second inf-sup condition that provides additional control over the time derivative of the discrete solution. We derive hp-a priori error bounds based on these inf-sup conditions, which we use to prove convergence rates for standard, tensor-product, and quasi-Trefftz polynomial spaces. Numerical experiments validate our theoretical results.

Neural network-based approaches have recently shown significant promise in solving partial differential equations (PDEs) in science and engineering, especially in scenarios featuring complex domains or incorporation of empirical data. One advantage of the neural network methods for PDEs lies in its automatic differentiation (AD), which necessitates only the sample points themselves, unlike traditional finite difference (FD) approximations that require nearby local points to compute derivatives. In this paper, we quantitatively demonstrate the advantage of AD in training neural networks. The concept of truncated entropy is introduced to characterize the training property. Specifically, through comprehensive experimental and theoretical analyses conducted on random feature models and two-layer neural networks, we discover that the defined truncated entropy serves as a reliable metric for quantifying the residual loss of random feature models and the training speed of neural networks for both AD and FD methods. Our experimental and theoretical analyses demonstrate that, from a training perspective, AD outperforms FD in solving PDEs.

In this paper, we present a novel training approach called the Homotopy Relaxation Training Algorithm (HRTA), aimed at accelerating the training process in contrast to traditional methods. Our algorithm incorporates two key mechanisms: one involves building a homotopy activation function that seamlessly connects the linear activation function with the $ReLU$ activation function; the other technique entails relaxing the homotopy parameter to enhance the training refinement process. We have conducted an in-depth analysis of this novel method within the context of the neural tangent kernel (NTK), revealing significantly improved convergence rates. Our experimental results, especially when considering networks with larger widths, validate the theoretical conclusions. This proposed HRTA exhibits the potential for other activation functions and deep neural networks.

Free-stream preservation is an essential property for numerical solvers on curvilinear grids. Key to this property is that the metric terms of the curvilinear mapping satisfy discrete metric identities, i.e., have zero divergence. Divergence-free metric terms are furthermore essential for entropy stability on curvilinear grids. We present a new way to compute the metric terms for discontinuous Galerkin spectral element methods (DGSEMs) that guarantees they are divergence-free. The proposed mimetic approach uses projections that fit within the de Rham Cohomology.

We develop a sparse hierarchical hp-finite element method (hp-FEM) for the Helmholtz equation with variable coefficients posed on a two-dimensional disk or annulus. The mesh is an inner disk cell (omitted if on an annulus domain) and concentric annuli cells. The discretization preserves the Fourier mode decoupling of rotationally invariant operators, such as the Laplacian, which manifests as block diagonal mass and stiffness matrices. Moreover, the matrices have a sparsity pattern independent of the order of the discretization and admit an optimal complexity factorization. The sparse hp-FEM can handle radial discontinuities in the right-hand side and in rotationally invariant Helmholtz coefficients. Rotationally anisotropic coefficients that are approximated by low-degree polynomials in Cartesian coordinates also result in sparse linear systems. e consider examples such as a high-frequency Helmholtz equation with radial discontinuities and rotationally anisotropic coefficients, singular source terms, țhe time-dependent Schrödinger equation, and an extension to a three-dimensional cylinder domain, with a quasi-optimal solve, via the Alternating Direction Implicit (ADI) algorithm.

We present and analyze a discontinuous Galerkin method for the numerical modeling of a Kelvin-Voigt thermo/poro-viscoelastic problem. We present the derivation of the model and we develop a stability analysis in the continuous setting that holds both for the full inertial and quasi-static problems and that is robust with respect to most of the physical parameters of the problem. For spatial discretization, we propose an arbitrary-order weighted symmetric interior penalty scheme that supports general polytopal grids and is robust with respect to strong heterogeneities in the model coefficients. For the semi-discrete problem, we prove the extension of the stability result demonstrated in the continuous setting and we provide an a-priori error estimate. A wide set of numerical simulations is presented to assess the convergence and robustness properties of the proposed method. Moreover, we test the scheme with literature and physically sound test cases for proof-of-concept applications in the geophysical context.

We present a new approach for nonlinear dimensionality reduction, specifically designed for computationally expensive mathematical models. We leverage autoencoders to discover a one-dimensional neural active manifold (NeurAM) capturing the model output variability, through the aid of a simultaneously learnt surrogate model with inputs on this manifold. Our method only relies on model evaluations and does not require the knowledge of gradients. The proposed dimensionality reduction framework can then be applied to assist outer loop many-query tasks in scientific computing, like sensitivity analysis and multifidelity uncertainty propagation. In particular, we prove, both theoretically under idealized conditions, and numerically in challenging test cases, how NeurAM can be used to obtain multifidelity sampling estimators with reduced variance by sampling the models on the discovered low-dimensional and shared manifold among models. Several numerical examples illustrate the main features of the proposed dimensionality reduction strategy and highlight its advantages with respect to existing approaches in the literature.

The computational effort entailed in the discretization of fluid-poromechanics systems is typically highly demanding. This is particularly true for models of multiphysics flows in the brain, due to the geometrical complexity of the cerebral anatomy-requiring a very fine computational mesh for finite element discretization-and to the high number of variables involved. Indeed, this kind of problems can be modeled by a coupled system encompassing the Stokes equations for the cerebrospinal fluid in the brain ventricles and Multiple-network Poro-Elasticity (MPE) equations describing the brain tissue, the interstitial fluid, and the blood vascular networks at different space scales. The present work aims to rigorously derive a posteriori error estimates for the coupled Stokes-MPE problem, as a first step towards the design of adaptive refinement strategies or reduced order models to decrease the computational demand of the problem. Through numerical experiments, we verify the reliability and optimal efficiency of the proposed a posteriori estimator and identify the role of the different solution variables in its composition.

This paper describes novel algorithms for the identification of (almost-)resonant behavior in scattering problems. Our methods, relying on rational approximation, aim at building surrogate models of what we call "field amplification", defined as the norm of the solution operator of the scattering problem, which we express through boundary-integral equations. To provide our techniques with theoretical foundations, we first derive results linking the field amplification to the spectral properties of the operator that defines the scattering problem. Such results are then used to justify the use of rational approximation in the surrogate-modeling task. Some of our proposed methods apply rational approximation in a "standard" way, building a rational approximant for either the solution operator directly or, in the interest of computational efficiency, for a randomly "sketched" version of it. Our other "hybrid" approaches are more innovative, combining rational-approximation-assisted root-finding with approximation using radial basis functions. Three key features of our methods are that (i) they are agnostic of the strategy used to discretize the scattering problem, (ii) they do not require any computations involving non-real wavenumbers, and (iii) they can adjust to different settings through the use of adaptive sampling strategies. We carry out some numerical experiments involving 2D scatterers to compare our approaches. In our tests, two of our approaches (one standard, one hybrid) emerge as the best performers, with one or the other being preferable, depending on whether emphasis is placed on accuracy or efficiency.