Engineering Papers⌕ Search

SEARCH · Engineering Papers

Results for “finite difference methods”

Search indexed NASA NTRS and DOE OSTI research on propulsion, heat transfer, battery materials and energy systems. Follow report and document links to the original sources.

Quote a phrase for an exact phrase match. Source license links do not imply unrestricted reuse.

At least 19 records
All sourcesNASA NTRSDOE OSTI
DOE OSTI2026

A variational mimetic finite difference method for elliptic interface problems on non-matching polytopal meshes with geometric interface inconsistencies

A new variational mimetic finite difference method for elliptic interface problems with perfect and imperfect thermal contacts on non-matching polytopal meshes with geometric interface inconsistencies is developed and analyzed theoretically and numerically. The method is defined on multiple non-matching submeshes with gaps and overlaps along their interfaces. The discrete equations are derived from a minimization problem for the augmented Dirichlet functional. For a perfect thermal contact, the functional uses a modified mimetic gradient with extended stencil which couples unknowns from both sides of an interface, as well as penalty terms to enforce weak continuity of temperature across the interface. The method leads to a symmetric positive definite matrix for any scaling of the penalty terms. For an imperfect thermal contact, the Dirichlet functional is supplemented with a quadratic jump term along the interface related to the interface thermal resistance. We prove that the method conserves the total heat flux across each interface. In conclusion, the obtained results are verified with numerical experiments showing convergence in the discrete L 2 and L ∞ norms.

97 MATHEMATICS AND COMPUTING↗
DOE OSTI2025

An implementation of a high-order generalized finite difference method for solving the time-harmonic cold plasma wave equation in toroidal geometry

A high-order physics-informed meshless finite difference numerical technique is introduced for solving the time-harmonic cold plasma wave equation in toroidal geometries, presenting a novel application of the generalized finite difference (GFD) method to plasma wave simulations. The algorithm employs an irregular distribution of computational points, with local point density informed by the shortest wavelength derived from the cold plasma dispersion relation. Numerical stability and robustness are addressed using regularization techniques. The algorithm, implemented for two spatial dimensions, solves for the wave electric field and is demonstrated to achieve convergence rates of $\mathcal{O}$($\mathcal{h}$ $\mathcal{P}$ )⁠. Verification tests reproduce plane wave solutions, and example simulations of ion cyclotron resonance heating and electron cyclotron resonance heating demonstrate its capability, approaching realistic tokamak plasma scenarios. This work contributes to laying a foundation for the GFD method to be used in more sophisticated, optimized, and physically realistic full-wave simulations in time-harmonic plasma wave research.

70 PLASMA PHYSICS AND FUSION TECHNOLOGY↗
DOE OSTI2025

Directional finite difference method for directly solving 3D gyrokinetic field equations with enhanced accuracy

The gyrokinetic (GK) field equation is a three-dimensional (3D) elliptic equation, but it is often simplified to a set of two-dimensional (2D) equations by assuming that the field does not vary along a specific direction. However, this simplification can introduce inevitable 0th-order numerical errors, as nonlinear mode coupling in toroidal geometry can produce undesirable harmonic modes that violate the assumption. In this work, we propose a novel directional finite difference method (FDM) with a local coordinate transformation to better resolve the target field of interest. The directional FDM can accurately solve 3D GK field equations without simplifications, which can overcome the limitations of conventional methods. The accuracy and efficiency of different FDMs are analyzed in great detail for a variety of geometries, from simple 2D Cartesian coordinates to realistic 3D curvilinear coordinates. The 0th-order numerical errors of simplified 2D GK equations were found to be more problematic for low-harmonic modes and low aspect ratio geometries such as spherical tokamaks. On the other hand, the directional 3D FDM can accurately resolve a much wider range of harmonic modes aligned to the direction of interest, including the low-harmonic modes. In conclusion, we demonstrate that the directional 3D FDM is a highly effective algorithm for solving the 3D GK field equations, achieving accuracy improvements of 10 to 100 times or more, particularly for low-harmonic modes in spherical tokamaks.

70 PLASMA PHYSICS AND FUSION TECHNOLOGY↗
DOE OSTI2024

2-D seismic wave propagation using the distributional finite-difference method: further developments and potential for global seismology

