Differential Equation Solver Research Articles

A Sparse Differential Algebraic Equation (DAE) and Stiff Ordinary Differential Equation (ODE) Solver in Maple

This paper implements efficient numerical methods in Maple to solve index-1 nonlinear Differential Algebraic Equations (DAEs) and stiff Ordinary Differential Equations (ODEs) systems. Single-step methods (like Trapezoid (TR), Implicit-mid point (IMP), Euler-backward (EB), Radau IIA (Rad) methods, TRBDF2, TRX2) and backward-difference formula of order 2 are implemented with adaptive time-stepping methods in Maple to solve index-1 nonlinear DAEs. Maple’s robust and efficient ability to search within a list/set is exploited to identify the sparsity pattern and automatically calculate the analytic Jacobian. The algorithm and implementation are robust and efficient for index-1 DAE problems and scale well for finite difference/finite element discretization of two-dimensional models with system size up to 10,000 nonlinear DAEs and solve the same in a few seconds. The computational efficiency of the proposed algorithm (provided as an open-access code) compares favorably with the commercial solver gPROMs, one of the most commonly used sparse DAE solvers in the industry.

Solving difference equations using fourier transform method

Solving difference equations using fourier transform method

Simulation and control of interactions in multi-physics, a Python package for port-Hamiltonian systems

The Python package SCRIMP (Simulation and ContRol of Interactions in Multi-Physics) is presented through a collection of port-Hamiltonian systems (pHs) of increasing complexity, stemming from mechanics and thermodynamics. A focus is made on the syntax of SCRIMP allowing the user to easily describe a distributed pHs and its discretization method using the Partitioned Finite Element Method (PFEM) in space, together with the Differential Algebraic Equation (DAE) solver to use. A Graphical User Interface (GUI) is presented.

Numerical modeling of time dependent Diffusive Shock Acceleration

Motivated by cosmic ray (CR) re-acceleration at a potential Galactic Wind Termination Shock (GWTS), we present a numerical model for time-dependent Diffusive Shock Acceleration (DSA). We use the stochastic differential equation solver (DiffusionSDE) of the cosmic ray propagation framework CRPropa3.2 with two modifications: An importance sampling module is introduced to improve statistics at high energies in order to keep the simulation time short. An adaptive time step is implemented in the DiffusionSDE module. This ensures to efficiently meet constraints on the time and diffusion step, which is crucial to obtain the correct shock spectra. The time evolution of the spectrum at a one-dimensional planar shock is verified against the solution obtained by the grid-based solver VLUGR3 for both energy-independent and energy-dependent diffusion. We show that the injection of pre-accelerated particles can lead to a broken power law spectrum in momentum if the incoming spectrum of CRs is harder than the re-accelerated spectrum. If the injected spectrum is steeper, the shock spectrum dominates at all energies. We finally apply the developed model to the GWTS by considering a spherically symmetric shock, a spiral Galactic magnetic field, and anisotropic diffusion. The time-dependent spectrum at the shock is modeled as a basis for further studies.

Deciphering and integrating invariants for neural operator learning with various physical mechanisms.

Neural operators have been explored as surrogate models for simulating physical systems to overcome the limitations of traditional partial differential equation(PDE) solvers. However, most existing operator learning methods assume that the data originate from a single physical mechanism, limiting their applicability and performance in more realistic scenarios. To this end, we propose the physical invariant attention neural operator (PIANO) to decipher and integrate the physical invariants for operator learning from the PDE series with various physical mechanisms. PIANO employs self-supervised learning to extract physical knowledge and attention mechanisms to integrate them into dynamic convolutional layers. Compared to existing techniques, PIANO can reduce the relative error by 13.6%-82.2% on PDE forecasting tasks across varying coefficients, forces or boundary conditions. Additionally, varied downstream tasks reveal that the PI embeddings deciphered by PIANO align well with the underlying invariants in the PDE systems, verifying the physical significance of PIANO.

