In quantum mechanics, the Rosen-Zener model represents a two-level quantum system. Its generalization to multiple degenerate sets of states leads to larger non-autonomous linear system of ordinary differential equations (ODEs). We propose a new method for computing the solution operator of this system of ODEs. This new method is based on a recently introduced expression of the solution in terms of an infinite matrix equation, which can be efficiently approximated by combining truncation, fixed point iterations, and low-rank approximation. This expression is possible thanks to the so-called $\star$-product approach for linear ODEs. In the numerical experiments, the new method's computing time scales linearly with the model's size. We provide a first partial explanation of this linear behavior.
The Gerber-Shiu function is a classical research topic in actuarial science.However, exact solutions are only available in the literature for very specific cases where the claim amounts follow distributions such as the exponential distribution. This presents a longstanding challenge, particularly from a computational perspective. For the classical risk process in continuous time, the Gerber-Shiu discounted penalty function satisfies a class of Volterra integral equations. In this paper, we use the collocation method to compute the Gerber-Shiu function for risk model with interest. Our methodology demonstrates that the function can be expressed as a linear algebraic system, which is straightforward to implement. One major advantage of our approach is that it does not require any specific distributional assumptions on the claim amounts, except for mild differentiability and continuity conditions that can be easily verified. We also examine the convergence orders of the collocation method. Finally, we present several numerical examples to illustrate the desirable performance of our proposed method.
Phase field models are gradient flows with their energy naturally dissipating in time. In order to preserve this property, many numerical schemes have been well-studied. In this paper we consider a well-known method, namely the exponential integrator method (EI). In the literature a few works studied several EI schemes for various phase field models and proved the energy dissipation by either requiring a strong Lipschitz condition on the nonlinear source term or certain $L^\infty$ bounds on the numerical solutions (maximum principle). However for phase field models such as the (non-local) Cahn-Hilliard equation, the maximum principle no longer exists. As a result, solving such models via EI schemes remains open for a long time. In this paper we aim to give a systematic approach on applying EI-type schemes to such models by solving the Cahn-Hilliard equation with a first order EI scheme and showing the energy dissipation. In fact second order EI schemes can be handled similarly and we leave the discussion in a subsequent paper. To our best knowledge, this is the first work to handle phase field models without assuming any strong Lipschitz condition or $L^\infty$ boundedness. Furthermore, we will analyze the $L^2$ error and present some numerical simulations to demonstrate the dynamics.
Magnetization dynamics in ferromagnetic materials is modeled by the Landau-Lifshitz (LL) equation, a nonlinear system of partial differential equations. Among the numerical approaches, semi-implicit schemes are widely used in the micromagnetics simulation, due to a nice compromise between accuracy and efficiency. At each time step, only a linear system needs to be solved and a projection is then applied to preserve the length of magnetization. However, this linear system contains variable coefficients and a non-symmetric structure, and thus an efficient linear solver is highly desired. If the damping parameter becomes large, it has been realized that efficient solvers are only available to a linear system with constant, symmetric, and positive definite (SPD) structure. In this work, based on the implicit-explicit Runge-Kutta (IMEX-RK) time discretization, we introduce an artificial damping term, which is treated implicitly. The remaining terms are treated explicitly. This strategy leads to a semi-implicit scheme with the following properties: (1) only a few linear system with constant and SPD structure needs to be solved at each time step; (2) it works for the LL equation with arbitrary damping parameter; (3) high-order accuracy can be obtained with high-order IMEX-RK time discretization. Numerically, second-order and third-order IMEX-RK methods are designed in both the 1-D and 3-D domains. A comparison with the backward differentiation formula scheme is undertaken, in terms of accuracy and efficiency. The robustness of both numerical methods is tested on the first benchmark problem from National Institute of Standards and Technology. The linearized stability estimate and optimal rate convergence analysis are provided for an alternate IMEX-RK2 numerical scheme as well.
This work delves into the exponential time differencing (ETD) schemes for the matrix-valued Allen-Cahn equation. In fact, the maximum bound principle (MBP) for the first- and second-order ETD schemes is presented in a prior publication [SIAM Review, 63(2), 2021], assuming a symmetric initial matrix field. Noteworthy is our novel contribution, demonstrating that the first- and second-order ETD schemes for the matrix-valued Allen-Cahn equation -- both being linear schemes -- unconditionally preserve the MBP, even in instances of nonsymmetric initial conditions. Additionally, we prove that these two ETD schemes preserve the energy dissipation law unconditionally for the matrix-valued Allen-Cahn equation. Some numerical examples are presented to verify our theoretical results and to simulate the evolution of corresponding matrix fields.
In PDE-constrained optimization, one aims to find design parameters that minimize some objective, subject to the satisfaction of a partial differential equation. A major challenges is computing gradients of the objective to the design parameters, as applying the chain rule requires computing the Jacobian of the design parameters to the PDE's state. The adjoint method avoids this Jacobian by computing partial derivatives of a Lagrangian. Evaluating these derivatives requires the solution of a second PDE with the adjoint differential operator to the constraint, resulting in a backwards-in-time simulation. Particle-based Monte Carlo solvers are often used to compute the solution to high-dimensional PDEs. However, such solvers have the drawback of introducing noise to the computed results, thus requiring stochastic optimization methods. To guarantee convergence in this setting, both the constraint and adjoint Monte Carlo simulations should simulate the same particle trajectories. For large simulations, storing full paths from the constraint equation for re-use in the adjoint equation becomes infeasible due to memory limitations. In this paper, we provide a reversible extension to the family of permuted congruential pseudorandom number generators (PCG). We then use such a generator to recompute these time-reversed paths for the heat equation, avoiding these memory issues.
In this paper, we formulate and analyse a symmetric low-regularity integrator for solving the nonlinear Klein-Gordon equation in the $d$-dimensional space with $d=1,2,3$. The integrator is constructed based on the two-step trigonometric method and the proposed integrator has a simple form. Error estimates are rigorously presented to show that the integrator can achieve second-order time accuracy in the energy space under the regularity requirement in $H^{1+\frac{d}{4}}\times H^{\frac{d}{4}}$. Moreover, the time symmetry of the scheme ensures the good long-time energy conservation which is rigorously proved by the technique of modulated Fourier expansions. A numerical test is presented and the numerical results demonstrate the superiorities of the new integrator over some existing methods.
This paper proposes a new Helmholtz decomposition based windowed Green function (HD-WGF) method for solving the time-harmonic elastic scattering problems on a half-space with Dirichlet boundary conditions in both 2D and 3D. The Helmholtz decomposition is applied to separate the pressure and shear waves, which satisfy the Helmholtz and Helmholtz/Maxwell equations, respectively, and the corresponding boundary integral equations of type $(\mathbb{I}+\mathbb{T})\bs\phi=\bs f$, that couple these two waves on the unbounded surface, are derived based on the free-space fundamental solution of Helmholtz equation. This approach avoids the treatment of the complex elastic displacement tensor and traction operator that involved in the classical integral equation method for elastic problems. Then a smooth ``slow-rise'' windowing function is introduced to truncate the boundary integral equations and a ``correction'' strategy is proposed to ensure the uniformly fast convergence for all incident angles of plane incidence. Numerical experiments for both two and three dimensional problems are presented to demonstrate the accuracy and efficiency of the proposed method.
A posteriori reduced-order models, e.g. proper orthogonal decomposition, are essential to affordably tackle realistic parametric problems. They rely on a trustful training set, that is a family of full-order solutions (snapshots) representative of all possible outcomes of the parametric problem. Having such a rich collection of snapshots is not, in many cases, computationally viable. A strategy for data augmentation, designed for parametric laminar incompressible flows, is proposed to enrich poorly populated training sets. The goal is to include in the new, artificial snapshots emerging features, not present in the original basis, that do enhance the quality of the reduced-order solution. The methodologies devised are based on exploiting basic physical principles, such as mass and momentum conservation, to devise physically-relevant, artificial snapshots at a fraction of the cost of additional full-order solutions. Interestingly, the numerical results show that the ideas exploiting only mass conservation (i.e., incompressibility) are not producing significant added value with respect to the standard linear combinations of snapshots. Conversely, accounting for the linearized momentum balance via the Oseen equation does improve the quality of the resulting approximation and therefore is an effective data augmentation strategy in the framework of viscous incompressible laminar flows.
Robotic systems are complex cyber-physical systems (CPS) commonly equipped with multiple sensors and effectors. Recent simulation methods enable the Digital Twin (DT) concept realisation. However, DT employment in robotic system development, e.g. in-development testing, is unclear. During the system development, its parts evolve from simulated mockups to physical parts which run software deployed on the actual hardware. Therefore, a design tool and a flexible development procedure ensuring the integrity of the simulated and physical parts are required. We aim to maximise the integration between a CPS's simulated and physical parts in various setups. The better integration, the better simulation-based testing coverage of the physical part (hardware and software). We propose a Domain Specification Language (DSL) based on Systems Modeling Language (SysML) that we refer to as SPSysML (Simulation-Physical System Modeling Language). SPSysML defines the taxonomy of a Simulation-Physical System (SPSys), being a CPS consisting of at least a physical or simulated part. In particular, the simulated ones can be DTs. We propose a SPSys Development Procedure (SPSysDP) that enables the maximisation of the simulation-physical integrity of SPSys by evaluating the proposed factors. SPSysDP is used to develop a complex robotic system for the INCARE project. In subsequent iterations of SPSysDP, the simulation-physical integrity of the system is maximised. As a result, the system model consists of fewer components, and a greater fraction of the system components are shared between various system setups. We implement and test the system with popular frameworks, Robot Operating System (ROS) and Gazebo simulator. SPSysML with SPSysDP enables the design of SPSys (including DT and CPS), multi-setup system development featuring maximised integrity between simulation and physical parts in its setups.
Partitioned neural network functions are used to approximate the solution of partial differential equations. The problem domain is partitioned into non-overlapping subdomains and the partitioned neural network functions are defined on the given non-overlapping subdomains. Each neural network function then approximates the solution in each subdomain. To obtain the convergent neural network solution, certain continuity conditions on the partitioned neural network functions across the subdomain interface need to be included in the loss function, that is used to train the parameters in the neural network functions. In our work, by introducing suitable interface values, the loss function is reformulated into a sum of localized loss functions and each localized loss function is used to train the corresponding local neural network parameters. In addition, to accelerate the neural network solution convergence, the localized loss function is enriched with an augmented Lagrangian term, where the interface condition and the boundary condition are enforced as constraints on the local solutions by using Lagrange multipliers. The local neural network parameters and Lagrange multipliers are then found by optimizing the localized loss function. To take the advantage of the localized loss function for the parallel computation, an iterative algorithm is also proposed. For the proposed algorithms, their training performance and convergence are numerically studied for various test examples.