SUMMARY We present a time-domain distributional finite-difference scheme based on the Lebedev staggered grid for the numerical simulation of wave propagation in acoustic and elastic media. The central aspect of the proposed method is the representation of the stresses and displacements with different sets of B-splines functions organized according to the staggered grid. The distributional finite-difference approach allows domain-decomposition, heterogeneity of the medium, curvilinear mesh, anisotropy, non-conformal interfaces, discontinuous grid and fluid–solid interfaces. Numerical examples show that the proposed scheme is suitable to model wave propagation through the Earth, where sharp interfaces separate large, relatively homogeneous layers. A few domains or elements are sufficient to represent the Earth’s internal structure without relying on advanced meshing techniques. We compare seismograms obtained with the proposed scheme and the spectral element method, and we show that our approach offers superior accuracy, reduced memory usage, and comparable efficiency.

Geochemistry & Geophysics↗
DOE OSTI2023

Generalized moving least squares vs. radial basis function finite difference methods for approximating surface derivatives

Approximating differential operators defined on two-dimensional surfaces is an important problem that arises in many areas of science and engineering. Over the past ten years, localized meshfree methods based on generalized moving least squares (GMLS) and radial basis function finite differences (RBF-FD) have been shown to be effective for this task as they can give high orders of accuracy at low computational cost, and they can be applied to surfaces defined only by point clouds. However, there have yet to be any studies that perform a direct comparison of these methods for approximating surface differential operators (SDOs). The first purpose of this work is to fill that gap. For this comparison, we focus on an RBF-FD method based on polyharmonic spline kernels and polynomials (PHS+Poly) since they are most closely related to the GMLS method. Additionally, we use a relatively new technique for approximating SDOs with RBF-FD called the tangent plane method since it is simpler than previous techniques and natural to use with PHS+Poly RBF-FD. Further, the second purpose of this work is to relate the tangent plane formulation of SDOs to the local coordinate formulation used in GMLS and to show that they are equivalent when the tangent space to the surface is known exactly. The final purpose is to use ideas from the GMLS SDO formulation to derive a new RBF-FD method for approximating the tangent space for a point cloud surface when it is unknown. For the numerical comparisons of the methods, we examine their convergence rates for approximating the surface gradient, divergence, and Laplacian as the point clouds are refined for various parameter choices. We also compare their efficiency in terms of accuracy per computational cost, both when including and excluding setup costs.

97 MATHEMATICS AND COMPUTING↗
DOE OSTI2024

A high-order finite difference method for moving immersed domain boundaries and material interfaces

Here, we present a high-order sharp treatment of immersed moving domain boundaries and material interfaces, and apply it to the advection-diffusion equation in two and three dimensions. The spatial discretization combines dimension-split finite difference schemes with an immersed boundary treatment based on a weighted least-squares reconstruction of the solution, providing stable discretizations with up to sixth order accuracy for diffusion terms and third order accuracy for advection terms. The temporal discretization relies on a novel strategy for maintaining high-order temporal accuracy in problems with moving boundaries that minimizes implementation complexity and allows arbitrary explicit or diagonally-implicit Runge-Kutta schemes. The approach is broadly compatible with popular PDE-specialized Runge-Kutta time integrators, including low-storage, strong stability preserving, and diagonally implicit schemes. Through numerical experiments we demonstrate that the full discretization maintains high-order spatial and temporal accuracy in the presence of complex 3D geometries and for a range of boundary conditions, including Dirichlet, Neumann, and flux conditions with large jumps in coefficients.

97 MATHEMATICS AND COMPUTING↗
DOE OSTI2025

Parallel diffusion operator for magnetized plasmas with improved spectral fidelity

Diffusive transport processes in magnetized plasmas are highly anisotropic, with fast parallel transport along the magnetic field lines sometimes faster than perpendicular transport by orders of magnitude. This constitutes a major challenge for describing non-grid-aligned magnetic structures in Eulerian (grid-based) simulations. Here, the present paper describes and validates a new method for parallel diffusion in magnetized plasmas based on the anti-symmetry representation [Halpern and Waltz, Phys. Plasmas 25, 060703 (2018)]. In the anti-symmetry formalism, diffusion manifests as a flow operator involving the logarithmic derivative of the transported quantity. Qualitative plane wave analysis shows that the new operator naturally yields better discrete spectral resolution compared to its conventional counterpart. Numerical simulations comparing the new method against existing finite difference methods are carried out, showing significant improvement. In particular, we find that combining anti-symmetry with finite differences in diagonally staggered grids essentially eliminates the so-called “artificial numerical diffusion” that affects conventional finite difference and finite volume methods.

Anisotropic diffusion↗
DOE OSTI2026

Enhancing photoionization rate calculations in low-temperature plasmas using spectral methods