Applications of the Z transformation in calculating sum of a series

The Z transformation has many applications in Mathematics, especially in solving difference equations to help handle discrete data models. This article provides one more application of the Z transformation, which is the summation of a series. The author through the research method of theory development from some properties of the Z transformation to obtain the result which is a theorem about the formula for the sum of a series, with the proof attached. From there, apply the theorem together with the Z transformation to calculate the sum of some series in specific problems. Thus, the problem of calculating the sum of the series has one more method of solving, as well as expanding the application of the Z transformation in the field of Mathematics, especially analysis.

Accurate error estimation for model reduction of nonlinear dynamical systems via data-enhanced error closure

Accurate error estimation is crucial in model order reduction, to obtain small reduced-order models as well as to certify their accuracy when deployed in downstream applications such as digital twins. In existing a posteriori error estimation approaches, knowledge about the time integration scheme is mandatory, e.g., the residual-based error estimators proposed for the reduced basis method. This poses a challenge when automatic ordinary differential equation solver libraries are used to perform the time integration. To address this, we present a data-enhanced approach for a posteriori error estimation. Our new formulation enables residual-based error estimators to be independent of any time integration method. To achieve this, we introduce a corrected reduced-order model which takes into account a data-driven closure term for improved accuracy. The closure term, subject to mild assumptions, is related to the local truncation error of the corresponding time integration scheme. We propose efficient computational schemes for approximating the closure term, at the cost of a modest amount of training data. Furthermore, the new error estimator is incorporated within a greedy process to obtain parametric reduced-order models. Numerical results on three different systems show the accuracy of the proposed error estimation approach and its ability to produce ROMs that generalize well.

Hyperbolic modeling of gradient damage and one-dimensional finite volume simulations

This study provides a new formulation of gradient damage model which allows an efficient explicit numerical solution of dynamics problems. The proposed methodology is based on an ”extended Lagrangian approach” developed by one of the authors for the nondissipative and dispersive shallow water equation. By using this strategy, the global minimization problem commonly derived for gradient damage models is recast as a purely local hyperbolic one with source terms that can be easily solved using finite volumes. The numerical solution of the governing system is then based on a fractional-step method consisting of a classical Godunov-type scheme and an implicit Ordinary Differential Equation solver for the local source terms. Numerical results are presented on one-dimensional multi-fragmentation and spalling tests for illustration purpose.

Unsteady magneto-hydrodynamic behavior of TiO2-kerosene nanofluid flow in wavy octagonal cavity

The aim of this study is to mathematically model the behavior of magneto-hydrodynamics (MHD) in a two-dimensional unsteady flow, specifically focusing on natural convection (NC) and heat transfer (HT) within a wavy octagonal cavity. The study investigates how heat is transferred and fluid flows within this cavity under specific conditions. In particular, it examines the impact of various parameters, such as the Hartmann number (Ha), Rayleigh number (Ra), and the volume fraction of nanoparticles (ϕ) on the flow and HT patterns. This cavity includes a rectangular vertical wall (RVW) at the center of its bottom wall, filled with kerosene-TiO2 nanofluid of spherical shape. Within this setup, the RVW is maintained at a high temperature (T = Th), while the wavy wall is kept at a lower temperature (T = Tc, where Tc <Th). All other boundaries of the domain are assumed to be adiabatic. The finite element method (FEM) is employed as the solver for the relevant partial differential equations in numerical simulations. The results show excellent agreement with previously published research papers. The numerical solution transitions from an unsteady state to a steady state in approximately 0.68 dimensionless time units during the HT process. Throughout the study, various parameters are explored, including Ha, ranging from 0 to 100, Ra, ranging from 103 to 106, and ϕ, ranging from 0 to 0.05. The findings reveal distinct patterns in streamlines and isotherms. Specifically, a reduction in the rate of HT is observed as the Lorentz force increases, while the rate of HT is enhanced with increasing buoyancy force. Additionally, an expansion in ϕ led to an increase in the rate of HT. The study's findings could contribute to our understanding of fluid dynamics, heat transfer, and the behavior of nanofluids in complex geometries, potentially leading to improvements in various engineering and industrial applications.

