Approximation of Systems with Delay and their Application Mykhaylo Petryk 1, Igor Cherevko 2 and Svitlana Ilika 2 1 Ternopil Ivan Puluj National Technical University, Ruska St., 56. Ternopil, 46001, Ukraine 2 Yuriy Fedkovych Chernivtsi National University, Kotsyubynsky St., 2. Chernivtsi, 58012, Ukraine Abstract In this paper, the schemes of approximation of systems with delays by special systems of ordinary differential equations are considered and the connections between their solutions are investigated. The equivalence of the exponential stability of systems with delay and of the proposed system of ordinary differential equations is established. We consider using approximation schemes to approximate location the asymptotic roots of quasipolynomials of linear systems of differential-difference equations with many delays through the roots of the characteristic equations of the corresponding approximating systems of ordinary differential equations. An algorithm for studying the stability of systems of linear differential-difference equations with many delays is proposed. Numerical algorithms for studying the stability of linear stationary systems with delay are constructed and their coefficient regions of stability are modeled. The implementation of the proposed algorithms for studying the stability of solutions of linear differential equations with delay is demonstrated on model test examples. The analysis of numerical experiments carried confirms the theoretical results presented in the paper. Keywords 1 differential-difference equations, dynamic systems, initial value problems, approximation schemes, stability and asymptotic stability, stability region, computer modeling. 1. Introduction Mathematical models of many applied processes in areas of complex structure are initial and boundary value problems for differential equations with a deviating argument. Such equations have found wide application due to the fact that they well describe many phenomena with consequences in technological processes, radio technical and electrical devices, economic and ecological systems. In many models, delay is introduced as a characteristic of poorly studied processes, which are not taken into account at this stage of model construction. All this became the reason why systems of differential- functional equations are a relevant object for research [1-5]. An important task is to study the stability of solutions of differential-difference equations. In connection with numerous practical applications, considerable attention is paid to obtaining stability conditions for linear differential equations with a delay [6-8]. Of particular interest are studies that allow using the methods of the theory of ordinary differential equations for the analysis of differential functional equations. Schemes of approximation of differential-difference equations by special systems of ordinary differential equations, in particular, make it possible to use approaches to the study of ordinary dynamic systems and to solve a number of problems from the theory of stability of linear systems with a delay. In this work, approximation schemes for initial problems for systems of linear differential-difference equations with many delays are investigated by the sequence of systems of ordinary differential equations, and algorithms for approximate finding of non-asymptotic roots of quasi-polynomials are constructed. Using these algorithms, a method of modeling the stability of solutions of linear systems Dynamical System Modeling and Stability Investigation (DSMSI-2023), December 19-21, 2023, Kyiv, Ukraine EMAIL: mykhaylo_petryk@tntu.edu.ua (A. 1); i.cherevko@chnu.edu.ua (A. 2); s.ilika@chnu.edu.ua (A. 3) ORCID: 0000-0001-6612-7213 (A. 1); 0000-0002-2690-2091 (A. 2); 0000-0002-7773-7655 (A. 3) ©️ 2023 Copyright for this paper by its authors. Use permitted under Creative Commons License Attribution 4.0 International (CC BY 4.0). CEUR Workshop Proceedings (CEUR-WS.org) CEUR Workshop ceur-ws.org ISSN 1613-0073 107 Proceedings with a delay was developed, and constructive computational algorithms were proposed for constructing coefficient regions of stability of linear systems with many delays [8-10]. 2. Approximation schemes Consider the initial problem for a system of differential equations with a delay 𝑑π‘₯(𝑑) = 𝐹(𝑑, π‘₯ (𝑑), π‘₯(𝑑 βˆ’ 𝜏)), 𝑑 β‰₯ 𝑑0 , (1) 𝑑𝑑 π‘₯ (𝑑) = πœ‘(𝑑), 𝑑 ∈ [𝑑0 βˆ’ 𝜏, 𝑇0 ], (2) where 𝐹 (𝑑, 𝑒, 𝑣) is a continuous function that satisfies the Lipshitz condition for 𝑒 + π‘Žπ‘£, πœ‘(𝑑) is a given continuous initial function π‘₯ ∈ 𝑅𝑛 , 𝜏 > 0. An approximate replacement of the initial problem (1)-(2) with the Cauchy problem for the system of ordinary differential equations is proposed in works [11-13]. It is based on the study of a sequence of connected π‘š delay elements and their replacement by a sequence of aperiodic links that carry out a 𝜏 shift by the amount π‘š according to the expansion in the Taylor series. 𝜏 𝜏 π‘§π‘–βˆ’1 (𝑑) = 𝑧𝑖 (𝑑 + ) β‰ˆ 𝑧𝑖 (𝑑) + 𝑧 β€² 𝑖 (𝑑) + β‹― (3) π‘š π‘š The initial problem (1)-(2) is matched by the Cauchy problem for the system of ordinary differential equations 𝑑𝑧0 (𝑑) = 𝐹(𝑑, 𝑧0 , π‘§π‘š ), 𝑑𝑑 (4) 𝑑𝑧𝑗 (𝑑) π‘š = (π‘§π‘—βˆ’1 (𝑑) βˆ’ 𝑧𝑗 (𝑑)) , 𝑗 = Μ…Μ…Μ…Μ…Μ…Μ… 1, π‘š , 𝑑𝑑 𝜏 π‘—πœ 𝑧𝑗 (𝑑0 ) = πœ‘ (𝑑 βˆ’ π‘š ) , 𝑗 = Μ…Μ…Μ…Μ…Μ…Μ… 0, π‘š. (5) Theorem 1 [9]. If the solution of the problem (1)-(2) π‘₯(𝑑) ∈ 𝐢[𝑑0 βˆ’ 𝜏, 𝑇], then π‘—πœ 𝜏 |π‘₯ (𝑑 βˆ’ ) βˆ’ 𝑧𝑗 (𝑑)| ≀ 𝛽 (πœ” (π‘₯, )) , 𝑗 = 0, Μ…Μ…Μ…Μ…Μ…Μ… π‘š , 𝑑 ∈ [𝑑0 , 𝑇], (6) π‘š π‘š 𝜏 where 𝛽(𝛿) β†’ 0 at 𝛿 β†’ 0 πœ” (π‘₯, π‘š) the continuity modulus of the function π‘₯(𝑑) on [𝑑0 βˆ’ 𝜏, 𝑇]. Remark 1. Since the solution of the initial problem (1)-(2) π‘₯(𝑑) ∈ 𝐢[𝑑0 βˆ’ 𝜏, 𝑇], then according to 𝜏 Cantor's theorem on uniform continuity πœ” (π‘₯, ) β†’ 0 when π‘š β†’ 0, which means that the solution of π‘š the Cauchy problem (4)-(5) approximates the solution of the initial problem for a system of differential equations with a delay (1)-(2). 2.1. Approximation nonasymptotic roots of quasipolynomials Consider a linear system of differential-difference equations π‘˜ 𝑑π‘₯ = βˆ‘ 𝐴𝑖 π‘₯ (𝑑 βˆ’ πœπ‘– ), (7) 𝑑𝑑 𝑖=0 where 𝐴𝑖 , 𝑖 = Μ…Μ…Μ…Μ…Μ… 0, π‘˜ A, Bi , i=1,k fixet 𝑛 Γ— 𝑛 matrices, π‘₯ ∈ 𝑅𝑛 , 0 = 𝜏0 < 𝜏1 < β‹― < πœπ‘˜ = 𝜏. Qualitative properties of system (7) depend on the location of the zeros of its quasi-polynomial Ξ¦(πœ†) = 𝑑𝑒𝑑(πœ†πΈ βˆ’ βˆ‘π‘˜π‘–=0 𝐴𝑖 𝑒 βˆ’πœ†πœπ‘– ) = 0. (8) The zeros of the quasi-polynomial (8) are divided into two groups: asymptotic and non-asymptotic. Asymptotic formulas can be obtained for asymptotic roots. Non-asymptotic roots located near the origin of coordinates, determine the main properties of solutions of system (7) and can be found only with the help of approximate algorithms. 108 The study of the approximation schemes of linear systems with a delay showed that the non- asymptotic roots of their quasi-polynomials can be effectively approximated by the roots of the characteristic polynomials of the corresponding approximating systems of ordinary differential equations. According to scheme (4)-(5), we correspond to the system with delay (7) by the system of ordinary differential equations π‘˜ 𝑑𝑧0 (𝑑) πœπ‘– π‘š = βˆ‘ 𝐴𝑖 𝑧𝑙𝑖 (𝑑), 𝑙𝑖 = [ ], 𝑑𝑑 𝜏 𝑖=0 (9) 𝑑𝑧𝑖 (𝑑) π‘š = πœ‡(π‘§π‘–βˆ’1 (𝑑) βˆ’ 𝑧𝑖 (𝑑)), 𝑖 = Μ…Μ…Μ…Μ…Μ…Μ… 1, π‘š , πœ‡ = . 𝑑𝑑 𝜏 Lemma 1 [9,10]. Equality holds for the characteristic equation of the approximating system (9) πœ‡ 𝑙𝑖 Ξ¨π‘š (πœ†) = 𝑑𝑒𝑑 (πœ†πΈ βˆ’ βˆ‘π‘˜π‘–=0 𝐴𝑖 (πœ‡+πœ† ) ) (πœ‡ + πœ†)π‘šπ‘› = 0 (10) and sequence of functions Ξ¨π‘š (πœ†) π»π‘š (πœ†) = (11) (πœ‡ + πœ†)π‘šπ‘› converges as π‘š β†’ ∞ to the quasipolynomial (8). Remark 2. [9, 10]. The zeros of the functions Ξ¨π‘š (π‘₯) and π»π‘š (πœ†) coincide, so the zeros of the characteristic polynomial (10) approximate the non-asymptotic roots of the quasi-polynomial (8). To increase the accuracy of the approximation of non-asymptotic roots, a scheme of increased accuracy is proposed π‘˜ 𝑑𝑧0 (𝑑) = βˆ‘ 𝐴𝑖 𝑧𝑙𝑖 (𝑑) , 𝑑𝑑 𝑖=0 (12) 𝑑𝑧𝑗 (𝑑) = π‘§π‘šβˆ’π‘— (𝑑), 𝑑𝑑 π‘‘π‘§π‘š+𝑗 (𝑑) = 2πœ‡2 (π‘§π‘—βˆ’1 (𝑑) βˆ’ 𝑧𝑗 (𝑑)) βˆ’ 2πœ‡π‘§π‘š+𝑗 (𝑑), Μ…Μ…Μ…Μ…Μ…Μ… 𝑗 = 1, π‘š. 𝑑𝑑 Lemma 2 For the characteristic equation of the approximating system of increased accuracy (2), equality holds π‘š π‘˜ π‘šβˆ’π‘™π‘– πœ†πœ πœ†πœ πœ†πœ πœ†πœ Ξ¨2π‘š+1 (πœ†) = 𝑑𝑒𝑑 [(𝐴0 βˆ’ πœ†πΈ ) (1 + (1 + )) + βˆ‘ 𝐴𝑖 (1 + (1 + )) ]Γ— π‘š 2π‘š π‘š 2π‘š 𝑖=1 π‘šπ‘› πœ†πœ πœ†πœ Γ— (1 + (1 + )) (13) π‘š 2π‘š and sequence of functions Ξ¨2π‘š+1 (πœ†) 𝐻2π‘š+1 (πœ†) = π‘šπ‘› (14) πœ†πœ πœ†πœ (1 + π‘š (1 + 2π‘š )) converges as π‘š β†’ ∞ to the quasipolynomial (8). 109 Example 1. Consider a scalar differential equation with a delay 𝑑π‘₯(𝑑) = π‘Žπ‘₯(𝑑) + 𝑏π‘₯ (𝑑 βˆ’ 𝜏), (15) 𝑑𝑑 where π‘Ž, 𝑏, 𝜏 ∈ 𝑅, 𝜏 > 0. The quasipolynomial for equation (15) has the form Ξ¦(πœ†) = π‘Ž βˆ’ πœ† + 𝑏𝑒 βˆ’πœ†πœ = 0. (16) The approximate system of ordinary differential equations corresponding to scheme (9) for equation (15) 𝑑𝑧0 = π‘Žπ‘§0 + π‘π‘§π‘š , 𝑑𝑑 (17) 𝑑𝑧𝑖 = πœ‡ (π‘§π‘–βˆ’1 βˆ’ 𝑧𝑖 ), Μ…Μ…Μ…Μ…Μ…Μ… 𝑖 = 1, π‘š. 𝑑𝑑 For the characteristic equation of system (17), we obtain equality πœ†πœ π‘š Ξ¨π‘š (πœ†) = (π‘Ž βˆ’ πœ†) (1 + ) + 𝑏 = 0. π‘š πœ†πœ After making the substitution 𝑆 = 1 + π‘š , we obtain the equation. 𝜏 𝜏 𝑆 π‘š+1 βˆ’ (1 + π‘Ž) 𝑆 π‘š βˆ’ 𝑏 = 0, (18) π‘š π‘š which is convenient for calculating its roots. The high-accuracy approximating system of ordinary differential equations corresponding to scheme (12) for equation (15) has the form 𝑑𝑧0 = π‘Žπ‘§0 + π‘π‘§π‘š , 𝑑𝑑 (19) 𝑑𝑧𝑗 = π‘§π‘š+𝑗 , 𝑑𝑑 π‘‘π‘§π‘š+𝑗 = 2πœ‡2 (π‘§π‘—βˆ’1 βˆ’ 𝑧𝑗 ) βˆ’ 2πœ‡π‘§π‘š+𝑗 , 𝑗 = Μ…Μ…Μ…Μ…Μ…Μ… 1, π‘š . 𝑑𝑑 For the characteristic equation of the system (19), the representation is valid π‘š πœ†πœ πœ†πœ Ξ¨2π‘š+1 (πœ†) = (π‘Ž βˆ’ πœ†) (1 + (1 + )) + 𝑏 = 0. (20) π‘š 2π‘š 2π‘š 1 By replacing πœ† = 𝜏 (𝑆 βˆ’ 2) in equation (20), the equation can be obtained in the standard form 𝛼0 𝑆 2π‘š+1 + 𝛼1 𝑆 2π‘š + 𝛼2 𝑆 2π‘šβˆ’1 + β‹― + 𝛼2π‘š 𝑆 + 𝛼2π‘š+1 = 0, where 2π‘š 𝑖 π‘šβˆ’2𝑖 𝛼2𝑖 = βˆ’ 𝐢 2 , Μ…Μ…Μ…Μ…Μ…Μ… 𝑖 = 0, π‘š, 𝜏 π‘š π‘š 𝑖 2π‘šβˆ’2𝑖 𝛼2𝑖+1 = (π‘Ž + ) πΆπ‘š 2 , 𝑖 = Μ…Μ…Μ…Μ…Μ…Μ…Μ…Μ…Μ…Μ…Μ… 0, π‘š βˆ’ 1, 𝜏 π‘š 𝛼2π‘š+1 = (π‘Ž + ) 2βˆ’π‘š + 𝑏. 𝜏 Consider the model equation with a delay 𝑑π‘₯(𝑑) = 2π‘₯ (𝑑) + π‘₯ (𝑑 βˆ’ 1), (22) 𝑑𝑑 the characteristic quasipolynomial of which πœ† = 2 + 𝑒 βˆ’πœ† . (23) The real root of the quasi-polynomial (23) with the largest real part, found by dividing the segment in half, is equal to Ξ»=0.120028. 110 We find the approximation of this root of the quasi-polynomial (23) using the roots of the characteristic polynomials of the approximating systems (18) and (21). The results of numerical experiments for different π‘š are shown in Table 1, where πœ†π‘–π‘ π‘‘. is obtained according to the standard scheme (17), and πœ†π‘–π‘Ž. 𝑖 is obtained according to the scheme of increased accuracy (19). Table 1 Table title π‘š πœ†π‘–π‘ π‘‘. Δ𝑠𝑑. 𝑖 πœ†π‘–π‘Ž. 𝑖 Ξ”π‘–π‘Ž. 𝑖 10 2.143414 0.023386 2.121493 0.001465 20 2.131895 0.011867 2.120422 0.000394 50 2.124817 0.004789 2.120094 0.000066 From Table 1, we see that the scheme of increased approximation accuracy is significantly more efficient than the classical approximation scheme. 2.2. Stability analysis of linear systems with a delay The theory of stability for systems of differential-difference equations (7) is currently sufficiently well researched. Theorem 2 [6, 14-16]. In order for the zero solution of system (7) to be exponentially stable, it is necessary and sufficient that all slopes of its quasipolynomial (8) lie in the left half-plane 𝑅𝑒 πœ† < 0. (24) Thus, the location of the zeros of quasi-polynomials characterizes the stability of linear stationary systems with a delay. In practice, this check is possible only in the simplest cases. We will use Lemma 1 and Lemma 2 to analyze the location of the roots of quasi-polynomials using the characteristic equations of the corresponding approximating systems of ordinary differential equations. Theorem 3 [10]. If the zero solution of the system with delay (1) is exponentially stable (not stable), then there is π‘š0 > 0 such that for all π‘š β‰₯ π‘š0 , the zero solution of the approximating system (3) is also exponentially stable (not stable). If for all π‘š β‰₯ π‘š0 the zero solution of the approximation system (3) is exponentially stable (not stable) then the zero solution of the system with a delay (1) is exponentially stable (not stable). Remark 3. From Theorem 3, we have that there exists a number π‘š0 such that π‘š β‰₯ π‘š0 the asymptotic stability of the zero solution of the system with a delay (7) is equivalent to the asymptotic stability of the zero solution of the approximating system of linear differential equations (9). Consider equality Ξ¦(πœ†) = π»π‘š (πœ†) + π‘…π‘š (πœ†), (25) where π»π‘š (πœ†) – the function defined by equality (11), which at π‘š β†’ ∞ converges to a quasi-polynomial Ξ¦(πœ†), Π° π‘…π‘š (πœ†) = Ξ¦(πœ†) βˆ’ π»π‘š (πœ†). The number π‘š0 should be chosen so that the inequality holds for π‘š β‰₯ π‘š0 min|Ξ¦(πœ†)| > max|π‘…π‘š (πœ†)|. πœ† πœ† Theorem 3 can be used to study the stability of solutions of linear differential-difference equations with a delay. Algorithms of such research are given in works [10, 17]. 111 2.3. Estimation of the effect of delay on the stability of solutions of linear equations with delay In engineering practice, systems with a delay are replaced by systems without a delay, if it is small. Let us consider the method of finding the upper limit of the delay, in which the stability of systems with a delay is equivalent to the stability of the corresponding systems without a delay. We will consider the mathematical justification of the possibility of replacing differential-difference equations with delay by ordinary differential equations is considered, and consider on finding upper bounds of delay, for which the stability regime of systems with delay is analogous to the stability regime of the corresponding systems without delay, is carried out. Consider the linear system with a delay (7) and its corresponding system without a delay π‘˜ 𝑑π‘₯(𝑑) = βˆ‘ 𝐴𝑖 π‘₯ (𝑑). (26) 𝑑𝑑 𝑖=0 Theorem 4 [18,19]. If the zero solution of the system (26) is asymptotically stable, then there exists a constant Ξ” > 0 such that for 0 < 𝜏 < Ξ” the zero solution of the system with a delay is also asymptotically stable. Consider the algorithm for applying the approximation scheme of differential-difference equations (7) by the system of ordinary differential equations (9) to find the upper limit of the delay in system (7), in which it is asymptotically stable [20]: 1. We match the system with a delay (7) to the approximating system (9). 2. We reduce the characteristic polynomial of the approximating system (10) to a form that is convenient for finding its roots. 3. We set the initial limit delay Ξ”0 and a step β„Ž (example, β„Ž = 1). 4. We find the approximate root of the quasi-polynomial 𝛼 with the largest real part 𝛼 at 𝜏 = Ξ”0 . 5. If 𝛼 < 0, then for the given Ξ”0 the zero solution of system with delay (7) is asymptotically stable. We put Δ𝑖+1 = Δ𝑖 + β„Ž, 𝑖 = 0,1, …, and we go to point 4 until this property is fulfilled. 6. If received 𝛼 > 0, then for this upper bound Δ𝑖+1 the zero solution of the delay system (7) is asymptotically unstable. 7. We obtained the interval [Δ𝑖 , Δ𝑖+1 ] of the values of the upper limit of the delay, at the ends of which the asymptotic stability changes. Using the method of division in half, we narrow this interval to the required accuracy and select the left end of the obtained interval for the required value. Example 2. Find the maximum value of Ο„ at which the zero solution of the differential-difference equation 𝑑π‘₯ = 0.25π‘₯ (𝑑) βˆ’ 0.5π‘₯ (𝑑 βˆ’ 𝜏) (27) 𝑑𝑑 will be exponentially stable. At 𝜏 = 0 we obtain a differential equation 𝑑π‘₯ = βˆ’0.25π‘₯ (𝑑), 𝑑𝑑 which is exponentially stable. Using the standard approximation scheme (17) and the increased accuracy approximation scheme (19) for the differential-difference equation (27), we will approximate equation (27) by quasipolynomial (16) with characteristic polynomials of the corresponding approximating systems of ordinary differential equations. Table 2 shows the approximate values of the upper limits of the delay βˆ†, which are found by the algorithm proposed above for different values of the dimensions of the approximating systems m according to the standard approximation scheme and the approximation scheme of increased accuracy. 112 Table 2 Table title M Δ𝑠𝑑. 𝑖 Ξ”π‘–π‘Ž. 𝑖 10 2.5264 2.4142 20 2.4713 2.4173 50 2.4393 2.4182 70 2.4332 2.4183 The exact value of the quantity βˆ† for which the exponential stability of equation (15) is preserved, found by the D-partition method [17] is equal to π‘Ž cos βˆ’1 (βˆ’ 𝑏 ) βˆ— 𝜏 = π‘Ž . 𝑏 sin (cos βˆ’1 (βˆ’ 𝑏 )) In the case of equation (27), we have 𝜏 βˆ— = 2.4181. It follows from Table 2 that the application of approximation schemes of differential-difference equations makes it possible to find an approximate value of the upper limit of the delay for which exponential stability is preserved. At the same time, the scheme of increased accuracy provides better accuracy of approximation. 3. Conclusions The aims to study the schemes of approximation of initial value problems for systems of linear differential equations with many delays by a sequence of systems of ordinary differential equations and their application to the study of the stability of systems of linear differential-difference equations with many delays and finding the upper bound of delay for which the stability of the system with delay is preserved. The study of these problems is reduced to checking the conditions for the negativity of the real parts of all zeros of the corresponding quasipolynomials. Since a direct verification of this condition in practice is possible only in the simplest cases, we analyze the roots of characteristic polynomials of the corresponding approximating systems of ordinary differential equations to solve it. The classical approximation scheme results in algorithms for finding the asymptotic roots of quasipolynomials that are convenient for computer implementation but require a high dimensionality of the corresponding approximating system. To improve the approximation accuracy of the nonasymptotic roots of quasipolynomials, a scheme for increased approximation accuracy is proposed and compared on test model examples. The implementation of the proposed algorithms for studying the stability of solutions of linear differential equations is demonstrated on model test examples. The analysis of numerical experiments carried out for the test model examples confirms the theoretical results presented in the paper. 4. References [1] Corduneanu Π‘., Li Y., Mahdavi M. Functional Differential Equations: Advances and Applications. – John Wiley & Sons, 2016. – 368 p. [2] Schiesser W.E.Time Delay ODE/PDE Models. Applications in Biomedical Science and Engineering. – Boca Rona, 2019. – 250 p. 113 [3] Fathalla A. Rihan. Delay Differential Equations and Applications to Biology. Springer, 2021. – 303 p. [4] Gopalsamy K. Stability and Oscillation in Delay Differential Equations of Population Dynamics. Dordrecht: Kluwer Academic Publishers, 1992. Vol. 74. 501 p. [5] Yang Kuang. Delay differential equations: with applications in population dynamics. New York : Academic Press, 1993. 398 p. [6] Rezvan V. Absolute Stability Automatic System with Delay. Nauka, 1983. 360 p. [7] Kolmanovskii V.B., Myshkis A.D. Introduction to the theory and applications of functional differential equations, Kluwer Acad. Publ. 1999. 664 p [8] Cherevko I., Klevchuk I., Pernay S. Building of stability regions of linear differential-difference equations. Adv. Academy of Sciences of Ukraine. - 2012. - β„– 7. - P. 28-34. [9] Cherevko I., Piddubna L. Approximations of differential difference equations and calculation of nonasymptotic roots of quasipolynomials. Revue D’Analyse numerique et de theorie de l’approximation.- 1999.-28, β„–1 .-P. 15-21. [10] Matviy O.V., Cherevko I.M. About approximation of system with delay and them stability. Nonlinear oscilations .-2004.-7, β„–2.-P. 208-216. [11] Krasovskii N. N. The approximation of a problem of analytic design of controls in a system with time-lag. J. Appl. Math. Mech. 1964.- 28, β„–. 4, P. 876-885. [12] Halanay A. Approximations of delays by ordinary differential equations. Recent advances in differential equations. New York: Academic Press, 1981. P. 155 – 197. [13] Repin Yu. M. On the approximate replacement of systems with lag by ordinary dynamical systems. J. Appl. Math. Mech. 1965. Vol. 29, No. 2, P. 254- 264. [14] Bellman, R., Cooke K. Differential Difference Equations, Academic Press, New York, 1963. 474 p. [15] Halanay A. Differential equations. Stability. Oscillations. Time Lags. New York; London: Acad. Press, 1968. 528 p. [16] Jack K. Hale, Sjoerd M. Verduyn Lunel. Introduction to functional differential equations. Springer New York, 1993. 458 p. [17] Vyzinska I.I. Modeling the stability of differential-difference equations with a delay. Bukovinian Mathematical Journal. 2023. Vol. 11, No. 1. P. 71-79. [18] Qin Yuan-Xun, Liou Iong-Qing, Uang Lian. Effect of time-lags on stability of dynamical systems. 1st International IFAC Congress on automation and Remote Control. 1960. Vol. 1, β„– 1. P. 79 - 94. [19] Khusainov D.Ya., Shatyrko A.V. The Lyapunov function method in studying the stability of differential functional systems. Kyiv: Kyiv University Publishing House, 1987. 236 p. [20] Cherevko I., Tuzyk I., Ilika S., Pertsov A. Approximation of Systems with Delay and Algorithms for Modeling Their Stability. 11th International Conference on Advanced Computer Information Technologies ACIT’2021. Deggendorf, Germany, 15-17 September 2021. P. 49-52. 114