High-dimensional Partial Differential Equations Research Articles

Deep finite volume method for partial differential equations

In this paper, we introduce the Deep Finite Volume Method (DFVM), an innovative deep learning framework tailored for solving high-order (order ≥2) partial differential equations (PDEs). Our approach centers on a novel loss function crafted from local conservation laws derived from the original PDE, distinguishing DFVM from traditional deep learning methods. By formulating DFVM in the weak form of the PDE rather than the strong form, we enhance accuracy, particularly beneficial for PDEs with less smooth solutions compared to strong-form-based methods like Physics-Informed Neural Networks (PINNs). A key technique of DFVM lies in its transformation of all second-order or higher derivatives of neural networks into first-order derivatives which can be computed directly using Automatic Differentiation (AD). This adaptation significantly reduces computational overhead, particularly advantageous for solving high-dimensional PDEs. Numerical experiments demonstrate that DFVM achieves equal or superior solution accuracy compared to existing deep learning methods such as PINN, Deep Ritz Method (DRM), and Weak Adversarial Networks (WAN), while drastically reducing computational costs. Notably, for PDEs with nonsmooth solutions, DFVM yields approximate solutions with relative errors up to two orders of magnitude lower than those obtained by PINN. The implementation of DFVM is available on GitHub at https://github.com/Sysuzqs/DFVM.

Approximate averaged derivative moments and stochastic response of high-DOF quasi-non-integrable Hamiltonian systems under Poisson white noises

Stochastic response prediction is one of the principal topics in diverse branches of engineering, such as aerospace engineering, mechanical engineering, and vehicle engineering. However, the challenge of addressing high-dimensional partial differential equations or evaluating high-dimensional domain integrals in determining averaged derivative moments has limited the research. This study proposes a method to predict the stochastic response of high-degree-of-freedom (DOF) quasi-non-integrable Hamiltonian systems under Poisson white noises. The original finite DOF system is reduced to a one-dimensional one governing the Hamiltonian by applying the stochastic averaging method, and the associated reduced generalized Fokker–Planck–Kolmogorov (GFPK) equation is derived. Then, the two-step generalized elliptical coordinate transformation is extended to evaluate averaged derivative moments that are the key to obtaining the stationary probability density function (SPDF) for probabilistic information prediction. Moreover, the perturbation method is introduced to address the high-order GFPK equation. Finally, an 8-DOF strongly nonlinear system under purely additive Poisson white noises is given as an example to highlight the accuracy of the proposed method. Results from the proposed method show good agreement with those from Monte Carlo simulations, which indicates that the proposed method can be applied to assess the response of high-DOF quasi-non-integrable Hamiltonian systems under Poisson white noises. In addition, if the potential energy of an arbitrary finite n-DOF system is the same form as that in Appendix B, pertinent high-dimensional domain integrals are converted into explicit expressions when n is even.

Inf-sup neural networks for high-dimensional elliptic PDE problems

Solving high dimensional partial differential equations (PDEs) has historically posed a considerable challenge when utilizing conventional numerical methods, such as those involving domain meshes. Recent advancements in the field have seen the emergence of neural PDE solvers, leveraging deep networks to effectively tackle high dimensional PDE problems. This study introduces Inf-SupNet, a model-based unsupervised learning approach designed to acquire solutions for a specific category of elliptic PDEs. The fundamental concept behind Inf-SupNet involves incorporating the inf-sup formulation of the underlying PDE into the loss function. The analysis reveals that the global solution error can be bounded by the sum of three distinct errors: the numerical integration error, the duality gap of the loss function (training error), and the neural network approximation error for functions within Sobolev spaces. To validate the efficacy of the proposed method, numerical experiments conducted in high dimensions demonstrate its stability and accuracy across various boundary conditions, as well as for both semi-linear and nonlinear PDEs.

Space-Time Deep Neural Network Approximations for High-Dimensional Partial Differential Equations

Space-Time Deep Neural Network Approximations for High-Dimensional Partial Differential Equations

Tackling the curse of dimensionality with physics-informed neural networks