A fast and efficient computation method for reflective diffraction simulations

We present a new computation method for simulating reflection high-energy electron diffraction and the total-reflection high-energy positron diffraction experiments. The two experiments are used commonly for the structural analysis of material surface. The present paper improves the conventional numerical method, the multi-slice method, for faster computation, since the present method avoids the matrix-eigenvalue solver for the computation of matrix exponentials and can adopt higher-order ordinary differential equation solvers. Moreover, we propose a high-performance implementation based on multi-thread parallelization and cache-reusable subroutines. In our tests, this new method performs up to 2,000 times faster than the conventional method without compromising accuracy.

Automated translation and accelerated solving of differential equations on multiple GPU platforms

We demonstrate a high-performance vendor-agnostic method for massively parallel solving of ensembles of ordinary differential equations (ODEs) and stochastic differential equations (SDEs) on GPUs. The method is integrated with a widely used differential equation solver library in a high-level language (Julia’s DifferentialEquations.jl) and enables GPU acceleration without requiring code changes by the user. Our approach achieves state-of-the-art performance compared to hand-optimized CUDA-C++ kernels while performing 20–100× faster than the vectorizing map (vmap) approach implemented in JAX and PyTorch. Performance evaluation on NVIDIA, AMD, Intel, and Apple GPUs demonstrates performance portability and vendor-agnosticism. We show composability with MPI to enable distributed multi-GPU workflows. The implemented solvers are fully featured – supporting event handling, automatic differentiation, and incorporation of datasets via the GPU’s texture memory – allowing scientists to take advantage of GPU acceleration on all major current architectures without changing their model code and without loss of performance. We distribute the software as an open-source library, DiffEqGPU.jl.

Structural reliability analysis of aircraft wing rib fatigue cracking using surrogate dynamic Bayesian network

AbstractDynamic Bayesian networks (DBNs) are crucial for evaluating time‐sensitive structural risks in aircraft, yet their inherent complexity can create substantial computational demands. This research introduces a surrogate model, specifically a two‐layer artificial neural network (ANN), to replace the most computationally intensive node in the DBN, typically associated with an ordinary differential equation solver. The surrogate model was chosen after comparing 24 regression machine learning models and optimizing hyperparameters, leading to a significant reduction in computational time and an improvement over traditional methods such as Monte Carlo simulations. The surrogate model demonstrates practical significance with its conservative risk estimations and successful application to a real‐world structural fatigue issue in an F‐4E aircraft intermediate rib. The unique aspect of this research lies in the strategic application of ANN to mitigate computational challenges, thereby enhancing the performance of DBN in fatigue risk analysis.

Direct Estimation of Parameters in ODE Models Using WENDy: Weak-Form Estimation of Nonlinear Dynamics

We introduce the Weak-form Estimation of Nonlinear Dynamics (WENDy) method for estimating model parameters for non-linear systems of ODEs. Without relying on any numerical differential equation solvers, WENDy computes accurate estimates and is robust to large (biologically relevant) levels of measurement noise. For low dimensional systems with modest amounts of data, WENDy is competitive with conventional forward solver-based nonlinear least squares methods in terms of speed and accuracy. For both higher dimensional systems and stiff systems, WENDy is typically both faster (often by orders of magnitude) and more accurate than forward solver-based approaches. The core mathematical idea involves an efficient conversion of the strong form representation of a model to its weak form, and then solving a regression problem to perform parameter inference. The core statistical idea rests on the Errors-In-Variables framework, which necessitates the use of the iteratively reweighted least squares algorithm. Further improvements are obtained by using orthonormal test functions, created from a set of C^{infty } bump functions of varying support sizes.We demonstrate the high robustness and computational efficiency by applying WENDy to estimate parameters in some common models from population biology, neuroscience, and biochemistry, including logistic growth, Lotka-Volterra, FitzHugh-Nagumo, Hindmarsh-Rose, and a Protein Transduction Benchmark model. Software and code for reproducing the examples is available at https://github.com/MathBioCU/WENDy.