Photoionization plays a central role in the development of streamer discharges and other non-equilibrium plasma phenomena. It creates seed electrons, which are essential for positive streamer propagation, allowing the ionization front to move forward. Because of this, accurate modeling of photoionization is very important for predicting streamer behavior and plasma evolution. The photoionization process in air (N 2 – O 2 mixture) is often described by the Zheleznyak model (1982). This model is usually solved through Helmholtz-type equations that approximate the Zheleznyak photoionization model (Zheleznyak et al. 1982) as Partial Differential Equations (PDEs). Conventional numerical methods, such as the Finite Difference Method (FDM) or Finite Volume Method (FVM), are widely used to solve these equations. Although they are prevalent, the computational cost of these methods due to their need for matrix operations and iterative solver is demanding. To address this challenge, this work develops a spectral solver based on the Fast Fourier Transform (FFT) combined with Discrete Cosine Transform (DCT) and Discrete Sine Transform (DST) to calculate the photoionization rate efficiently in an axisymmetric cylindrical domain. This method naturally satisfies the boundary conditions used in the model and converts the PDE into algebraic ones in spectral space. Thus, avoids the need for iterative matrix solvers. When compared with FDM results, it is demonstrated that the new solver not only maintains accuracy, but also reduces the computational cost, showing a performance increase of approximately 100 compared to FDM over a wide range of problem sizes. The method is parallelized using Message Passing Interface (MPI) and has been integrated into a fluid plasma model for streamer simulation. Here, this FFT-based approach provides a fast and reliable alternative for calculating photoionization in fluid models, helping large-scale plasma simulations run faster and efficiently, and allows higher-resolution simulation without extra computational cost.

Axisymmetric system↗
DOE OSTI2025

An Early Investigation of the HHL Quantum Linear Solver for Scientific Applications

In this paper, we explore using the Harrow–Hassidim–Lloyd (HHL) algorithm to address scientific and engineering problems through quantum computing, utilizing the NWQSim simulation package on a high-performance computing platform. Focusing on domains such as power-grid management and climate projection, we demonstrate the correlations of the accuracy of quantum phase estimation, along with various properties of coefficient matrices, on the final solution and quantum resource cost in iterative and non-iterative numerical methods such as the Newton–Raphson method and finite difference method, as well as their impacts on quantum error correction costs using the Microsoft Azure Quantum resource estimator. We summarize the exponential resource cost from quantum phase estimation before and after quantum error correction and illustrate a potential way to reduce the demands on physical qubits. This work lays down a preliminary step for future investigations, urging a closer examination of quantum algorithms’ scalability and efficiency in domain applications.

hybrid software for QC-HPC↗
DOE OSTI2024

A finite difference informed random walker (FDiRW) solver for strongly inhomogeneous diffusion problems

In nature, many complex multi-physics coupling problems exhibit strong diffusivity inhomogeneity. For instance, in the context of radionuclide absorption by porous wasteform materials within a flowing waste stream, the difference of species’ diffusivity in solid and liquid phases spans by 3~8 orders of magnitude. To solve the diffusion equations with strongly inhomogeneous diffusivity, traditional discretization-based methods, such as the Finite Difference Method (FDM), require infinitesimally small time steps (<10 -10 ) as high spatial resolutions are employed in most microstructure evolution processes, leading to prohibitively high computational costs. Here, this work developed an integrated numerical approach (FDiRW: Finite Difference informed Random Walk) to tackle this challenge. The idea is that utilizing the Random Walk concept, the fast diffusion is modeled as a superposition of point source’s solution for a concentration distribution while FDM is used to obtain the point source’s solution at each node. A mesh-coarsening algorithm is developed to generate an exclusive coarse mesh for FDiRW approach to maximize its efficiency. The effectiveness of the coarse mesh-based FDiRW approach is validated by benchmarking Finite Difference solutions. Numerical results demonstrated that FDiRW achieves a remarkable 1000x computational efficiency improvement over FDM while preserving desired accuracy for a medium-sized model of 192 × 192 × 192 grids. Finally, as models scale up, a floating-point operations (PLOPs) analysis of the FDiRW algorithm reveals that its computational complexity grows quadratically in terms of the number of nodes employed in computation.

36 MATERIALS SCIENCE↗
DOE OSTI2023

Acoustic Codes in 2D Spherical Coordinate