The curse-of-dimensionality taxes computational resources heavily with exponentially increasing computational cost as the dimension increases. This poses great challenges in solving high-dimensional partial differential equations (PDEs), as Richard E. Bellman first pointed out over 60 years ago. While there has been some recent success in solving numerical PDEs in high dimensions, such computations are prohibitively expensive, and true scaling of general nonlinear PDEs to high dimensions has never been achieved. We develop a new method of scaling up physics-informed neural networks (PINNs) to solve arbitrary high-dimensional PDEs. The new method, called Stochastic Dimension Gradient Descent (SDGD), decomposes a gradient of PDEs’ and PINNs’ residual into pieces corresponding to different dimensions and randomly samples a subset of these dimensional pieces in each iteration of training PINNs. We prove theoretically the convergence and other desired properties of the proposed method. We demonstrate in various diverse tests that the proposed method can solve many notoriously hard high-dimensional PDEs, including the Hamilton–Jacobi-Bellman (HJB) and the Schrödinger equations in tens of thousands of dimensions very fast on a single GPU using the PINNs mesh-free approach. Notably, we solve nonlinear PDEs with nontrivial, anisotropic, and inseparable solutions in less than one hour for 1000 dimensions and in 12 h for 100,000 dimensions on a single GPU using SDGD with PINNs. Since SDGD is a general training methodology of PINNs, it can be applied to any current and future variants of PINNs to scale them up for arbitrary high-dimensional PDEs.

A numerical B-spline Galerkin method with proper generalized decomposition for reduced order modeling of partial differential equations

A numerical study is presented on reduced order modeling of low and high-dimensional partial differential equations using a new B-spline Galerkin proper generalized decomposition (PGD) method. Galerkin-PGD schemes based on B-splines are numerically investigated for various benchmark problems, including the Poisson equation, linear/nonlinear advection–diffusion equations, and the Navier–Stokes equations. The treatment of inhomogeneous boundary conditions, non-separable terms and nonlinearities is numerically investigated. Numerical errors for high-fidelity B-spline solutions and reduced order solutions are provided to demonstrate the effectiveness and accuracy of the proposed method. Computational costs are reported to illustrate the efficiency of the developed approaches.

Low-dose CT reconstruction based on high-dimensional partial differential equation projection recovery

We propose a low-dose CT reconstruction method using partial differential equation (PDE) denoising under high-dimensional constraints. The projection data were mapped into a high-dimensional space to construct a high-dimensional representation of the data, which were updated by moving the points in the high-dimensional space. The data were denoised using partial differential equations and the CT image was reconstructed using the FBP algorithm. Compared with those by FBP, PWLS-QM and TGV-WLS methods, the relative root mean square error of the Shepp-Logan image reconstructed by the proposed method were reduced by 68.87%, 50.15% and 27.36%, the structural similarity values were increased by 23.50%, 8.83% and 1.62%, and the feature similarity values were increased by 17.30%, 2.71% and 2.82%, respectively. For clinical image reconstruction, the proposed method, as compared with FBP, PWLS-QM and TGV-WLS methods, resulted in reduction of the relative root mean square error by 42.09%, 31.04% and 21.93%, increased the structural similarity values by 18.33%, 13.45% and 4.63%, and increased the feature similarity values by 3.13%, 1.46% and 1.10%, respectively. The new method can effectively reduce the streak artifacts and noises while maintaining the spatial resolution in reconstructed low-dose CT images.

Physics-Informed Neural Networks for Solving High-Index Differential-Algebraic Equation Systems Based on Radau Methods

As is well known, differential algebraic equations (DAEs), which are able to describe dynamic changes and underlying constraints, have been widely applied in engineering fields such as fluid dynamics, multi-body dynamics, mechanical systems, and control theory. In practical physical modeling within these domains, the systems often generate high-index DAEs. Classical implicit numerical methods typically result in varying order reduction of numerical accuracy when solving high-index systems. Recently, the physics-informed neural networks (PINNs) have gained attention for solving DAE systems. However, it faces challenges like the inability to directly solve high-index systems, lower predictive accuracy, and weaker generalization capabilities. In this paper, we propose a PINN computational framework, combined Radau IIA numerical method with an improved fully connected neural network structure, to directly solve high-index DAEs. Furthermore, we employ a domain decomposition strategy to enhance solution accuracy. We conduct numerical experiments with two classical high-index systems as illustrative examples, investigating how different orders and time-step sizes of the Radau IIA method affect the accuracy of neural network solutions. For different time-step sizes, the experimental results indicate that utilizing a 5th-order Radau IIA method in the PINN achieves a high level of system accuracy and stability. Specifically, the absolute errors for all differential variables remain as low as 10−6, and the absolute errors for algebraic variables are maintained at 10−5. Therefore, our method exhibits excellent computational accuracy and strong generalization capabilities, providing a feasible approach for the high-precision solution of larger-scale DAEs with higher indices or challenging high-dimensional partial differential algebraic equation systems.

