Partial Differential Equation Solver Research Articles

Optimal Dual Schemes for Adaptive Grid Based Hexmeshing

Hexahedral meshes are a ubiquitous domain for the numerical resolution of partial differential equations. Computing a pure hexahedral mesh from an adaptively refined grid is a prominent approach to automatic hexmeshing, and requires the ability to restore the all hex property around the hanging nodes that arise at the interface between cells having different size. The most advanced tools to accomplish this task are based on mesh dualization. These approaches use topological schemes to regularize the valence of inner vertices and edges, such that dualizing the grid yields a pure hexahedral mesh. In this article, we study in detail the dual approach, and propose four main contributions to it: (i) We enumerate all the possible transitions that dual methods must be able to handle, showing that prior schemes do not natively cover all of them; (ii) We show that schemes are internally asymmetric, therefore not only their construction is ambiguous, but different implementative choices lead to hexahedral meshes with different singular structure; (iii) We explore the combinatorial space of dual schemes, selecting the minimum set that covers all the possible configurations and also yields the simplest singular structure in the output hexmesh; (iv) We enlarge the class of adaptive grids that can be transformed into pure hexahedral meshes, relaxing one of the tight topological requirements imposed by previous approaches. Our extensive experiments show that our transition schemes consistently outperform prior art in terms of ability to converge to a valid solution, amount and distribution of singular mesh edges, and element count. Last but not least, we publicly release our code and reveal a conspicuous amount of technical details that were overlooked in previous literature, lowering an entry barrier that was hard to overcome for practitioners in the field.

Integrating deep neural networks with full-waveform inversion: Reparameterization, regularization, and uncertainty quantification

Full-waveform inversion (FWI) is an accurate imaging approach for modeling the velocity structure by minimizing the misfit between recorded and predicted seismic waveforms. However, the strong nonlinearity of FWI resulting from fitting oscillatory waveforms can trap the optimization in local minima. We have adopted a neural-network-based full-waveform inversion (NNFWI) method that integrates deep neural networks with FWI by representing the velocity model with a generative neural network. Neural networks can naturally introduce spatial correlations as regularization to the generated velocity model, which suppresses noise in the gradients and mitigates local minima. The velocity model generated by neural networks is input to the same partial differential equation (PDE) solvers used in conventional FWI. The gradients of the neural networks and PDEs are calculated using automatic differentiation, which back propagates gradients through the acoustic PDEs and neural network layers to update the weights of the generative neural network. Experiments on 1D velocity models, the Marmousi model, and the 2004 BP model determine that NNFWI can mitigate local minima, especially for imaging high-contrast features such as salt bodies, and it significantly improves the inversion in the presence of noise. Adding dropout layers to the neural network model also allows analyzing the uncertainty of the inversion results through Monte Carlo dropout. NNFWI opens a new pathway to combine deep learning and FWI for exploiting the characteristics of deep neural networks and the high accuracy of PDE solvers. Because NNFWI does not require extra training data and optimization loops, it provides an attractive and straightforward alternative to conventional FWI.

Numerical Simulations of a Fluidized Granular Flow Entry Into Water: Insights Into Modeling Tsunami Generation by Pyroclastic Density Currents

AbstractThe tsunami generation potential of pyroclastic density currents (PDCs) entering the sea is poorly understood, due to limited data and observations. Thus far, tsunami generation by PDCs has been modeled in a similar manner to tsunami generation associated with landslides or debris flows, using two‐layer depth‐averaged approaches. Using the adaptive partial differential equation solver Basilisk and benchmarking with published laboratory experiments, this work explores some of the important parameters not yet accounted for in numerical models of PDC‐generated tsunamis. We use assumptions derived from experimental literature to approximate the granular, basal flow component of a PDC as a dense Newtonian fluid flowing down an inclined plane. This modeling provides insight into how the boundary condition of the slope and the viscosity of the dense granular‐fluid influence the characteristics of the waves generated. It is shown that the boundary condition of the slope has a first‐order impact on the interaction dynamics between the fluidized granular flow and water, as well as the energy transfer from the flow to the generated wave. The experimental physics is captured well in the numerical model, which confirms the underlying assumption of Newtonian fluid‐like behavior in the context of wave generation. The results from this study suggest the importance of considering vertical density and velocity stratification in wave generation models. Furthermore, we demonstrate that granular‐fluids denser than water are capable of shearing the water surface and generating significant amplitude waves, despite vigorous overturning.

Learning the solution operator of parametric partial differential equations with physics-informed DeepONets.