Finite-difference methods are widely used to simulate infrasound propagation in the atmosphere. Flexibility of finite-difference scheme allows implementation of highly heterogeneous media for sound propagation as well as complex source models for sound generation. While full 3-D finite-difference methods have been utilized for local infrasound propagation with pronounced topography, 2-D modeling approach has been preferred for regional and global propagation as full 3-D methods generally require enormous computational resources. Infrasound propagation is often simulated with a second-order finite difference scheme. This lowest-order finite-difference scheme is robust and straightforward to implement complex boundary conditions, but the solution includes large error with numerical dispersion and dissipation. This large numerical error may make the second-order finite-difference not optimal for long range infrasound propagation modeling as the numerical dispersion degrades the accuracy of the solution unacceptably. Here, we developed a high-order finite-difference solver for long-range infrasound simulation. The high-order scheme is particularly popular for linear wave modeling in aeroacoustics owing to its low-dispersive and low-dissipative behavior. We develop and evaluate a high-order finite difference scheme in 2-D axisymmetric coordinates. The axisymmetry allows to approximate 3-D spherical sound propagation and amplitude attenuation by a 2-D method. AC2Dr is developed to simulate infrasound propagation in realistic atmosphere, but can be used for linear acoustic waves in general materials with background flow. AC2Dr in an axisymmetric coordinates allows for spherical radiation of acoustic waves from compact sources.

Sjogreen, Bjorn↗
DOE OSTI2023

Finite-difference time-domain methods

The finite-difference time-domain (FDTD) method is a widespread numerical tool for full-wave analysis of electromagnetic fields in complex media and for detailed geometries. Applications of the FDTD method cover a range of time and spatial scales, extending from subatomic to galactic lengths and from classical to quantum physics. Technology areas that benefit from the FDTD method include biomedicine — bioimaging, biophotonics, bioelectronics and biosensors; geophysics — remote sensing, communications, space weather hazards and geolocation; metamaterials — sub-wavelength focusing lenses, electromagnetic cloaks and continuously scanning leaky-wave antennas; optics — diffractive optical elements, photonic bandgap structures, photonic crystal waveguides and ring-resonator devices; plasmonics — plasmonic waveguides and antennas; and quantum applications — quantum devices and quantum radar. This Primer summarizes the main features of the FDTD method, along with key extensions that enable accurate solutions to be obtained for different research questions. Additionally, hardware considerations are discussed, plus examples of how to extract magnitude and phase data, Brillouin diagrams and scattering parameters from the output of an FDTD model. Furthermore, the Primer ends with a discussion of ongoing challenges and opportunities to further enhance the FDTD method for current and future applications.

42 ENGINEERING↗
DOE OSTI2024

A Finite Difference informed Random Walk solver for simulating radiation defect evolution in polycrystalline structures with strongly inhomogeneous diffusivity

Diffusivity of species and defects on grain boundaries is usually several orders of magnitude larger than that inside grains. Such strongly inhomogeneous diffusivity requires prohibitively high computational demands for modeling microstructural evolution. Here, this paper presents a highly-efficient numerical solver, combining the Finite Difference method and Random Walk model, designed for accurately modeling strongly inhomogeneous diffusion within polycrystalline structures. The proposed solver, termed Finite Difference informed Random Walk (FDiRW), integrates a customized Finite Difference (cFD) scheme tailored for fast diffusion along thin grain boundaries represented by a single-layer of nodes. Numerical experiments demonstrate that the FDiRW solver achieves an impressive efficiency gain of 1560x compared to traditional Finite Difference methods while maintaining accuracy, making it feasible for personal computer machines to handle diffusional systems with strongly inhomogeneous diffusivity across static polycrystalline microstructures. The model has been successfully applied to simulate radiation defect evolution, showcasing its scalability to engineering scales in both length and time dimensions.

36 MATERIALS SCIENCE↗
DOE OSTI2022

Deep Learning-Assisted Real-Time Forward Modeling of Electromagnetic Logging in Complex Formations

Higher dimensional (i.e., 2-D and 3-D) modeling is indispensable to correctly evaluate the responses of electromagnetic logging tools in complex formation environments. However, limited by the high computational cost of rigorous modeling such as finite difference method and the finite element method, the real-time applications in the well logging industry primarily rely on the 1-D forward solver, which would result in erroneous formation evaluation for complex scenarios. As a result, aiming at realizing fast modeling for electromagnetic logging tools in complex formations, this paper proposes a general framework assisted by deep neural networks (DNNs). The framework consists of three modules: earth model classification, parameter extraction, and surrogate construction. Separate DNNs are trained and tested for different modules. The accuracy and efficiency of the DNN assisted fast modeling are validated by several experiments. Here, this study finds that the fast modeling assisted by DNNs is able to calculate the tool responses and reconstruct the subsurface formations in real-time.