An Interpretable Approach to the Solutions of High-Dimensional Partial Differential Equations

In recent years, machine learning algorithms, especially deep learning, have shown promising prospects in solving Partial Differential Equations (PDEs). However, as the dimension increases, the relationship and interaction between variables become more complex, and existing methods are difficult to provide fast and interpretable solutions for high-dimensional PDEs. To address this issue, we propose a genetic programming symbolic regression algorithm based on transfer learning and automatic differentiation to solve PDEs. This method uses genetic programming to search for a mathematically understandable expression and combines automatic differentiation to determine whether the search result satisfies the PDE and boundary conditions to be solved. To overcome the problem of slow solution speed caused by large search space, we propose a transfer learning mechanism that transfers the structure of one-dimensional PDE analytical solution to the form of high-dimensional PDE solution. We tested three representative types of PDEs, and the results showed that our proposed method can obtain reliable and human-understandable real solutions or algebraic equivalent solutions of PDEs, and the convergence speed is better than the compared methods. Code of this project is at https://github.com/grassdeerdeer/HD-TLGP.

Tensor networks for solving the time-independent Boltzmann neutron transport equation

Tensor network techniques, known for their low-rank approximation ability that breaks the curse of dimensionality, are emerging as a foundation of new mathematical methods for ultra-fast numerical solutions of high-dimensional Partial Differential Equations (PDEs). Here, we present a mixed Tensor Train (TT)/Quantized Tensor Train (QTT) approach for the numerical solution of time-independent Boltzmann Neutron Transport equations (BNTEs) in Cartesian geometry. Discretizing a realistic three-dimensional (3D) BNTE by (i) diamond differencing in space, (ii) multigroup-in-energy, and (iii) discrete ordinates collocation in angle leads to large generalized eigenvalue problems that generally require a matrix-free approach and large computer clusters. Starting from this discretization, we construct a TT representation of the PDE fields and discrete operators, followed by a QTT representation of the TT cores. We then solve the tensorized generalized eigenvalue problem using a fixed-point scheme with tensor network optimization techniques. We validate our approach by applying the method to two examples of 3D neutron transport problems, currently solved by the Los Alamos National Laboratory PARallel TIme-dependent SN (PARTISN) solver.1 We demonstrate that our TT/QTT method, executed on a standard desktop computer, leads to large compression. This allows for the storage of terrabyte-sized neutron angular flux eigenvectors in megabytes. Additionally, we create megabyte-sized full access TT representations of yottabyte-sized transport matrix operators. By leveraging the TT operators and solution methods, we obtain a 7500 times speedup when compared to the PARTISN solution time with an error of less than 10−5.

DF-ParPINN: parallel PINN based on velocity potential field division and single time slice focus

Modified Benney-Luke equation (mBL equation) is a three-dimensional temporal-spatial equation with complex structures, that is a high-dimensional partial differential equation (PDE), it is also a new equation of the physical ocean field, and its solution is important for studying the internal wave-wave interaction of inclined seafloor. For conventional PDE solvers such as the pseudo-spectral method, it is difficult to solve mBL equation with both accuracy and speed. Physics-informed neural network (PINN) incorporates physical prior knowledge in deep neural networks, which can solve PDE with relative accuracy and speed. However, PINN is only suitable for solving low-dimensional PDE with simple structures, and not suitable for solving high-dimensional PDE with complex structures. This is mainly because high-dimensional PDEs usually have complex structures and high-order derivatives and are likely to be high-dimensional non-convex functions, and the high-dimensional non-convex optimization problem is an NP-hard problem, resulting in the PINN easily falling into inaccurate local optimal solutions when solving high-dimensional PDEs. Therefore, we improve the PINN for the characteristics of mBL equation and propose “DF-ParPINN: parallel PINN based on velocity potential field division and single time slice focus” to solve mBL equation with large amounts of data. DF-ParPINN consists of three modules: temporal-spatial division module of overall velocity potential field, data rational selection module of multiple time slices, and parallel computation module of high-velocity fields and low-velocity fields. The experimental results show that the solution time of DF-ParPINN is no more than 0.5s, and its accuracy is much higher than that of PINN, PIRNN, cPINN, and DeepONet. Moreover, the relative error of DF-ParPINN after deep training 1000000 epochs can be reduced to less than 0.1. The validity of DF-ParPINN proves that the improved PINN also can solve high dimensional PDE with complex structures and large amounts of data quickly and accurately, which is of great significance to the deep learning of the physical ocean field.