Partial differential equations (PDEs) play a central role in the mathematical analysis and modeling of complex dynamic processes across all corners of science and engineering. Their solution often requires laborious analytical or computational tools, associated with a cost that is markedly amplified when different scenarios need to be investigated, for example, corresponding to different initial or boundary conditions, different inputs, etc. In this work, we introduce physics-informed DeepONets, a deep learning framework for learning the solution operator of arbitrary PDEs, even in the absence of any paired input-output training data. We illustrate the effectiveness of the proposed framework in rapidly predicting the solution of various types of parametric PDEs up to three orders of magnitude faster compared to conventional PDE solvers, setting a previously unexplored paradigm for modeling and simulation of nonlinear and nonequilibrium processes in science and engineering.

Coral: a parallel spectral solver for fluid dynamics and partial differential equations

Coral: a parallel spectral solver for fluid dynamics and partial differential equations

High-performance symbolic-numerics via multiple dispatch

As mathematical computing becomes more democratized in high-level languages, high-performance symbolic-numeric systems are necessary for domain scientists and engineers to get the best performance out of their machine without deep knowledge of code optimization. Naturally, users need different term types either to have different algebraic properties for them, or to use efficient data structures. To this end, we developed Symbolics.jl, an extendable symbolic system which uses dynamic multiple dispatch to change behavior depending on the domain needs. In this work we detail an underlying abstract term interface which allows for speed without sacrificing generality. We show that by formalizing a generic API on actions independent of implementation, we can retroactively add optimized data structures to our system without changing the pre-existing term rewriters. We showcase how this can be used to optimize term construction and give a 113x acceleration on general symbolic transformations. Further, we show that such a generic API allows for complementary term-rewriting implementations. Exploiting this feature, we demonstrate the ability to swap between classical term-rewriting simplifiers and e-graph-based term-rewriting simplifiers. We illustrate how this symbolic system improves numerical computing tasks by showcasing an e-graph ruleset which minimizes the number of CPU cycles during expression evaluation, and demonstrate how it simplifies a real-world reaction-network simulation to halve the runtime. Additionally, we show a reaction-diffusion partial differential equation solver which is able to be automatically converted into symbolic expressions via multiple dispatch tracing, which is subsequently accelerated and parallelized to give a 157x simulation speedup. Together, this presents Symbolics.jl as a next-generation symbolic-numeric computing environment geared towards modeling and simulation.

A Fast and Fully Parallel Analog CMOS Solver for Nonlinear PDEs

A general-purpose analog computing method is proposed to compute the continuous-time solutions of nonlinear partial differential equations (PDEs). The discrete-time difference operator in the standard finite difference time domain (FDTD) method is replaced by continuous-time delay operators that can be realized using analog all-pass filters. The resulting spatially discrete time-continuous (SDTC) update equations are realized using analog circuits which compute continuous-time solutions of the PDE with prescribed initial and boundary conditions. The proposed concept is demonstrated in simulation via an integrated circuit (IC) design of a nonlinear acoustic wave equation solver in 180 nm CMOS technology. Analog arithmetic operations (multiply, scale, and add) are realized in parallel using fully differential op-amps and analog multipliers. The proposed IC computes the PDE solution in parallel at 33 discrete spatial points and has a simulated bandwidth and power consumption of approximately 2 MHz and 3 W, respectively. The performance of the IC is simulated using foundry-supplied device models and quantified using i) the mean squared difference between the circuit simulation results and FDTD simulations, and ii) the noise to signal energy ratio. Acceptable accuracy is obtained, with error metric values varying between -7 and -30 dB for various configurations of the problem. Comparison of the custom analog IC simulations with MATLAB- and C-based FDTD code running on a modern workstation shows an expected average speedup of 205× and 140×, respectively.

Kelvin transformations for simulations on infinite domains

Solving partial differential equations (PDEs) on infinite domains has been a challenging task in physical simulations and geometry processing. We introduce a general technique to transform a PDE problem on an unbounded domain to a PDE problem on a bounded domain. Our method uses the Kelvin Transform, which essentially inverts the distance from the origin. However, naive application of this coordinate mapping can still result in a singularity at the origin in the transformed domain. We show that by factoring the desired solution into the product of an analytically known (asymptotic) component and another function to solve for, the problem can be made continuous and compact, with solutions significantly more efficient and well-conditioned than traditional finite element and Monte Carlo numerical PDE methods on stretched coordinates. Specifically, we show that every Poisson or Laplace equation on an infinite domain is transformed to another Poisson (Laplace) equation on a compact region. In other words, any existing Poisson solver on a bounded domain is readily an infinite domain Poisson solver after being wrapped by our transformation. We demonstrate the integration of our method with finite difference and Monte Carlo PDE solvers, with applications in the fluid pressure solve and simulating electromagnetism, including visualizations of the solar magnetic field. Our transformation technique also applies to the Helmholtz equation whose solutions oscillate out to infinity. After the transformation, the Helmholtz equation becomes a tractable equation on a bounded domain without infinite oscillation. To our knowledge, this is the first time that the Helmholtz equation on an infinite domain is solved on a bounded grid without requiring an artificial absorbing boundary condition.