97 MATHEMATICS AND COMPUTING↗
DOE OSTI2025

Solving high-dimensional partial integral differential equations: The finite expression method

Partial integro-differential equations (PIDEs) have broad applications in the sciences, from electro-magnetism to options pricing. Here, in this paper, we introduce a new finite expression method (FEX) to solve PIDEs. This approach builds upon the original FEX and its inherent advantages with new advances: 1) A novel method of parameter grouping is proposed to reduce the number of coefficients in high-dimensional function approximation; 2) A Taylor series approximation method is implemented to significantly improve the computational efficiency and accuracy of the evaluation of the integral terms of PIDEs. The new FEX based method, denoted FEX-PG to indicate the addition of the parameter grouping (PG) step to the algorithm, provides both high accuracy and interpretable numerical solutions, with the outcome being an explicit equation that facilitates intuitive understanding of the underlying solution structures. These features are often absent in traditional methods, such as finite element methods (FEM) and finite difference methods, as well as in deep learning-based approaches. To benchmark our method against recent advances, we apply the new FEX-PG to solve benchmark PIDEs in the literature. In high-dimensional settings, FEX-PG exhibits strong and robust performance, achieving relative errors on the order of single precision machine epsilon, significantly outperforming existing approaches based on neural networks.

Combinatorial optimization↗
DOE OSTI2025

Constrained or unconstrained? Neural-network-based equation discovery from data

Throughout many fields, practitioners often rely on differential equations to model systems. Yet, for many applications, the theoretical derivation of such equations and/or the accurate resolution of their solutions may be intractable. Instead, recently developed methods, including those based on parameter estimation, operator subset selection, and neural networks, allow for the data-driven discovery of both ordinary and partial differential equations (PDEs), on a spectrum of interpretability. The success of these strategies is often contingent upon the correct identification of representative equations from noisy observations of state variables and, as importantly and intertwined with that, the mathematical strategies utilized to enforce those equations. Specifically, the latter has been commonly addressed via unconstrained optimization strategies. Representing the PDE as a neural network, we propose to discover the PDE (or the associated operator) by solving a constrained optimization problem and using an intermediate state representation similar to a physics-informed neural network (PINN). The objective function of this constrained optimization problem promotes matching the data, while the constraints require that the discovered PDE is satisfied at a number of spatial collocation points. We present a penalty method and a widely used trust-region barrier method to solve this constrained optimization problem, and we compare these methods on numerical examples. Our results on several example problems demonstrate that the latter constrained method outperforms the penalty method, particularly for higher noise levels or fewer collocation points. This work motivates further exploration into using sophisticated constrained optimization methods in scientific machine learning, as opposed to their commonly used, penalty-method or unconstrained counterparts. For both of these methods, we solve these discovered neural network PDEs with classical methods, such as finite difference methods, as opposed to PINNs-type methods relying on automatic differentiation. Here, we briefly highlight how simultaneously fitting the data while discovering the PDE improves the robustness to noise and other small, yet crucial, implementation details.

Data-driven discovery↗
DOE OSTI2025

Explicit block encodings of boundary value problems for many-body elliptic operators

Simulation of physical systems is one of the most promising use cases of future digital quantum computers. In this work we systematically analyze the quantum circuit complexities of block encoding the discretized elliptic operators that arise extensively in numerical simulations for partial differential equations, including high-dimensional instances for many-body simulations. When restricted to rectangular domains with separable boundary conditions, we provide explicit circuits to block encode the many-body Laplacian with separable periodic, Dirichlet, Neumann, and Robin boundary conditions, using standard discretization techniques from low-order finite difference methods. To obtain high-precision, we introduce a scheme based on periodic extensions to solve Dirichlet and Neumann boundary value problems using a high-order finite difference method, with only a constant increase in total circuit depth and subnormalization factor. We then present a scheme to implement block encodings of differential operators acting on more arbitrary domains, inspired by Cartesian immersed boundary methods. We then block encode the many-body convective operator, which describes interacting particles experiencing a force generated by a pair-wise potential given as an inverse power law of the interparticle distance. This work provides concrete recipes that are readily translated into quantum circuits, with depth logarithmic in the total Hilbert space dimension, that block encode operators arising broadly in applications involving the quantum simulation of quantum and classical many-body mechanics.

Kharazi, Tyler [University of California, Berkeley↗