Numerical investigation of high-dimensional option pricing PDEs by utilizing a hybrid radial basis function - finite difference procedure

The target of this research is to resolve high-dimensional partial differential equations (PDEs) for multi-asset options, modeled as parabolic time-dependent PDEs. We present a hybrid radial basis function - finite difference (RBF-FD) solver, which combines the advantages of Gaussian and multiquadric functions. Additionally, we employ the Krylov subspace method on the resulting system of ordinary differential equations, reducing the computation load for finding the numerical solution. Computational tests, up to seven dimensions, and comparisons support the superiority of our hybrid solver.

Hutchinson Trace Estimation for high-dimensional and high-order Physics-Informed Neural Networks

Physics-Informed Neural Networks (PINNs) have proven effective in solving partial differential equations (PDEs), especially when some data are available by seamlessly blending data and physics. However, extending PINNs to high-dimensional and even high-order PDEs encounters significant challenges due to the computational cost associated with automatic differentiation in the residual loss function calculation. Herein, we address the limitations of PINNs in handling high-dimensional and high-order PDEs by introducing the Hutchinson Trace Estimation (HTE) method. Starting with the second-order high-dimensional PDEs, which are ubiquitous in scientific computing, HTE is applied to transform the calculation of the entire Hessian matrix into a Hessian vector product (HVP). This approach not only alleviates the computational bottleneck via Taylor-mode automatic differentiation but also significantly reduces memory consumption from the Hessian matrix to an HVP’s scalar output. We further showcase HTE’s convergence to the original PINN loss and its unbiased behavior under specific conditions. Comparisons with the Stochastic Dimension Gradient Descent (SDGD) highlight the distinct advantages of HTE, particularly in scenarios with significant variability and variance among dimensions. We further extend the application of HTE to higher-order and higher-dimensional PDEs, specifically addressing the biharmonic equation. By employing tensor-vector products (TVP), HTE efficiently computes the colossal tensor associated with the fourth-order high-dimensional biharmonic equation, saving memory and enabling rapid computation. The effectiveness of HTE is illustrated through experimental setups, demonstrating comparable convergence rates with SDGD under memory and speed constraints. Additionally, HTE proves valuable in accelerating the Gradient-Enhanced PINN (gPINN) version as well as the Biharmonic equation. Overall, HTE opens up a new capability in scientific machine learning for tackling high-order and high-dimensional PDEs.

A DEIM Tucker tensor cross algorithm and its application to dynamical low-rank approximation

We introduce a Tucker tensor cross approximation method that constructs a low-rank representation of a d-dimensional tensor by sparsely sampling its fibers. These fibers are selected using the discrete empirical interpolation method (DEIM). Our proposed algorithm is referred to as DEIM fiber sampling (DEIM-FS). For a rank-r approximation of an O(Nd) tensor, DEIM-FS requires access to only dNrd−1 tensor entries, a requirement that scales linearly with the tensor size along each mode. We demonstrate that DEIM-FS achieves an approximation accuracy close to the Tucker-tensor approximation obtained via higher-order singular value decomposition at a significantly reduced cost. We also present DEIM-FS (iterative) that does not require access to singular vectors of the target tensor unfolding and can be viewed as a black-box Tucker tensor algorithm. We employ DEIM-FS to reduce the computational cost associated with solving nonlinear tensor differential equations (TDEs) using dynamical low-rank approximation (DLRA). The computational cost of solving DLRA equations can become prohibitive when the exact rank of the right-hand side tensor is large. This issue arises in many TDEs, especially in cases involving non-polynomial nonlinearities, where the right-hand side tensor has full rank. This necessitates the storage and computation of tensors of size O(Nd). We show that DEIM-FS results in significant computational savings for DLRA by constructing a low-rank Tucker approximation of the right-hand side tensor on the fly. Another advantage of using DEIM-FS is to significantly simplify the implementation of DLRA equations, irrespective of the type of TDEs. We demonstrate the efficiency of the algorithm through several examples including solving high-dimensional partial differential equations.