Kelvin transformations for simulations on infinite domains

Solving partial differential equations (PDEs) on infinite domains has been a challenging task in physical simulations and geometry processing. We introduce a general technique to transform a PDE problem on an unbounded domain to a PDE problem on a bounded domain. Our method uses the Kelvin Transform, which essentially inverts the distance from the origin. However, naive application of this coordinate mapping can still result in a singularity at the origin in the transformed domain. We show that by factoring the desired solution into the product of an analytically known (asymptotic) component and another function to solve for, the problem can be made continuous and compact, with solutions significantly more efficient and well-conditioned than traditional finite element and Monte Carlo numerical PDE methods on stretched coordinates. Specifically, we show that every Poisson or Laplace equation on an infinite domain is transformed to another Poisson (Laplace) equation on a compact region. In other words, any existing Poisson solver on a bounded domain is readily an infinite domain Poisson solver after being wrapped by our transformation. We demonstrate the integration of our method with finite difference and Monte Carlo PDE solvers, with applications in the fluid pressure solve and simulating electromagnetism, including visualizations of the solar magnetic field. Our transformation technique also applies to the Helmholtz equation whose solutions oscillate out to infinity. After the transformation, the Helmholtz equation becomes a tractable equation on a bounded domain without infinite oscillation. To our knowledge, this is the first time that the Helmholtz equation on an infinite domain is solved on a bounded grid without requiring an artificial absorbing boundary condition.

Numerical Analysis on the Storage of Nuclear Waste in Gas-Saturated Deep Coal Seam

Nuclear power has contributed humanity a lot since its successful usage in electricity power generation. According to the global statistics, nuclear power accounts for 16% of the total electricity generation in 2020. However, the rapid development of nuclear power also brings up some problems, in which the storage of nuclear waste is the thorny one. This work carries out a series of modeling and simulation analysis on the geological storage of nuclear waste in a gas-saturated deep coal seam. As the first step, a coupled heat-solid-gas model with three constitutional fields of heat transfer, coal deformation, and gas seepage that based on three governing conservation equations is proposed. The approved mechanical model covers series of interactive influences among temperature change, dual permeability of coal, thermal stress, and gas sorption. As the second step, a finite element numerical model and numerical simulation are developed to analyze the storage of nuclear waste in a gas-saturated deep coal seam based on the partial differential equations (PDE) solver of COMSOL Multiphysics with MATLAB. The numerical simulation is implemented and solved then to draw the following conclusions as the nuclear waste chamber heats up the surrounding coal seam firstly in the initial storage stage of 400 years and then be heated by the far-field reservoir. The initial velocity of gas flow decreases gradually with the increment of distance from the storage chamber. Coal gas flows outward from the central storage chamber to the outer area in the first 100 years when the gas pressure in the region nearby the central storage chamber is higher than that in the far region and flows back then while the temperature in the outer region is higher. The modeling and simulation studies are expected to provide a deep understanding on the geological storage of nuclear waste.

Solving of partial differential equations with distributed order in time using fractional-order Bernoulli-Legendre functions

In this paper, an efficient numerical method is used to provide the approximate solution of distributed-order fractional partial differential equations (DFPDEs). The proposed method is based on the fractional integral operator of fractional-order Bernoulli-Legendre functions and the collocation scheme. In our technique, by approximating functions that appear in the DFPDEs by fractional-order Bernoulli functions in space and fractional-order Legendre functions in time using Gauss-Legendre numerical integration, the under study problem is converted to a system of algebraic equations. This system is solved by using Newton's iterative scheme, and the numerical solution of DFPDEs is obtained. Finally, some numerical experiments are included to show the accuracy, efficiency, and applicability of the proposed method.

A robust solver for elliptic PDEs in 3D complex geometries

