Standard Runge-Kutta Method Research Articles

Finite Elements with Switch Detection for numerical optimal control of nonsmooth dynamical systems with set-valued heaviside step functions

This paper develops high-accuracy methods for the numerical solution of optimal control problems subject to nonsmooth differential equations with set-valued Heaviside step functions. An important subclass of these systems are Filippov systems. By writing the Heaviside step function as the solution map of a linear program and using its optimality conditions, the initial nonsmooth system is rewritten into an equivalent Dynamic Complementarity System (DCS). The Finite Elements with Switch Detection (FESD) method (Nurkanović et al., 2024) was originally developed for Filippov systems transformed via Stewart’s reformulation into DCS (Stewart, 1990). This paper extends it to the above mentioned class of nonsmooth systems. The key ideas are to start with a standard Runge–Kutta method for the DCS and to let the integration step sizes to be degrees of freedom. Then, additional conditions are introduced to allow implicit but accurate switch detection and to remove possible spurious degrees of freedom when no switches occur. The theoretical properties of the FESD method are studied. The motivation for these developments is to obtain a computationally tractable formulation of nonsmooth optimal control problems. Numerical simulations and optimal control examples are used to illustrate the favorable properties of the proposed approach. All methods introduced in this paper are implemented in the open-source software package nosnoc (Nurkanović and Diehl, 2022).

A generalized model for the population dynamics of a two stage species with recruitment and capture using a nonstandard finite difference scheme

The aim of this work is to formulate and analyze a new and generalized discrete-time population dynamics model for a two-stage species with recruitment and capture factors. This model is derived from a well-known continuous-time population dynamics model of a two-stage species with recruitment and capture developed by Ladino and Valverde and the nonstandard finite difference (NSFD) methodology proposed by Mickens. We establish positivity and asymptotic stability of the proposed discrete-time population dynamics model. As an important consequence, the population dynamics of the new discrete-time model is determined fully. Also, a set of numerical examples is conducted to illustrate the theoretical results and to demonstrate advantages of the new model. The theoretical results and numerical examples show that the proposed discrete-time model not only preserves correctly the population dynamics of the continuous one but is also easy to be implemented. However, some discrete-time models based on the standard Runge–Kutta methods fail to preserve the population dynamics of the continuous-time model. As a result, they generate numerical approximations which are not only non-negative but also unstable.

A modified Runge–Kutta method for increasing stability properties

In control theory, we are interested mainly in investigating asymptotic stability of systems. In this paper explicit Runge–Kutta methods are investigated for numerical solutions of nonlinear ODEs with known strict Lyapunov functions. In order to increase stability properties, a modified version of explicit Runge–Kutta methods is presented. The modified Runge–Kutta method is explicit. It is different from the projection methods in literature which are implicit. Comparing with the standard Runge–Kutta method, the modified Runge–Kutta one has better stability properties. Numerical examples are provided to illustrate the effectiveness of the modified Runge–Kutta method.

Explicit two-step peer methods with reused stages

Two-step peer methods for the numerical solution of Initial Value Problems (IVP) combine the advantages of Runge-Kutta (RK) and multistep methods to obtain high stage order and provide in a natural way a dense output. In general, explicit s-stage peer methods require s evaluations of the vector field at each step. Nevertheless, Klinge and coworkers (BIT Numer Math, 2018) have shown that some methods use less function calls se<s, here called effective stages, by re-using sr=s−se previously computed stages (shifted stages) from the previous steps in the current one.In this paper we propose a new approach, different from the one used by Klinge and coworkers, to re-use previously computed stages, that we call peer methods with reused stages, showing that methods with reused stages and se effective stages are equivalent to three-step peer methods with se stages. Then, we analyze all the families of methods with two effective stages, obtaining methods with s=3 and orders 4 and 5 in which the free parameters of the families have been used to minimize the coefficient of the leading error term as well as to maximize the absolute stability interval. We have also studied one family of peer methods with s=4 and three effective stages, obtaining a method with order 6, superconvergent of order 7, and optimized leading error term as well as absolute stability interval. Some numerical experiments show the performance of the obtained methods by comparing them with other previously obtained peer methods as well as other standard Runge-Kutta and multistep methods.