The One Step Malliavin scheme: new discretization of BSDEs implemented with deep learning regressions

Abstract A novel discretization is presented for decoupled forward–backward stochastic differential equations (FBSDE) with differentiable coefficients, simultaneously solving the BSDE and its Malliavin sensitivity problem. The control process is estimated by the corresponding linear BSDE driving the trajectories of the Malliavin derivatives of the solution pair, which implies the need to provide accurate $\varGamma $ estimates. The approximation is based on a merged formulation given by the Feynman–Kac formulae and the Malliavin chain rule. The continuous time dynamics is discretized with a theta-scheme. In order to allow for an efficient numerical solution of the arising semidiscrete conditional expectations in possibly high dimensions, it is fundamental that the chosen approach admits to differentiable estimates. Two fully-implementable schemes are considered: the BCOS method as a reference in the one-dimensional framework and neural network Monte Carlo regressions in case of high-dimensional problems, similarly to the recently emerging class of Deep BSDE methods (Han et al. (2018 Solving high-dimensional partial differential equations using deep learning. Proc. Natl. Acad. Sci., 115, 8505–8510); Huré et al. (2020 Deep backward schemes for high-dimensional nonlinear PDEs. Math. Comp., 89, 1547–1579)). An error analysis is carried out to show $\mathbb{L}^2$ convergence of order $1/2$, under standard Lipschitz assumptions and additive noise in the forward diffusion. Numerical experiments are provided for a range of different semilinear equations up to $50$ dimensions, demonstrating that the proposed scheme yields a significant improvement in the control estimations.

Adaptive deep neural networks for solving corner singular problems

Deep neural networks (DNNs) for numerical solutions to partial differential equations (PDEs) have exhibited their remarkable merits of meshless methods, dimensionless features, and nonlinear approximation powers. The majority of conventional DNNs are mainly concentrated on the low- or high-dimensional PDE with smooth solutions. In singular problems the DNN does not behave so efficiently as in the smooth problem because of low regularities of the solution. In this paper, we propose a novel adaptive DNN for the PDE with corner singular solutions. In the corner singular problem gradients of solution vary rapidly around some particular vertices or lines. This poses a difficulty in development of DNN solvers. We design adaptive DNN techniques for the corner singularity, which incorporates (a) an adaptive loss function that assigns radial weights to local losses, (b) an adaptive activation function that varies the choice of activation functions from neuron to neuron, and (c) an adaptive sampling strategy that selectively focuses on sampling points near the singular locations. The adaptive DNN with these adaptive schemes greatly enhances approximation accuracy and computation efficiency. This is verified by a comprehensive series of numerical experiments of 2D and 3D corner singular problems. The comparison analyzes with the existing DNN solvers are also made to demonstrate effectiveness of the new method. The L2 relative errors of the proposed adaptive DNN are reduced by an order of magnitude with respect to several conventional DNN methods.

A deep-genetic algorithm (deep-GA) approach for high-dimensional nonlinear parabolic partial differential equations

We propose a new method, called a deep-genetic algorithm (deep-GA), to accelerate the performance of the so-called deep-BSDE method, which is a deep learning algorithm to solve high dimensional partial differential equations through their corresponding backward stochastic differential equations (BSDEs). Recognizing the sensitivity of the solver to the initial guess selection, we embed a genetic algorithm (GA) into the solver to optimize the selection. We aim to achieve faster convergence for the nonlinear PDEs on a broader interval than deep-BSDE. Our proposed method is applied to two nonlinear parabolic PDEs, i.e., the Black-Scholes (BS) equation with default risk and the Hamilton-Jacobi-Bellman (HJB) equation. We compare the results of our method with those of the deep-BSDE and show that our method provides comparable accuracy with significantly improved computational efficiency.