We develop a boundary integral equation solver for elliptic partial differential equations on complex 3D geometries. Our method is efficient, high-order accurate and robustly handles complex geometries. A key component is our singular and near-singular layer potential evaluation scheme, hedgehog: a simple extrapolation of the solution along a line to the boundary. We present a series of geometry-processing algorithms required for hedgehog to run efficiently with accuracy guarantees on arbitrary geometries and an adaptive upsampling scheme based on a iteration-free heuristic for quadrature error. We validate the accuracy and performance with a series of numerical tests and compare our approach to a competing local evaluation method.

Solving inverse-PDE problems with physics-aware neural networks

We propose a novel composite framework to find unknown fields in the context of inverse problems for partial differential equations (PDEs). We blend the high expressibility of deep neural networks as universal function estimators with the accuracy and reliability of existing numerical algorithms for partial differential equations as custom layers in semantic autoencoders. Our design brings together techniques of computational mathematics, machine learning and pattern recognition under one umbrella to incorporate domain-specific knowledge and physical constraints to discover the underlying hidden fields. The network is explicitly aware of the governing physics through a hard-coded PDE solver layer in contrast to most existing methods that incorporate the governing equations in the loss function or rely on trainable convolutional layers to discover proper discretizations from data. This subsequently focuses the computational load to only the discovery of the hidden fields and therefore is more data efficient. We call this architecture Blended inverse-PDE networks (hereby dubbed BiPDE networks) and demonstrate its applicability for recovering the variable diffusion coefficient in Poisson problems in one and two spatial dimensions, as well as the diffusion coefficient in the time-dependent and nonlinear Burgers' equation in one dimension. We also show that the learned hidden parameters are robust to added noise on input data.

Efficient training of physics‐informed neural networks via importance sampling

AbstractPhysics‐informed neural networks (PINNs) are a class of deep neural networks that are trained, using automatic differentiation, to compute the response of systems governed by partial differential equations (PDEs). The training of PINNs is simulation free, and does not require any training data set to be obtained from numerical PDE solvers. Instead, it only requires the physical problem description, including the governing laws of physics, domain geometry, initial/boundary conditions, and the material properties. This training usually involves solving a nonconvex optimization problem using variants of the stochastic gradient descent method, with the gradient of the loss function approximated on a batch of collocation points, selected randomly in each iteration according to a uniform distribution. Despite the success of PINNs in accurately solving a wide variety of PDEs, the method still requires improvements in terms of computational efficiency. To this end, in this paper, we study the performance of an importance sampling approach for efficient training of PINNs. Using numerical examples together with theoretical evidences, we show that in each training iteration, sampling the collocation points according to a distribution proportional to the loss function will improve the convergence behavior of the PINNs training. Additionally, we show that providing a piecewise constant approximation to the loss function for faster importance sampling can further improve the training efficiency. This importance sampling approach is straightforward and easy to implement in the existing PINN codes, and also does not introduce any new hyperparameter to calibrate. The numerical examples include elasticity, diffusion, and plane stress problems, through which we numerically verify the accuracy and efficiency of the importance sampling approach compared to the predominant uniform sampling approach.

DiscretizationNet: A machine-learning based solver for Navier–Stokes equations using finite volume discretization

Over the last few decades, existing Partial Differential Equation (PDE) solvers have demonstrated a tremendous success in solving complex, non-linear PDEs. Although accurate, these PDE solvers are computationally costly. With the advances in Machine Learning (ML) technologies, there has been a significant increase in the research of using ML to solve PDEs. The goal of this work is to develop an ML-based PDE solver, that couples’ important characteristics of existing PDE solvers with ML technologies. The two solver characteristics that have been adopted in this work are: (1) the use of discretization-based schemes to approximate spatio-temporal partial derivatives and (2) the use of iterative algorithms to solve linearized PDEs in their discrete form. In the presence of highly non-linear, coupled PDE solutions, these strategies can be very important in achieving good accuracy, better stability and faster convergence. Our ML-solver, DiscretizationNet, employs a generative CNN-based encoder–decoder model with PDE variables as both input and output features. During training, the discretization schemes are implemented inside the computational graph to enable faster GPU computation of PDE residuals, which are used to update network weights that result into converged solutions. A novel iterative capability is implemented during the network training to improve the stability and convergence of the ML-solver. The ML-Solver is demonstrated to solve the steady, incompressible Navier–Stokes equations in 3-D for several cases such as, lid-driven cavity, flow past a cylinder and conjugate heat transfer.

Deep Nitsche Method: Deep Ritz Method with Essential Boundary Conditions

