The basic set of test problems for ODE system solvers

Known mathematical packages, MATLAB, Maple, Mathematica, MathCAD and others may get wrong, often plausible, the result of numerical solution of ODE systems with low, given by default, the requirements for mathematical accuracy of the results of the numerical solution of ODE systems. Since ODE system parameters obtained usually experimental, with a low mathematical precision, so the requirements for the precision of the results of mathematical solutions of ODE systems is low (for example, in MATLAB package required precision is 0.001). The article offers a basic set of tests to assess the range of applicability of the relevant solvers. The basic set of test problems for ODE solvers systems include linear ODE systems with the known analytic solution and nonlinear systems with known graphics solution. Presented comparative results of solutions for proposed problems using MATLAB solvers and manzhuk program from a library of standard mathematical programs SADEL (Sets of Algebraic and Differential Equations solvers Library), which has been designed for reliable and accurate solving of systems of linear algebraic equations (LAE) and ODE systems. The results can be used in mathematical modeling of dynamic systems described by ODE systems.


Introduction
The programs for solving ODE systems use three classical problem statements for solving ODE systems: 1. The Cauchy normal form in the space of differential state variables (explicit form) [1]: ⁄ Where is the vector of the coordinate basis of differential variables of the state of dimension ; is the vector function of the right parts of the dimension ; is the independent variable (usually me . Initial conditions specified as , time of integration is , is the time of the beginning the integration, time of completion of integration. 2. Differential-algebraic form of different indices in the full coordinate basis (semi-explicit form) [2,3]: Where -vector of algebraic variables of the full coordinate basis of differential-algebraic variables of dimension , is the vector function of dimension . The agreed initial conditions are now specified as . This statement (1) of the task of ODE systems solutions are  2 widely used in the mathematical and numerical modeling automatic regulation and control systems, for example, MATLAB-SIMULINK [4], MVTU [5,6] and others. 3. The differential-algebraic form in the full coordinate basis, i.e. system of differential-algebraic equations (DAE) (implicit form): Wherevector function of dimension . The initial conditions are the same as in the previous case. Statement (2) is used in mathematical and numerical modeling of electronic and mechanical engineering systems, for example, the ODE systems solvers SPICE and ADAMS. The task (2) in the space of differential algebraic variables was set in the full coordinate basis by L. Petzold in 1982 and she developed the DASS (Differential Algebraic Systems Solver) program based on the Gear method with replacement of derivatives according to the Gear formula [7]:

Where
-Backward Differential Formula,integration step. The resulting system of nonlinear algebraic equations was solved by the Newton method at each integration step. The task of solving ODE systems in the extended space of differential-algebraic variables without transformations of the initial ODE systems is based on the statement (2) of the problem and has the form: (3) Here ⁄vector of derivatives of differential state variables with respect to time of dimension . Considering the fact that some matrix-vector transformations on a computer are not equivalent due to the limited bit grid of computers, the classical forms of representation of ODE systems should be used mainly for theoretical studies, provided that all transformations of the original DAE systems are performed analytically or in symbolic form. The stability domain of any explicit method is limited, therefore, implicit methods are recommended for solving hard problems [1,2]. But implicit methods also have their drawbacks. Using the implicit method for solving stiff and locally unstable problems may lead to qualitatively wrong solution [8]. It is undesirable to use implicit methods for modeling hybrid dynamic systems containing discrete and relay elements, since in this case the Jacobi matrix must be updated frequently. Finally, implementing implicit methods is much more complicated than explicit ones. Therefore, attempts to build explicit methods that are effective for hard problems do not stop [8]. This article discusses the basic set of test ODE systems, which includes both stiff and nonstiff ODE systems and comparing solvers of MATLAB with manzhuk program from the library of standard mathematical programs SADEL in this set of tests.

Methods for solving ODE-DAE systems
There are implemented 8 methods for solving ODE-DAE systems in the MATLAB. The choice of methods for solving ODE-DAE systems in the SADEL library is explained below. Since the explicit methods for solving systems of ODE and DAE are not universal and do not meet the requirements of the stability in the solution with variable step integration of the stiff ODE systems, then as basic integration methods for SADEL library were selected implicit methods of integration [2]. Implicit methods for solving ODE systems are reduced to iteratively solving nonlinear algebraic systems. In the software implementation of these methods, two main types of errors in the numerical solution of ODE systems should be distinguishedglobal error (or qualitative error) and local error (or quantitative error). In theoretical terms, the properties of ODE-DAE systems and methods for their solution are considered in relation to linear inhomogeneous ODE systems of the form: Where is a constant real matrix with size , is a known vector function of time. If the integration method is not suitable for solving the system (4), it is unsuitable for system (3) too. Denote the determinant of the matrix through . The spectrum of eigenvalues is the set of all eigenvalues of matrix , . In the general case is a complex value. It is known that if the eigenvalues are different, then the general solution of system (4) will be a linear combination of fundamental solutions and , is some constant. The analytical solution of system (4) can be represented as the sum of the individual fundamental solutions of a homogeneous system of ordinary differential equations with a matrix and the private solution of the heterogeneous system (4):  The classical control of the local integration error will ensure the accuracy of the solution, and to control the global integration error it is necessary to ensure the appropriate stability of the integration methods. To ensure obtaining a qualitatively correct solution for the above classes of problems, the integration method must be stable for stable ODE systems (I and II classes of tasks) and should be unstable for unstable ODE systems (IV and V classes of problems), therefore it is necessary that the integration methods be appropriately AL-stable (AL-stability -Absolute stability in Left half-plane of the complex stability plane of numerical integration), i.e. absolutely stable strictly in the left half-plane of the complex plane of stability of the integration methods [9] . With low requirements to the precision of solving ODE systems, this area of absolute stability guarantees the preservation of a qualitatively correct solution with increasing integration step in the plane , since with low integration accuracy the integration step can become so large that non-ALstable integration methods can obtain incorrect solution trajectories. Implicit one-step Runge-Kutta integration methods meet the AL-stability requirement [10]. Algorithms for the implementation of implicit one-step integration methods of Runge-Kutta type are based on the joint solution of systems of nonlinear algebraic equations for solved DAE system of the form (3) and systems of linear algebraic equations for integration methods , formed at the corresponding stages of the integration methods, relative to differential-algebraic variables , and derivatives of differential variables [11]. Were developed DAbc formulas (as the development of Ab T c Butcher methods [12]) for S-phasic implicit methods for solving DAE systems in an extended coordinate basis: Wherethe number of stages of the corresponding one-step method ; --th integration step; , , , method parameters. Based on these formulas, three integration methods were developed and implemented:  M1 -A-stable one-step implicit method of the first order of accuracy;  M2 -AL-stable one-step implicit method of the second order of accuracy;  M3 -AL-stable two-stage implicit method of the fourth order of accuracy. The main drawback of AL-stable methods is that they do not satisfy the requirements of Lstability and sometimes give "false fluctuations" in the solution, so these methods are not widely used in practice [2]. The solution to this problem is given in [10].

Test problems with a well-known analytical solution
It is proposed to develop 5 test linear homogeneous ODE systems 1) -5) with analytical solutions corresponding to the classes of tasks shown in Figure 1 as the basic set of test problems. Given the test linear system of ODE of the 2nd order:

⁄
Where is the vector of unknowns of system , is the matrix of constant coefficients of a system of size 2 2, is a constant vector.

Test problems with a graphical solution:
6) An example of the calculation of a tough second-order ODE system (MU is a stiffness parameter) -Van der Pol's test. Correct solution for this test shown in Figure 2.

7)
Example of electronic circuit simulation with a multi-period solution (5th order ODE system). Figure 3. Electronic circuit for test 7 The ODE system for this test is:   Waivers for problems with graphical solution are time of calculation more than 10 seconds or the wrong plot of the solution.

Test problems with a well-known analytical solution :
Failures for tasks with an analytical solution were recorded if the maximum relative error was more than 0.001 in a given integration interval. These tasks were used to test manzhuk solver and different MATLAB solvers. The final test results are shown in Table 1.

Test problems with a graphical solution
In view of the results published in [10,13,14], Table 2 shows the summary comparison solvers of stiff ODE systems of the aforementioned SADEL library, from well-known libraries and packages of mathematical programs. The real laser parameters manzhuk (М3 The sign (+) means the correct solution of the test using the method recommended for solving stiff ODE systems, with an accuracy . The sign (-) indicates an incorrect test solution using the method recommended for solving hard ODE systems with default parameters without warning the user and refusing to get the right solution. The sign (-+) means an incorrect test solution by the method recommended for solving hard ODE systems, with default parameters without warning the user and obtaining the correct solution by setting appropriate methods and integration parameters and increasing the calculation time. Sign (-++) means an incorrect test solution by the method recommended for solving hard ODE systems, with default parameters without warning the user and obtaining the correct solution by setting only the parameter specified accuracy of the method of integration and increasing the calculation time. NDno data. Figure 6 shows as an example the result of a qualitatively incorrect solution for the 1st variable of test 7) by the implicit method ode15s and correct solution by method M3 of program manzhuk. There were no warnings about the incorrect decision in the ode15s solver.

Conclusion
Most of the problems of solving systems of ordinary differential equations arising in practice cannot be solved analytically [15]. Numerical methods descend to the seventeenth centuryto Newton, who proposed using them to calculate the trajectories of comets. Manual calculations continued to be used until the 60s of the twentieth century, before the appearance of the first solver programs written in FORTRAN for the first computers [16]. Since that time, development of test problems for solving ODE systems has been carried out. Classical sets of test considered tests from books [1,2,17]. The first test comparison of the MATLAB ODE solvers is given in [18]. In MATLAB we can always choose the method and adjust its parameters to obtain qualitatively correct results, but most users are not experts in numerical methods, so the comparison should be made for the recommended default parameters of integration ODE solvers. For specific classes of problems to be solved, it is necessary to have sets of test practical and mathematical problems with a known deliberately reliable and accurate solution for the correct selection and adjustment of the corresponding solvers of ODE systems and the correct evaluation of the obtained results of numerical simulation. This article presents the results obtained in the continuation of the above works. Further studies focused on developing effective methods and algorithms for the automatic integration step selection for AL-stable implicit methods for the numerical solution of high-and ultra-high-dimensional ODE-DAE systems and testing these new solvers. The relevance and importance of these works is indicated in the [19,20].