Mathematical model to investigate transmission dynamics of COVID-19 with vaccinated class

&lt;abstract&gt;&lt;p&gt;The susceptible, exposed, infected, quarantined and vaccinated (SEIQV) population is accounted for in a mathematical model of COVID-19. This model covers the therapy for diseased people as well as therapeutic measures like immunization for susceptible people to enable understanding of the dynamics of the disease's propagation. Each of the equilibrium points, i.e., disease-free and endemic, has been proven to be globally asymptotically stable under the assumption that $ \mathscr{R}_0 $ is smaller or larger than unity, respectively. Although vaccination coverage is high, the basic reproduction number depends on the vaccine's effectiveness in preventing disease when $ \mathscr{R}_0 &amp;gt; 0 $. The Jacobian matrix and the Routh-Hurwitz theorem are used to derive the aforementioned analysis techniques. The results are further examined numerically by using the standard second-order Runge-Kutta (RK2) method. In order to visualize the global dynamics of the aforementioned model, the proposed model is expanded to examine some piecewise fractional order derivatives. We may comprehend the crossover behavior in the suggested model's illness dynamics by using the relevant derivative. To numerical present the results, we use RK2 method.&lt;/p&gt;&lt;/abstract&gt;

Numerical Methods That Preserve a Lyapunov Function for Ordinary Differential Equations

The paper studies numerical methods that preserve a Lyapunov function of a dynamical system, i.e., numerical approximations whose energy decreases, just like in the original differential equation. With this aim, a discrete gradient method is implemented for the numerical integration of a system of ordinary differential equations. In principle, this procedure yields first-order methods, but the analysis paves the way for the design of higher-order methods. As a case in point, the proposed method is applied to the Duffing equation without external forcing, considering that, in this case, preserving the Lyapunov function is more important than the accuracy of particular trajectories. Results are validated by means of numerical experiments, where the discrete gradient method is compared to standard Runge–Kutta methods. As predicted by the theory, discrete gradient methods preserve the Lyapunov function, whereas conventional methods fail to do so, since either periodic solutions appear or the energy does not decrease. Moreover, the discrete gradient method outperforms conventional schemes when these do preserve the Lyapunov function, in terms of computational cost; thus, the proposed method is promising.

Two-dimensional Darcy–Bénard convection evolving in Fourier space

Nonlinear transient convection in a porous rectangle heated from below is studied by an analytically based method. Fourier series for the temperature and streamfunction are applied, where each Fourier coefficient evolves in time according to a coupled set of ordinary differential equations. The mathematical method can be considered as a recursive nonlinear mapping in Fourier space from a given state at a time t to an updated state at time t + dt, with an infinitesimal time increment. This nonlinear evolution in Fourier space requires normal-mode compatible boundary conditions along the entire boundary. Each Fourier coefficient gets three contributions during its updating in time: (i) one decay term due to thermal diffusion, governed by linear theory; (ii) one growth term due to buoyancy, governed by linear theory; and (iii) quadratic nonlinearities from the convection term in the heat equation, involving all pairwise interactions between the Fourier modes. We present numerical computations with a standard Runge–Kutta method, with Rayleigh numbers up to ten times the well-known critical value 4π2. Our plots are produced with Fourier series truncated to include 15 normal modes in both the horizontal and vertical directions. For validation of our method, we present tables for the Nusselt number of steady convection, with a higher number (up to 20) of normal modes included in the truncated system of equations. Our computations of transient nonlinear convection lead to steady states. A final steady state is not unique for a given geometry, but depends on the initial state and the Rayleigh number. This ambiguity of steady states is depicted by a hysteresis loop. The Malkus hypothesis of maximal heat transfer is put into perspective. This hypothesis does not pick a preferred cell width, but it nevertheless constrains the hysteresis loop.