We propose a new method to deal with the essential boundary conditions encountered in the deep learning-based numerical solvers for partial differential equations. The trial functions representing by deep neural networks are non-interpolatory, which makes the enforcement of the essential boundary conditions a nontrivial matter. Our method resorts to Nitsche's variational formulation to deal with this difficulty, which is consistent, and does not require significant extra computational costs. We prove the error estimate in the energy norm and illustrate the method on several representative problems posed in at most 100 dimension.

A Low-Rank Schwarz Method for Radiative Transfer Equation With Heterogeneous Scattering Coefficient

Random sampling has been used to find low-rank structure and to build fast direct solvers for multiscale partial differential equations of various types. In this work, we design an accelerated Schw...

Mathematical Modeling of Unsteady Flow with Uniform/Non-Uniform Thermal and Magnetic Fields in a Half-Moon Shaped Domain

The mathematical modeling of two-dimensional unsteady free convective flow and thermal transport inside a half-moon shaped domain charged in the presence of uniform/non-uniform temperature and magnetic effects with Brownian motion of the nanoparticles has been conducted. Thirty-two types of nanofluids in a combination of various nanoparticles and base fluids having different sizes, shapes, and solid concentrations of nanoparticles are chosen to examine the better performance of heat transfer. The circular boundary is cooled while the diameter boundary is heated with uniform/non-uniform temperature. An external uniform/non-uniform/periodic magnetic field is imposed along diameter. The powerful partial differential equations solver, finite element method of Galerkin type, has been engaged in numerical simulation. The numerical solution's heat transfer mechanism reaches a steady-state from the unsteady situation within a very short dimensionless time of about 0.65. The thermal transport rate enhances for increasing buoyancy force whereas decreases with higher magnetic intensity. The uniform thermal condition along the diameter of half-moon gives a higher thermal transport rate compared to non-uniform heating conditions. The non-uniform magnetic field provides greater values of the mean Nusselt number than the uniform field. In addition, the outcomes also predict that a better rate of temperature transport for kerosene-based nanofluid than water-based, ethylene glycol-based, and engine oil-based nanofluid. The heat transfer rate is observed at about 67.86 and 23.78% using Co-Kerosene and Co-water nanofluids, respectively, with an additional 1% nanoparticles volume fraction. The blade shape nanoparticles provide a better heat transfer rate than spherical, cylindrical, brick, and platelet shapes. Small size nanoparticles confirm a higher value of average Nusselt number than big size. Mean Nusselt number increases 22.1 and 5.4% using 1% concentrated Cu-water and Cu-engine oil nanofluid, respectively, than base fluid.

Effect of probe geometry during measurement of >100 A Ga2O3 vertical rectifiers

The high breakdown voltage and low on-state resistance of Schottky rectifiers fabricated on β-Ga2O3 leads to low switching losses, making them attractive for power inverters. One of the main goals is to achieve high forward currents, requiring the fabrication of large area (>1 cm2) devices in order to keep the current density below the threshold for thermally driven failure. A problem encountered during the measurement of these larger area devices is the dependence of current spreading on the probe size, resistance, number, and geometry, which leads to lower currents than expected. We demonstrate how a multiprobe array (6 × 8 mm2) provides a means of mitigating this effect and measure a single sweep forward current up to 135 A on a 1.15 cm2 rectifier fabricated on a vertical Ga2O3 structure. Technology computer-aided design simulations using the floods code, a self-consistent partial differential equation solver, provide a systematic insight into the role of probe placement, size (40–4120 μm), number (1–5), and the sheet resistance of the metal contact on the resultant current-voltage characteristics of the rectifiers.

A nonlinear solver based on an adaptive neural network, introduction and application to porous media flow

Sensitivity to derivatives and a need for proper initial guesses are the main disadvantages of classic nonlinear solvers like Newton's method. To overcome the obstacles, a numerical solver for second-order nonlinear Partial Differential Equations (PDEs) based on an Adaptive Neural Network (AdNN) is introduced. While using Newton's method needs to form the Jacobian matrix and its inversion, which both of them are time-consuming and dependent on the nature of the problem, the proposed solver tries to find the roots of the algebraic equations by changing the weights of AdNN by the help adaptive laws. The proposed approach has been applied to solve the governing PDEs of the gas flow in shale resources and the immiscible two-phase flow of water and oil in hydrocarbon reservoirs as two highly nonlinear phenomena. The generated profiles of pressures and saturations show a satisfying match with the outputs of Newton's method. However, using the presented algorithm not only removes the former necessities but also helps the community to solve the relevant PDEs governing the critical elements of the energy market in the future with a higher level of confidence.