Fast and scalable computation of shape-morphing nonlinear solutions with application to evolutional neural networks

We develop fast and scalable methods for computing reduced-order nonlinear solutions (RONS). RONS was recently proposed as a framework for reduced-order modeling of time-dependent partial differential equations (PDEs), where the modes depend nonlinearly on a set of time-varying parameters. RONS uses a set of ordinary differential equations (ODEs) for the parameters to optimally evolve the shape of the modes to adapt to the PDE's solution. This method has already proven extremely effective in tackling challenging problems such as advection-dominated flows and high-dimensional PDEs. However, as the number of parameters grow, integrating the RONS equation and even its formation become computationally prohibitive. Here, we develop three separate methods to address these computational bottlenecks: symbolic RONS, collocation RONS and regularized RONS. We demonstrate the efficacy of these methods on two examples: Fokker–Planck equation in high dimensions and the Kuramoto–Sivashinsky equation. In both cases, we observe that the proposed methods lead to several orders of magnitude in speedup and accuracy. Our proposed methods extend the applicability of RONS beyond reduced-order modeling by making it possible to use RONS for accurate numerical solution of linear and nonlinear PDEs. Finally, as a special case of RONS, we discuss its application to problems where the PDE's solution is approximated by a neural network, with the time-dependent parameters being the weights and biases of the network. The RONS equations dictate the optimal evolution of the network's parameters without requiring any training.

An extreme learning machine-based method for computational PDEs in higher dimensions

We present two effective methods for solving high-dimensional partial differential equations (PDE) based on randomized neural networks. Motivated by the universal approximation property of this type of networks, both methods extend the extreme learning machine (ELM) approach from low to high dimensions. With the first method the unknown solution field in d dimensions is represented by a randomized feed-forward neural network, in which the hidden-layer parameters are randomly assigned and fixed while the output-layer parameters are trained. The PDE and the boundary/initial conditions, as well as the continuity conditions (for the local variant of the method), are enforced on a set of random interior/boundary collocation points. The resultant linear or nonlinear algebraic system, through its least squares solution, provides the trained values for the network parameters. With the second method the high-dimensional PDE problem is reformulated through a constrained expression based on an Approximate variant of the Theory of Functional Connections (A-TFC), which avoids the exponential growth in the number of terms of TFC as the dimension increases. The free field function in the A-TFC constrained expression is represented by a randomized neural network and is trained by a procedure analogous to the first method. We present ample numerical simulations for a number of high-dimensional linear/nonlinear stationary/dynamic PDEs to demonstrate their performance. These methods can produce accurate solutions to high-dimensional PDEs, in particular with their errors reaching levels not far from the machine accuracy for relatively lower dimensions. Compared with the physics-informed neural network (PINN) method, the current method is both cost-effective and more accurate for high-dimensional PDEs.

Spectro-spatial evolution of the CMB. Part I. Discretisation of the thermalisation Green's function

Spectral distortions of the cosmic microwave background (CMB) have been recognized as an important future probe of the early Universe. Existing theoretical studies primarily focused on describing the evolution and creation of average distortions, ignoring spatial perturbations in the plasma. One of the main reasons for this choice is that a treatment of the spectro-spatial evolution of the photon field deep into the primordial Universe requires solving a radiative transfer problem for the distortion signals, which in full detail is computationally challenging. Here we provide the first crucial step towards tackling this problem by formulating a new spectral discretisation of the underlying average thermalisation Green's function. Our approach allows us to convert the high-dimensional partial differential equation system (≃ 103–104 equations) into and set of ordinary differential equations of much lower dimension (≃ 10 equations). We demonstrate the precision of the approach and highlight how it may be further improved in the future. We also clarify the link of the observable spectral distortion parameters (e.g., μ and y) to the computational spectral basis that we use in our frequency discretisation. This reveals how several basis-dependent ambiguities can be interpreted in future CMB analysis. Even if not exact, the new Green's function discretisation can be used to formulate a generalised photon Boltzmann-hierarchy, which can then be solved with methods that are familiar from theoretical studies of the CMB temperature and polarisation anisotropies. We will carry this program out in a series of companion papers, thereby opening the path to full spectro-spatial exploration of the CMB with future CMB imagers and spectrometers.