Compact High-Order Gas-Kinetic Scheme for Three-Dimensional Flow Simulations

A third-order compact gas-kinetic scheme is proposed for three-dimensional compressible flow computations. The new scheme is based on three key ingredients: the time-accurate gas evolution model for interface flux, the Hermite weighted essentially nonoscillatory reconstruction, and the two-stage temporal discretization. In contrast to the Riemann-solver-based methods, due to the use of time-accurate flux function the proposed scheme can achieve a third-order temporal accuracy with two stages instead of three stages in the standard Runge–Kutta method. As an extension of the existing high-order reconstructions, a Hermite weighted essentially nonoscillatory reconstruction is specifically designed for the current scheme with the implementation of the constrained least-square technique, which subsequently improves the accuracy and robustness of the compact scheme. There is no trouble cell identification in the current scheme. A third-order accuracy can be achieved even with curved boundary. Numerical examples for both smooth and discontinuous flows show the robustness and high accuracy of the compact third-order scheme, which has a comparable performance as the fifth-order noncompact gas-kinetic scheme. A large Courant–Friedrichs–Lewy number around 0.5 can be used in the computations.

Revised version of exponentially fitted pseudo-Runge-Kutta method

In this paper, we have proposed the revised version of exponentially-fitted (ef) pseudo-Runge-Kutta method (ef-PRKM). The motivation behind the revision is to fill the leakage of error in the internal stages during the fitting process. Generally, the internal stage operator, in an ef-PRKM or ef-Runge-Kutta method (ef-RKM) integrates two exponential functions exp(± ωx), with unknown frequency ω e R (for trigonometric-fitting exp(± ωx), ω e iR). However, these internal stage operators produce some amount of errors in integrating other functions like xk exp(± ωx), k ≥ 0 . Here, we first measure the error expression and taking into account this error, we redefine the external stage (solution) operators and compute the revised weights of the ef-PRKM. The revised ef-PRKM is tested on two initial value problems (IVPs). The results are reported in tables and figure. The proposed method would be a good option to find the numerical solution of IVPs efficiently with less cost consumption in the form of slopes/function evaluations than the standard Runge-Kutta methods.

Fractal Density Matrix for the Case of Time Dependent Radio Frequency Pulses Defined in the First Rotating Frame

We treat in this effort the problem of time evolution during a time dependent radio-frequency pulse in the First Rotating Frame (FRF) with Fractal Time Derivatives for the Representation of the Schroedinger Equation and Fractal Density Matrix for Like Spins ½. The resultant time evolutions during a Sin-Cos Pulse are compared with the Time Dependence during the Fractal Bloch equations without relaxation solved using Standard Runge Kutta methods.

On the realization of explicit Runge-Kutta schemes preserving quadratic invariants of dynamical systems

We implement several explicit Runge-Kutta schemes that preserve quadratic invariants of autonomous dynamical systems in Sage. In this paper, we want to present our package ex.sage and the results of our numerical experiments. In the package, the functions rrk_solve, idt_solve and project_1 are constructed for the case when only one given quadratic invariant will be exactly preserved. The function phi_solve_1 allows us to preserve two specified quadratic invariants simultaneously. To solve the equations with respect to parameters determined by the conservation law we use the elimination technique based on Grbner basis implemented in Sage. An elliptic oscillator is used as a test example of the presented package. This dynamical system has two quadratic invariants. Numerical results of the comparing of standard explicit Runge-Kutta method RK(4,4) with rrk_solve are presented. In addition, for the functions rrk_solve and idt_solve, that preserve only one given invariant, we investigated the change of the second quadratic invariant of the elliptic oscillator. In conclusion, the drawbacks of using these schemes are discussed.

Stability analysis of a virulent code in a network of computers