Numerical investigations of activation energy on the peristaltic transport of Carreau nanofluid through a curved asymmetric channel

Nanometric particles with base liquids cause the production of nanofluids, which are distinguished by their outstanding thermally conductive fluid properties and the expansion of electrical and mechanical devices. Based on these considerations, we devised a study to investigate the effect of activation energy on the peristaltic motion of Carreau nanofluid inside a curved asymmetric channel under the influence of a magnetic field. The governing equations for the curved channel of non-Newtonian fluid flow are formulated. The nonlinear partial differential equations system has been reduced to ordinary differential equations by the assumptions of low Reynolds number and long wavelength approximations. The resulting nonlinear coupled differential equations are numerically solved directly using NDSolve (numerical differential equation solver) coding of computational mathematical software Mathematica, and velocity, temperature, concentration, and streamlines are plotted. With graphical demonstrations, the influence of essential parameters on velocity, temperature, concentration, and streamlines is explained in detail. The dimensionless temperature distribution grows as the activation energy parameter grows. In reality, the number of energetic particles (with energies equal to or greater than activation energy) increases, resulting in improved temperature distribution.

On-chip all-optical second-order ordinary differential equation solver based on a single microdisk resonator.

This paper proposes an all-optical second-order ordinary differential equation (SODE) solver based on a single microdisk resonator. We validate the feasibility of our structure for constant and complex coefficient SODE solutions for Gaussian and super-Gaussian pulses. The results demonstrate a good agreement between the solutions obtained with the designed structure and those obtained through mathematical calculations for both constant and complex coefficient SODEs. We also discuss the influence of input optical signal pulse width on solution result deviations. Furthermore, we validate the capability of the designed structure to achieve tunable solutions for complex-coefficient SODEs with a tuning power of less than 10 mW. The device footprint is approximately 20×30 μm2, and it is 3-4 times smaller than the current smallest solving unit. The maximum Q-factor reaches 9.8×104. The proposed device avoids the traditional approach of cascading two resonators for SODE solving. Moreover, achieving mode alignment within the same resonator reduces the process challenges associated with aligning multiple devices in a cascade. Furthermore, it offers wider applicability for solving SODEs, namely, the ability to solve both constant and complex coefficient SODEs with complete derivative terms.

Solving the pulsar equation using physics-informed neural networks

ABSTRACT In this study, Physics-Informed Neural Networks (PINNs) are skilfully applied to explore a diverse range of pulsar magnetospheric models, specifically focusing on axisymmetric cases. The study successfully reproduced various axisymmetric models found in the literature, including those with non-dipolar configurations, while effectively characterizing current sheet features. Energy losses in all studied models were found to exhibit reasonable similarity, differing by no more than a factor of three from the classical dipole case. This research lays the groundwork for a reliable elliptic Partial Differential Equation solver tailored for astrophysical problems. Based on these findings, we foresee that the utilization of PINNs will become the most efficient approach in modelling three-dimensional magnetospheres. This methodology shows significant potential and facilitates an effortless generalization, contributing to the advancement of our understanding of pulsar magnetospheres.

Wideband parametric baseband macromodeling of linear and passive photonic circuits via complex vector fitting