Network of computers are vulnerable to attack by virulent codes which can halt an organization’s activities and in the process resulting in loss of revenue. The need to understand their dynamics if utmost importance and hence there is a need to develop mathematical models. These models will help us understand the impact these virulent codes have on the network of computers. In this article, we develop and solve numerically a mathematical model which can be used to understand the dynamics of a virulent code in a network of computers. This model is called the Immune, susceptible, exposed, infectious, quarantine, and recovered (MSEIQR). The model is solved using a very robust spectral method called the piecewise pseudospectral relaxation method (PPRM). PPRM accuracy is validated by comparing the results with the standard Runge–Kutta method. Stability analysis is also performed on the modified MSEIQR model for malicious code. Results generated are in agreement with the stability analysis performed. Results showing effect of crucial parameters in the dynamics of the model are presented in graphical form.

An Analytical Study for Nonlinear Free and Forced Vibration of Electrostatically Actuated MEMS Resonators

This work aims to construct accurate and simple lower-order analytical approximation solutions for the free and forced vibration of electrostatically actuated micro-electro-mechanical system (MEMS) resonators, in which geometrical and material nonlinearities are induced by the mid-plane stretching, dynamic pull-in characteristics, electrostatic forces and other intrinsic properties. Due to the complexity of nonlinear MEMS systems, the quest of exact closed-form solutions for these problems is hardly obtained for system design and analysis, in particular for harmonically forced nonlinear systems. To examine the simplicity and effectiveness of the present analytical solutions, two illustrative cases are taken into consideration. First, the free vibration of a doubly clamped microbeam suspended on an electrode due to a suddenly applied DC voltage is considered. Based on the Euler–Bernoulli beam theory and the von Karman type nonlinear kinematics, the dynamic motion of the microbeam is further discretized by the Galerkin method to an autonomous system with general nonlinearity, which can be solved analytically by using the Newton harmonic balance method. In addition to large-amplitude free vibration, the primary resonance response of a doubly clamped microbeam driven by two symmetric electrodes is also investigated, in which the microbeam is actuated by a bias DC voltage and a harmonic AC voltage. Following the same decomposition approach, the governing equation of a harmonically forced beam model can be transformed to a nonautonomous system with odd nonlinearity only. Then, lower-order analytical approximation solutions are derived to analyze the steady-state resonance response of such a problem under a combination of various DC and AC voltage effects. Finally, the analytical approximation results of both cases are validated, and they are in good agreement with those obtained by the standard Runge–Kutta method.

A Class of Two-Derivative Two-Step Runge-Kutta Methods for Non-Stiff ODEs

In this paper, a new class of two-derivative two-step Runge-Kutta (TDTSRK) methods for the numerical solution of non-stiff initial value problems (IVPs) in ordinary differential equation (ODEs) is considered. The TDTSRK methods are a special case of multi-derivative Runge-Kutta methods proposed by Kastlunger and Wanner (1972). The methods considered herein incorporate only the first and second derivatives terms of ODEs. These methods possess large interval of stability when compared with other existing methods in the literature. The experiments have been performed on standard problems, and comparisons were made with some standard explicit Runge-Kutta methods in the literature.

A Modified Runge–Kutta Method for Nonlinear Dynamical Systems With Conserved Quantities

In this paper, explicit Runge–Kutta methods are investigated for numerical solutions of nonlinear dynamical systems with conserved quantities. The concept, ε-preserving is introduced to describe the conserved quantities being approximately retained. Then, a modified version of explicit Runge–Kutta methods based on the optimization technique is presented. With respect to the computational effort, the modified Runge–Kutta method is superior to implicit numerical methods in the literature. The order of the modified Runge–Kutta method is the same as the standard Runge–Kutta method, but it is superior in preserving the conserved quantities to the standard one. Numerical experiments are provided to illustrate the effectiveness of the modified Runge–Kutta method.

Factorized Runge–Kutta–Chebyshev Methods