A novel wideband parametric baseband macromodeling technique for passive photonic devices and circuits is presented. It allows to efficiently estimate the baseband scattering representations of a linear, passive photonic system as a function of a set of design variables, such as geometrical layout or substrate features. The proposed technique relies on the interpolation of macromodels computed via a complex vector fitting (CVF) algorithm, by adopting a methodology based on amplitude and frequency scaling that preserves, by construction, the physical properties of the system, such as causality, stability and passivity. For a specified combination of the design parameters, a rational CVF model is derived that can be simulated by a wide range of ordinary differential equation (ODE) solvers or circuit simulators. Additionally, time-domain simulations using the computed model can be performed at arbitrary optical carrier frequencies, thus allowing for the simulation of multi-wavelength systems. Two application examples are presented to demonstrate the flexibility and advantages of the proposed method.

Impact on Casson Nanofluid Flow Across an Inclined, Slanted Surface by Radiation, Energy and Mass Transfer

Casson nanofluid flow across a porous, slanted surface is investigated. In addition to the impact of a heat source, thermal energy, and a chemical reaction also are considered. Applying similarity transformations, the controlling nonlinear PDEs are converted to ODEs. To resolve the ordinary differential equations (ODEs), we used the universal numerical differential equation solver Wolfram Language function NDSolve. The distributions under heat source are affected by chemical processes and other variables are explored. While liquid velocity reduced due to inclination and the Casson factor, mass and energy transit rates rose. The findings are represented graphically and are consistent with past research. The findings would provide valuable insights into the flow patterns, temperature distribution, and concentration profiles within the system. This information can be crucial for designing and optimizing heat transfer systems, energy-efficient processes, and catalytic reactors involving Casson nanofluids and porous media.

Express packaging waste recycling: Stakeholders’ dynamic behavioral changes based on a computational evolutionary game theoretic approach

Identifying stakeholder behaviors is crucial in managing municipal express packaging waste. Prior studies have neglected the impacts of the advertising effect and environmental preference and thus fail to detect stakeholders’ behavioral changes under different policies. To address this limitation, we build a computational evolutionary game model that integrates the advertising effect, environmental preference benefits, and economic and information policies. A numerical simulation is then conducted via a differential equation solver embedded in MATLAB software to identify the behavioral changes of the government, express delivery enterprises, and consumers. The results indicate that subsidies induce enterprises and consumers to reach their stability strategy more quickly, with maximum behavioral change rates of 0.82% and 0.64%, respectively. The marginal impacts of information policy and environmental preference benefits are greatest for consumers, whose behavioral change rates can reach 1.33% and 1.52%, respectively. For enterprises, a 1.71% behavioral change rate can be achieved by increasing the government’s effort level. An increase in enterprises’ effort level causes both enterprises and consumers to reach a stability strategy more quickly and achieves behavioral change rates of 4.21% and 1.51%, respectively. Our model can be used to predict stakeholders’ behavioral changes in the waste disposal industry, and the findings of this study provide novel ideas for waste management practitioners to use in encouraging recycling behaviors.

An Extended Generalized Average Modeling Framework for Power Converters

The Generalized Averaged Modeling (GAM) technique is traditionally employed to capture the dynamic performance of power electronic converters. This paper proposes an improved version of it, named the Extended-GAM (EGAM) technique, which supports the multiplication of two Double Fourier Series (DFS) signals in the time domain. Multiplication of DFS signals in the time domain translates to the 2D-convolution of coefficients of the DFS terms of their equivalent Discrete Fourier Image (DFI) representations. Thus, the proposed EGAM technique, capable of capturing many harmonics present in the output of a power converter, effectively captures the dynamic behavior of power converters excited by two distinct frequencies. The proposed technique is then converted into an algorithm suitable for numerical platforms, which typically use Ordinary Differential Equation (ODE) solvers. The proposed algorithm is validated based on the observations of the effects of harmonic truncation. The efficacy of the proposed technique is assessed through a case study, wherein a single-phase inverter employs LC filters on both the dc-link and the ac-side. Finally, it is shown that the results obtained with the proposed method show an excellent congruence between simulation and hardware experimental models. Additionally, the proposed algorithm is packaged into a MATLAB toolbox and shared for future implementations.