The second-order extended stability Factorized Runge–Kutta–Chebyshev (FRKC2) explicit schemes for the integration of large systems of PDEs with diffusive terms are presented. The schemes are simple to implement through ordered sequences of forward Euler steps with complex stepsizes, and easily parallelised for large scale problems on distributed architectures. Preserving 7 digits for accuracy at 16 digit precision, the schemes are theoretically capable of maintaining internal stability for acceleration factors in excess of 6000 with respect to standard explicit Runge-Kutta methods. The extent of the stability domain is approximately the same as that of RKC schemes, and a third longer than in the case of RKL2 schemes. Extension of FRKC methods to fourth-order, by both complex splitting and Butcher composition techniques, is also discussed.A publicly available implementation of FRKC2 schemes may be obtained from maths.dit.ie/frkc

Order reduction phenomenon for general linear methods

The order reduction phenomenon for general linear methods (GLMs) for stiff differential equations is investigated. It turns out that, similarly as for standard Runge–Kutta methods, the effective order of convergence for a large class of GLMs applied to stiff differential systems, is equal to the stage order of the method. In particular, it is demonstrated that the global error ‖e[n]‖ of GLMs of order p and stage order q applied to the Prothero–Robinson test problem y′(t)=λ(y(t)−φ(t))+φ′(t), t∈[t0,T], y(t0)=φ(t0), is O(hq)+O(hp) as h→0 and hλ→−∞. Moreover, for GLMs with Runge–Kutta stability which are A(0)-stable and for which the stability function R(z) of the underlying Runge–Kutta methods, (i.e., the corresponding RK methods which have the same absolute stability properties as the GLMs), is such that R(∞)≠1, the global error satisfies ‖e[n]‖=O(hq+1)+O(hp) as h→0 and hλ→−∞. These results are confirmed by numerical experiments.

Benchmark of Delayed Equilibrium Model (DEM) and classic two-phase critical flow models against experimental data

The safety analysis of Pressurized Water Reactors, in the event of LOCA, strongly depends on the ability to evaluate the discharge rate of coolant inventory through the breach. Due to the huge pressure difference between the primary system and the reactor containment, the mass flow rate is choked at the break. Under such conditions, both mechanical and thermal equilibrium between phases are not ensured.A general theory to evaluate the two-phase critical mass flow rate is not yet available. However, some models are capable of providing accurate evaluations of either critical mass flux or critical pressure and such a model is the Delayed Equilibrium Model (DEM), which is examined in this article. Here we show how to integrate the DEM system of equations coupling a standard Runge-Kutta method with the Possible-Impossible Flow algorithm. Hence a simple procedure which does not require sophisticate computational schemes.The main objective of this work is to compare DEM, Homogeneous Equilibrium Model, Moody (1965) and Henry and Fauske (1971) models to experimental data. The four models were tested and the results from experimental data containing a sample range in excess of 450 conditions compared in determining an appropriate benchmark. Each of the chosen models is representative of a particular category of critical flow models. Furthermore, two-phase critical models provide good estimations depending on the configuration or set of conditions. Consequently, the models have been individually tested incorporating long tubes, short tubes and slits.This analysis has been carried out for both critical mass flux and critical pressure evaluations.

A Unified IMEX Runge--Kutta Approach for Hyperbolic Systems with Multiscale Relaxation

In this paper we consider the development of Implicit-Explicit (IMEX) Runge-Kutta schemes for hyperbolic systems with multiscale relaxation. In such systems the scaling depends on an additional parameter which modifies the nature of the asymptotic behavior which can be either hyperbolic or parabolic. Because of the multiple scalings, standard IMEX Runge-Kutta methods for hyperbolic systems with relaxation loose their efficiency and a different approach should be adopted to guarantee asymptotic preservation in stiff regimes. We show that the proposed approach is capable to capture the correct asymptotic limit of the system independently of the scaling used. Several numerical examples confirm our theoretical analysis.

Perturbation-Iteration Algorithm to Solve Fractional Giving Up Smoking Mathematical Model

In this paper, a numerical technique is applied to a five variable giving up smoking fractional mathematical model. This model is based on five types of smokers, i.e. potential, occasional, heavy, temporary quitters and permanent quitters. Efficacy of Perturbation Iteration Algorithm on fractional system of differential equations is shown graphically between standard Runge-Kutta method and PIA.