# Stability and Boundedness of Solutions to Some Multidimensional Time-Varying Nonlinear Systems

Authors: Mark A. Pinsky
Journal: Mathematical Problems in Engineering
 (2022)

---

## Abstract

Assessment of the degree of boundedness/stability of multidimensional nonlinear systems with time-dependent and nonperiodic coefficients is an important problem in various applied areas which has no adequate resolution yet. Most of the known techniques provide computationally intensive and conservative stability criteria in this field and frequently fail to estimate the regions of boundedness/stability of solutions to the corresponding systems. Recently, we outlined a new approach to this problem which is based on the analysis of solutions to a scalar auxiliary equation bounded from the above time histories of the norms of solutions to the original system. This paper develops a novel technique casting the auxiliary equation in a modified form which extends the application domain and reduces the computational burden of our prior approach. Consequently, we developed more general boundedness/stability criteria and estimated trapping/stability regions for some multidimensional nonlinear systems with nonperiodic time-dependent coefficients that are common in various application domains. This lets us to assess in target simulations the extent of boundedness/stability of some multidimensional, nonlinear, and time-varying systems which were intractable with our prior technique.


---

## Body

## 1. Introduction and Motivation Example

### 1.1. Background

Analysis of boundedness and stability of nonlinear systems with variable and nonperiodic coefficients plays a vital role in various engineering and natural science problems, which, for instance, are concerned with the design of controllers and observers. It appears that currently there are no confirmable necessary and sufficient conditions of local stability of the trivial solution to homogeneous systems of such kind (see, e.g., [1–7]). It was shown by Lyapunov [8] that under some additional conditions, the trivial solution to a homogeneous system with time-varying coefficients is asymptotically stable if the linearization of this system at zero is regular and its maximal Lyapunov exponent is negative (see also [1] for contemporary review of this subject). Yet, the confirmation of the former condition presents a considerable problem in applications. Subsequently, it was demonstrated by Perron [9] that arbitrary small perturbations can reverse sign of the Lyapunov exponents of a linearized system with time-varying coefficients and alter its stability. This attested that Lyapunov’s regularity condition is essential. Consequent works were focused on examining the Lyapunov stability conditions through a review of the stability of Lyapunov exponents for linearized nonautonomous systems. Finally, the necessary and sufficient conditions of stability of Lyapunov exponents for solutions to the linearized system were given in [10, 11]. However, it turns out that authentication of these conditions requires the apprehension of the fundamental set of solutions for the underlying systems, which is dubious in applications.Application of the concept of generalized exponents, that was introduced in [6], provides sufficient conditions of stability of the trivial solution to time-varying nonlinear systems which, in principle, can be checked with the aid of numerical simulations. Still, the upper generalized exponent is larger than or equal to the maximal Lyapunov exponent, which heightens the conservatism of this more robust approach.To retrieve a simple sufficient stability condition from the last approach, we firstly write the following equation with marked linearization:(1)x˙=Btx+f∗t,x+F∗t,∀t≥t0,xt,t0,x0∈ℝn,f∗t,0=0,xt0,t0,x0=x0,where functions, f∗:t0,∞×ℝn⟶ℝn, F∗:t0,∞⟶ℝn, and matrix, B:t0,∞⟶ℝn×n, are continuous, t0∈Τ:=ζ,∞,ζ∈ℝ, x0∈H⊂ℝn, and H is a neighborhood of x≡0.ℝn is a real n−dimensional Euclidean space, F∗t=F0ηt, supt≥t0ηt=1,F0>0, ⋅ stands for induced 2-norm of a matrix or 2-norm of a vector, and xt,t0,x0:t0,∞×Τ×H⟶ℝn is a solution to (1). To shorten the notation, we will write below that xt,t0,x0≡xt,x0 and assume that (1) possesses a unique solution ∀x0∈H and ∀t≥t0.Due to continuity of the right side of (1), the last assumption implies that xt,x0 is a continuous and continuously differential function that is bounded for any t∈t0,t∗,t∗<∞. Consequently, the remainder of this paper focuses on behavior of xt,x0 for t⟶∞.We also will assess stability of the trivial solution to the following homogeneous equation:(2)x˙=Btx+f∗t,x,xt0,x0=x0,and its linear counterpart:(3)x˙=Btx,xt0,x0=x0.Clearly, the solution to (3) can be written as xt,x0=Wt,t0x0, where Wt,t0=wtw−1t0 and wt are transition and fundamental solution matrices for (3), respectively [1, 4].Next, we write Lipschitz continuity condition forf∗t,x as follows:(4)f∗t,x≤l1t,t0x,∀x∈Ω1⊂ℝn,∀t≥t0,where Ω1 is a bounded neighborhood of x≡0 and l1t,t0≤l^1>0,∀t≥t0, is a continuous function. In turn, let us assume that(5)Wt,t0≤Ne−vt−t0,∀t≥t0,N,v>0,and(6)Nl^1−v<0.Then, the trivial solution to (2) is asymptotically stable if (4)–(6) hold [4, 6, 7].A more general but less tractable condition on the norm of the transition matrices was presented in [1, 2] in the form(7)Wt,t0≤Nexp∫t0tϒsds,∀t>t0,ϒ:t0,∞⟶ℝ,N>0,where it was shown that (4) and (7) and the following condition:(8)limt⟶∞sup1t−t0∫t0tϒsds+Nl^1<0,embrace the asymptotic stability of the trivial solution to (2). More compelling stability conditions of (2) were given in [12], yet verification of these conditions in applications can be challenging as well.In the control literature, the analysis of stability of nonautonomous nonlinear systems was frequently aided by application of the Lyapunov function method [15–23]. Nonetheless, the application of this methodology is strenuous in this area since adequate Lyapunov functions are rarely known for multidimensional and nonlinear systems with nonperiodic coefficients.The problem of estimating the stability regions of autonomous nonlinear systems has attracted numerous publications in the last few decades [24–38], but these techniques fail for systems with time-varying coefficients. To our knowledge, the problem of estimating the trapping regions for nonautonomous nonlinear systems has not been virtually addressed in the current literature.Assessment of stability of nonlinear systems based on the analysis of convergence of adjacent trajectories was developed in [39]. However, to our knowledge, this approach rarely leads to computationally sound stability conditions for multidimensional systems.In [13], we developed a novel technique for estimation of upper bounds of the norms of solutions to (1) or (2) and use it to distinguish the boundedness/stability criteria and estimate the trapping/stability regions for these systems. Consequently, under normalizing condition, wt0=1, we derive the following inequality: xt,x0≤Xt,X0,∀t≥t0, where xt,x0 is a solution to (1) and Xt,X0:t0,∞×ℝ≥0⟶ℝ≥0 (ℝ≥0 is a set of nonnegative real numbers) is a solution to a scalar equation we derived in [13]:(9)X˙=ptX+ctf∗t,xt,x0+F∗t,Xt0,X0=w−1t0x0=X0.Note that(10)pt=dlnwtdt,and(11)ct=wtw−1t=σmaxwσminw.In turn, Pinsky and Koblik [13] also introduced a nonlinear extension of the Lipschitz continuity condition which was presented as follows:(12)f∗t,x≤Lt,x,∀x∈Ω2x⊂ℝn,∀t≥t0,where L:=t0,∞×ℝ≥0⟶ℝ≥0 is continuous function in t and x and locally Lipschitz in x, Lt,0=0, and Ω2x is a bounded neighborhood of x≡0. Note that Pinsky and Koblik [13] defined (12) in a closed form if f∗ is either a piece-wise polynomial function in x or can be approximated by such function, e.g., bounded in Ω2x error term, which, for instance, can be written in Lagrange form. In the former case, Ω2x≡ℝn and in the latter case, Ω2x≡ℝn if the error term is also globally bounded for ∀x∈ℝn. Thus, condition Ω2x≡ℝn holds for a large set of nonlinear systems emerging in science and engineering applications and we are going to use it in the remainder of this paper.Let us illustrate how to defineLt,x for a simple but representative vector function. Assume that x=x1x2T and a vector function f∈ℝ2 is defined, e g., as follows: f=a1tx13x22a2tx22T. Then, f2≤f1≤a1x13x22+a2x22≤a1x5+a2x2, where we use that xin≤xn,i=1,2,n∈ℕ, and ℕ is a set of positive integers. Note that additional and more complex examples of this kind are provided in Section 5 (see also [13, 14]).Consequently, using (12), we reviewed (9) in the following more tractable form:(13)z˙=ptz+ctLt,z+F∗t,zt0,z0=z0=w−1t0x0.To reduce the notation, we adopt in (13) and throughout this paper that zt,t0,z0=zt,z0.In turn, it was shown in [13] that(14)xt,x0≤Xt,X0≤zt,z0,where xt,x0,Xt,X0, and zt,z0 are solutions of (1), (9), and (13), respectively.Thus, (14) can be used to bound from above the norm of solutions of (1) by matching solutions of the scalar equation (13).It turned out that (13) provides sound estimates of the trapping/stability regions under the assumption that only the bound of f∗t,x is known—a frequent pronouncement in the control literature [3–5]. Nonetheless, such estimates may become more conservative if f∗t,x is defined explicitly. Consequently, we refined this methodology in [14], where (13) was used to estimate the error of successive approximations of solutions to (1) or (2) stemming from the trapping/stability regions of these systems. This modified approach enhanced our boundedness/stability criteria and delivered approximations that increased successively the accuracy of estimation of the boundaries of the corresponding trapping/stability regions.Nonetheless, the methodologies developed in [13, 14] work under the condition that ct<∞,∀t≥t0, which considerably limits the scope of its applications. Nonetheless, it follows from (11) that frequently limt⟶∞ct=∞ even if A is a stable and time-invariant matrix.The current paper lifts this limitation for a practically important class of nonlinear systems with time-varying and nonperiodic coefficients. Its main contribution is in the development of a modified auxiliary equation withct≡1 and pt≡const under some conditions that are frequently met in various applications. This new technique prompts the development of novel criteria of bondedness/stability and estimation of bondedness/stability regions for a wide class of systems that were intractable to our former methodology [13] and voids elaborate simulations of wt, pt, and ct that were required previously.To make this paper more inclusive, we present the definitions of some standard principles of stability theory which are going to be used in the remainder of this paper. Note that the standard definitions of Lyapunov stability and asymptotic stability for time-varying nonlinear systems that are accepted below can be found, e g., in [1–5]. In the remainder of this paper, we will call these properties shortly either stability or asymptotic stability.We begin with a formal definition of Lyapunov exponents, where we adopt the exposition of these quantities made in [1, 5].For (2), the Lyapunov exponents ϕixt,x0=limt0<t⟶∞supt−t0−1lnxt,x0,i≤n, where xt,x0 is a solution to (2). The maximal Lyapunov exponent, ϕ=maxiϕixt,x0, measures the maximal rate of exponential growth/decay of the corresponding solutions as t⟶∞ and bears a pivoting role in stability theory. For a linear system (3), the Lyapunov exponents are defined as ϕl,i=limt0<t⟶∞supt−t0−1lnσiwt,i≤n, where σi are the singular values of the fundamental solution matrix of (3), and for linear systems, maxiϕl,i=ϕ [1]. This lets us to bound the transition matrix of (3) as follows:(15)Wt,t0≤Dϑe−χt−t0,∀t≥t0,D>0,−χ=ϕ+ϑ,where ϑ is a small positive number [6].In order, let us bring the definition of the comparison principle [4] which is frequently used below. Consider a scalar differential equation(16)u˙1=gt,u1,∀t≥t0,u1t,u10∈ℝ,u1t0,u10=u10,where function gt,u1 is continuous in t and locally Lipschitz in u1 for ∀u1∈℘⊂ℝ. Suppose that the solution to this equation is u1t,u10∈℘,∀t≥t0. Next, consider a differential inequality(17)D+u2≤gt,u2,∀t≥t0,u2t0,u20=u20≤u10,where D+u2 denotes the upper right-hand derivative in t of u2t,u20 [4] and u2t,u20∈℘,∀t≥t0. Then, u1t,u10≥u2t,u20,∀t≥t0.Now we present for convenience some conventional definitions of the trapping/stability regions as follows.Definition 1.
A connected and compact set of all initial vectors,ℑ1t0, is called a trapping region of equation (1) if condition x0∈ℑ1t0 implies that xt,x0∈ℑ1t0,∀t>t0.
Clearly, this definition acknowledges thatℑ1 is the invariant set of (1) containing zero.Definition 2.
A connected and open set of all initial vectors,ℑ2t0, that includes zero vector, is called a region of stability of the trivial solution to (2) if condition x0∈ℑ2t0 implies that xt,x0 is stable.Definition 3.
A connected and open set of all initial vectors,ℑ3t0, that includes zero vector, is called a region of asymptotic stability of the trivial solution to (2) if condition x0∈ℑ3t0 implies that limt⟶∞xt,x0=0.


### 1.2. Motivation Example

The first part of this section derives a novel auxiliary equation withct=1 and pt=const for a simple planar system. The second part infers stability and estimates the stability region of the derived auxiliary equation and extends these inferences to the stability assessment of our planar model (2).Let us firstly apply our former methodology [13] to a simple case, where B=diagλ1λ2, λ1,2∈ℝ,λ1>λ2. Obviously, in this case, wt=diageλ1t−t0eλ2t−t0,w−1t=diage−λ1t−t0e−λ2t−t0 which implies that wt=eλ1t−t0,w−1t=eλ2t−t0. Thus, ct=wtw−1t=eλ1−λ2t−t0 and limt⟶∞ct=∞ if λ1≠λ2 which makes our prior estimates of xt,x0, which is based on (13), overconservative for large values of t. Nonetheless, ct≡1 if λ1,2 are complex conjugates and our prior inferences hold in this case.Let us break out our current approach into a sequence of straightforward steps for a simple version of system (1) where x∈ℝ2, Bt=A+G∗t, A∈ℝ2×2,A=const,A≠0, G∗t∈ℝ2×2, and limt⟶∞t−t0−1∫t0tG∗sds=0. Assume also that matrix A is diagonalizable and eigenvalue and eigenvector matrices of A are Λ=diagλ,λ∗, V=vkj,k,j=1,2, where λ=α+iβ,λ∗=α−iβ, f∗=0εx23T, F∗=0F∗2T.(13) can be applied to the stability analysis of such a planar system if (11) implies that ct<∞,∀t≥t0 which, in general, is difficult to warranty in advance. Thus, we outline a different technique recasting the auxiliary equation in the form, where ct≡1 and pt≡const.(1)
Firstly, let us write our planar model of system (1) in eigenbasis of A as follows:(18)y˙=Λy+Gty+V−10−ε∑j=12v2jyj3T+F2t,wherey∈ℂ2, ℂ2 is a two-dimensional space of complex numbers, x=Vy, Gt=V−1G∗tV, and Ft=V−1F∗t.(2)
Next, we rewrite the last equation as follows:(19)y˙=λI+iby+a−λIy+Gty+V−10−ε∑j=12v2jyj3T+Ft,whereI is an identity matrix, b=diagβ−β, a=diagαα, and λ∈ℝ is going to be determined subsequently.In order, let us select a linear subsystem of (19) with an underlined diagonal matrix as(20)y˙=λI+iby.For this subsystem, a diagonal fundamental solution matrix can be written as follows:wt=eλI+ibt−t0=eλt−t0eibt−t0, which implies that wt=eλt−t0eibt−t0=eλt−t0 since eibt−t0=1. Same reasoning yields that w−1t=e−λt−t0 which shows that ct=wtw−1t=1, and due to (10), pt=λ=const.Thus, if we use (18) as the underlying linear system, then the first addition in right side of the auxiliary equation (13), which is derived for (17), is ptz=λz.(3)
Further application of (13) to (17) brings the second term in the right side of our modified auxiliary equation in the following form:(21)a−λIy+Gty+V−10−ε∑j=12v2jyj3T+Ft≤a−λI+Gty+V−10−ε∑j=12v2jyj3T+Ft.Then, the utility of standard inequality yields(22)0−ε∑j=12v2jyj3T≤absε∑j=12v2jyj3≤εabs∑j=12v2jyj3≤ε∑j=12absv2jabsyj3≤ε∑j=12absv2jy3≤ε∑j=12absv2j3y3,whereabsq=a2+b2 if q=a+ib and absq=q if Imq=0.This lets us to write (13) as follows:(23)z˙=λ+a−λIz+Gtz+σz3+Ft,z0=V−1x0,whereσ=εV−1∑j=12absv2j3.(4)
Lastly, we selectλ to maximize the degree of stability of a linear equation z˙=λ+a−λIz which minimizes conservatism of our estimates and reduces to the following condition.(24)minλλ+a−λI=minλλ+α−λ,since a is a diagonal matrix with equal eigenvalues. Resolving (24) yields that λ=α which casts (23) in the following form:(25)z˙=α+Gtz+σz3+Ft,∀t≥t0,zt0,z0=z0=V−1x0.For further references, we write a homogeneous counterpart to (25) as follows:(26)z˙=α+Gtz+σz3,∀t≥t0,z0=V−1x0.Let us recall now that, due to (14), xt,x0≤zt,z0,z0=V−1x0,∀t≥t0, where xt,x0 and zt,z0 are solutions to planar versions of either equation (1) or (2) and equation (25) or (26), respectively.To further abridge the stability analysis of the trivial solution to (26), we are going to develop for this equation its linear, scalar, and, thus, integrable upper bound. Due to the comparison principle [4], solutions to such equations resolve the essential boundedness/stability properties of equation (26). In general, Lipschitz continuity condition (4) can be used to develop such a linear equation. Nonetheless, in our case, we apply a less conservative bound, z3≤z^2z,z∈0,z^,z^>0. Application of the last inequality brings the following linear equation:(27)Z˙=st,z^Z,∀t≥t0,Zt0,z0,z^=z0=V−1x0,where st,z^=α+Gt+σz^2. Note that solutions to (27) can be used to estimate solutions to (26) if z∈0,z^.In turn, to simplify the assessment of stability of (27), we set that supt0≤tst,z^=−α∗+z^2σ+Gs=−s¯z^, where we assumed that α=−α∗, α∗≥0 and Gs=sup∀t≥t0Gt<∞.Subsequently, the solutions to (26) can be bounded as follows: Z=z0exp∫t0tsτ,z^dτ≤z0exp−∫t0ts¯z^dτ≤z0exp−s¯z^t−t0 which implies that (26) is asymptotically stable if(28)s¯z^>0,z^∈0,z^max,where z^max is the maximal value of z^ for which (28) holds.The last assessment can be drawn through the utility of a more general reasoning which is repeatedly applied in this paper. In fact, determine a comparison equation to (27) as follows:(29)Z¯˙=−s¯z^Z¯,∀t≥t0,Z¯t0,z0,z^=z0=V−1x0.Clearly, the right side of (29) bounds from above the right side of (27), whereas the initial values for both equations are equal. Thus, Zt,z0,z^≤Z¯t,z0,z^=z0exp−s¯z^t−t0,∀t≥t0, due to comparison principle.Next, we infer that (28) implies the asymptotic stability of the trivial solution to (26). In fact, due to the comparison principle,(30)zt,z0≤Zt,z0,z^,∀t≥t0,if zt,z0∈0,z^,∀t≥t0, where zt,z0 and Zt,z0,z^ are solutions to (26) and (27).Let us show thatzt,z0∈0,z^,∀t≥t0 if z0<z^ and (28) holds. Clearly, zt0,z0=z0∈0,z^ under our assumption. Assume, on the contrary, that t=t∗>t0 is the smallest value of t for which zt∗,z0=z^. Then, (30) holds for ∀t∈t0,t∗, and thus zt∗,z0≤Zt∗,z0,z^≤Z¯t∗,z0,z^=z0exp−s¯z^t∗−t0<z0<z^. Hence, zt,z0<z^,∀t≥t0 if (28) holds and z0<z^. Consequently, the trivial solutions to (25) and, in turn, to our planar model of equation (2) are asymptotically stable.Furthermore, solving forz^ equation s¯z^=0 trumps estimation of the range of values of z^ which assure asymptotic stability of (29), i.e., values of z^∈0,z^max, where(31)z^max=α∗−Gsσ12.The last formula leads to the estimation of the stability regions for (26) and the corresponding model of equation (2). Really, (28) and, subsequently, (30) hold if z0<z^∈0,z^max. Hence, the last inequality implies that (29) and, in turn, (27) are asymptotically stable under that condition which, in turn, yields that the trivial solution to (26) is asymptotically stable as well, and a solution to (26) approaches zero as t⟶∞ if z0∈0,z^max. In order, this condition ensures asymptotic stability of the trivial solutions to our simplified model of (2) which relates to (26). Hence, due to the second relation in (26), the region of asymptotic stability of the trivial solution to our planar model of (2) includes interior points of the ellipsoid, which is defined as follows: V−1x0=z^max.Section3 generalizes this last formula and derives the boundedness/stability criteria and estimates the trapping/stability regions for multidimensional (1) and (2).


## 1.1. Background

Analysis of boundedness and stability of nonlinear systems with variable and nonperiodic coefficients plays a vital role in various engineering and natural science problems, which, for instance, are concerned with the design of controllers and observers. It appears that currently there are no confirmable necessary and sufficient conditions of local stability of the trivial solution to homogeneous systems of such kind (see, e.g., [1–7]). It was shown by Lyapunov [8] that under some additional conditions, the trivial solution to a homogeneous system with time-varying coefficients is asymptotically stable if the linearization of this system at zero is regular and its maximal Lyapunov exponent is negative (see also [1] for contemporary review of this subject). Yet, the confirmation of the former condition presents a considerable problem in applications. Subsequently, it was demonstrated by Perron [9] that arbitrary small perturbations can reverse sign of the Lyapunov exponents of a linearized system with time-varying coefficients and alter its stability. This attested that Lyapunov’s regularity condition is essential. Consequent works were focused on examining the Lyapunov stability conditions through a review of the stability of Lyapunov exponents for linearized nonautonomous systems. Finally, the necessary and sufficient conditions of stability of Lyapunov exponents for solutions to the linearized system were given in [10, 11]. However, it turns out that authentication of these conditions requires the apprehension of the fundamental set of solutions for the underlying systems, which is dubious in applications.Application of the concept of generalized exponents, that was introduced in [6], provides sufficient conditions of stability of the trivial solution to time-varying nonlinear systems which, in principle, can be checked with the aid of numerical simulations. Still, the upper generalized exponent is larger than or equal to the maximal Lyapunov exponent, which heightens the conservatism of this more robust approach.To retrieve a simple sufficient stability condition from the last approach, we firstly write the following equation with marked linearization:(1)x˙=Btx+f∗t,x+F∗t,∀t≥t0,xt,t0,x0∈ℝn,f∗t,0=0,xt0,t0,x0=x0,where functions, f∗:t0,∞×ℝn⟶ℝn, F∗:t0,∞⟶ℝn, and matrix, B:t0,∞⟶ℝn×n, are continuous, t0∈Τ:=ζ,∞,ζ∈ℝ, x0∈H⊂ℝn, and H is a neighborhood of x≡0.ℝn is a real n−dimensional Euclidean space, F∗t=F0ηt, supt≥t0ηt=1,F0>0, ⋅ stands for induced 2-norm of a matrix or 2-norm of a vector, and xt,t0,x0:t0,∞×Τ×H⟶ℝn is a solution to (1). To shorten the notation, we will write below that xt,t0,x0≡xt,x0 and assume that (1) possesses a unique solution ∀x0∈H and ∀t≥t0.Due to continuity of the right side of (1), the last assumption implies that xt,x0 is a continuous and continuously differential function that is bounded for any t∈t0,t∗,t∗<∞. Consequently, the remainder of this paper focuses on behavior of xt,x0 for t⟶∞.We also will assess stability of the trivial solution to the following homogeneous equation:(2)x˙=Btx+f∗t,x,xt0,x0=x0,and its linear counterpart:(3)x˙=Btx,xt0,x0=x0.Clearly, the solution to (3) can be written as xt,x0=Wt,t0x0, where Wt,t0=wtw−1t0 and wt are transition and fundamental solution matrices for (3), respectively [1, 4].Next, we write Lipschitz continuity condition forf∗t,x as follows:(4)f∗t,x≤l1t,t0x,∀x∈Ω1⊂ℝn,∀t≥t0,where Ω1 is a bounded neighborhood of x≡0 and l1t,t0≤l^1>0,∀t≥t0, is a continuous function. In turn, let us assume that(5)Wt,t0≤Ne−vt−t0,∀t≥t0,N,v>0,and(6)Nl^1−v<0.Then, the trivial solution to (2) is asymptotically stable if (4)–(6) hold [4, 6, 7].A more general but less tractable condition on the norm of the transition matrices was presented in [1, 2] in the form(7)Wt,t0≤Nexp∫t0tϒsds,∀t>t0,ϒ:t0,∞⟶ℝ,N>0,where it was shown that (4) and (7) and the following condition:(8)limt⟶∞sup1t−t0∫t0tϒsds+Nl^1<0,embrace the asymptotic stability of the trivial solution to (2). More compelling stability conditions of (2) were given in [12], yet verification of these conditions in applications can be challenging as well.In the control literature, the analysis of stability of nonautonomous nonlinear systems was frequently aided by application of the Lyapunov function method [15–23]. Nonetheless, the application of this methodology is strenuous in this area since adequate Lyapunov functions are rarely known for multidimensional and nonlinear systems with nonperiodic coefficients.The problem of estimating the stability regions of autonomous nonlinear systems has attracted numerous publications in the last few decades [24–38], but these techniques fail for systems with time-varying coefficients. To our knowledge, the problem of estimating the trapping regions for nonautonomous nonlinear systems has not been virtually addressed in the current literature.Assessment of stability of nonlinear systems based on the analysis of convergence of adjacent trajectories was developed in [39]. However, to our knowledge, this approach rarely leads to computationally sound stability conditions for multidimensional systems.In [13], we developed a novel technique for estimation of upper bounds of the norms of solutions to (1) or (2) and use it to distinguish the boundedness/stability criteria and estimate the trapping/stability regions for these systems. Consequently, under normalizing condition, wt0=1, we derive the following inequality: xt,x0≤Xt,X0,∀t≥t0, where xt,x0 is a solution to (1) and Xt,X0:t0,∞×ℝ≥0⟶ℝ≥0 (ℝ≥0 is a set of nonnegative real numbers) is a solution to a scalar equation we derived in [13]:(9)X˙=ptX+ctf∗t,xt,x0+F∗t,Xt0,X0=w−1t0x0=X0.Note that(10)pt=dlnwtdt,and(11)ct=wtw−1t=σmaxwσminw.In turn, Pinsky and Koblik [13] also introduced a nonlinear extension of the Lipschitz continuity condition which was presented as follows:(12)f∗t,x≤Lt,x,∀x∈Ω2x⊂ℝn,∀t≥t0,where L:=t0,∞×ℝ≥0⟶ℝ≥0 is continuous function in t and x and locally Lipschitz in x, Lt,0=0, and Ω2x is a bounded neighborhood of x≡0. Note that Pinsky and Koblik [13] defined (12) in a closed form if f∗ is either a piece-wise polynomial function in x or can be approximated by such function, e.g., bounded in Ω2x error term, which, for instance, can be written in Lagrange form. In the former case, Ω2x≡ℝn and in the latter case, Ω2x≡ℝn if the error term is also globally bounded for ∀x∈ℝn. Thus, condition Ω2x≡ℝn holds for a large set of nonlinear systems emerging in science and engineering applications and we are going to use it in the remainder of this paper.Let us illustrate how to defineLt,x for a simple but representative vector function. Assume that x=x1x2T and a vector function f∈ℝ2 is defined, e g., as follows: f=a1tx13x22a2tx22T. Then, f2≤f1≤a1x13x22+a2x22≤a1x5+a2x2, where we use that xin≤xn,i=1,2,n∈ℕ, and ℕ is a set of positive integers. Note that additional and more complex examples of this kind are provided in Section 5 (see also [13, 14]).Consequently, using (12), we reviewed (9) in the following more tractable form:(13)z˙=ptz+ctLt,z+F∗t,zt0,z0=z0=w−1t0x0.To reduce the notation, we adopt in (13) and throughout this paper that zt,t0,z0=zt,z0.In turn, it was shown in [13] that(14)xt,x0≤Xt,X0≤zt,z0,where xt,x0,Xt,X0, and zt,z0 are solutions of (1), (9), and (13), respectively.Thus, (14) can be used to bound from above the norm of solutions of (1) by matching solutions of the scalar equation (13).It turned out that (13) provides sound estimates of the trapping/stability regions under the assumption that only the bound of f∗t,x is known—a frequent pronouncement in the control literature [3–5]. Nonetheless, such estimates may become more conservative if f∗t,x is defined explicitly. Consequently, we refined this methodology in [14], where (13) was used to estimate the error of successive approximations of solutions to (1) or (2) stemming from the trapping/stability regions of these systems. This modified approach enhanced our boundedness/stability criteria and delivered approximations that increased successively the accuracy of estimation of the boundaries of the corresponding trapping/stability regions.Nonetheless, the methodologies developed in [13, 14] work under the condition that ct<∞,∀t≥t0, which considerably limits the scope of its applications. Nonetheless, it follows from (11) that frequently limt⟶∞ct=∞ even if A is a stable and time-invariant matrix.The current paper lifts this limitation for a practically important class of nonlinear systems with time-varying and nonperiodic coefficients. Its main contribution is in the development of a modified auxiliary equation withct≡1 and pt≡const under some conditions that are frequently met in various applications. This new technique prompts the development of novel criteria of bondedness/stability and estimation of bondedness/stability regions for a wide class of systems that were intractable to our former methodology [13] and voids elaborate simulations of wt, pt, and ct that were required previously.To make this paper more inclusive, we present the definitions of some standard principles of stability theory which are going to be used in the remainder of this paper. Note that the standard definitions of Lyapunov stability and asymptotic stability for time-varying nonlinear systems that are accepted below can be found, e g., in [1–5]. In the remainder of this paper, we will call these properties shortly either stability or asymptotic stability.We begin with a formal definition of Lyapunov exponents, where we adopt the exposition of these quantities made in [1, 5].For (2), the Lyapunov exponents ϕixt,x0=limt0<t⟶∞supt−t0−1lnxt,x0,i≤n, where xt,x0 is a solution to (2). The maximal Lyapunov exponent, ϕ=maxiϕixt,x0, measures the maximal rate of exponential growth/decay of the corresponding solutions as t⟶∞ and bears a pivoting role in stability theory. For a linear system (3), the Lyapunov exponents are defined as ϕl,i=limt0<t⟶∞supt−t0−1lnσiwt,i≤n, where σi are the singular values of the fundamental solution matrix of (3), and for linear systems, maxiϕl,i=ϕ [1]. This lets us to bound the transition matrix of (3) as follows:(15)Wt,t0≤Dϑe−χt−t0,∀t≥t0,D>0,−χ=ϕ+ϑ,where ϑ is a small positive number [6].In order, let us bring the definition of the comparison principle [4] which is frequently used below. Consider a scalar differential equation(16)u˙1=gt,u1,∀t≥t0,u1t,u10∈ℝ,u1t0,u10=u10,where function gt,u1 is continuous in t and locally Lipschitz in u1 for ∀u1∈℘⊂ℝ. Suppose that the solution to this equation is u1t,u10∈℘,∀t≥t0. Next, consider a differential inequality(17)D+u2≤gt,u2,∀t≥t0,u2t0,u20=u20≤u10,where D+u2 denotes the upper right-hand derivative in t of u2t,u20 [4] and u2t,u20∈℘,∀t≥t0. Then, u1t,u10≥u2t,u20,∀t≥t0.Now we present for convenience some conventional definitions of the trapping/stability regions as follows.Definition 1.
A connected and compact set of all initial vectors,ℑ1t0, is called a trapping region of equation (1) if condition x0∈ℑ1t0 implies that xt,x0∈ℑ1t0,∀t>t0.
Clearly, this definition acknowledges thatℑ1 is the invariant set of (1) containing zero.Definition 2.
A connected and open set of all initial vectors,ℑ2t0, that includes zero vector, is called a region of stability of the trivial solution to (2) if condition x0∈ℑ2t0 implies that xt,x0 is stable.Definition 3.
A connected and open set of all initial vectors,ℑ3t0, that includes zero vector, is called a region of asymptotic stability of the trivial solution to (2) if condition x0∈ℑ3t0 implies that limt⟶∞xt,x0=0.


## 1.2. Motivation Example

The first part of this section derives a novel auxiliary equation withct=1 and pt=const for a simple planar system. The second part infers stability and estimates the stability region of the derived auxiliary equation and extends these inferences to the stability assessment of our planar model (2).Let us firstly apply our former methodology [13] to a simple case, where B=diagλ1λ2, λ1,2∈ℝ,λ1>λ2. Obviously, in this case, wt=diageλ1t−t0eλ2t−t0,w−1t=diage−λ1t−t0e−λ2t−t0 which implies that wt=eλ1t−t0,w−1t=eλ2t−t0. Thus, ct=wtw−1t=eλ1−λ2t−t0 and limt⟶∞ct=∞ if λ1≠λ2 which makes our prior estimates of xt,x0, which is based on (13), overconservative for large values of t. Nonetheless, ct≡1 if λ1,2 are complex conjugates and our prior inferences hold in this case.Let us break out our current approach into a sequence of straightforward steps for a simple version of system (1) where x∈ℝ2, Bt=A+G∗t, A∈ℝ2×2,A=const,A≠0, G∗t∈ℝ2×2, and limt⟶∞t−t0−1∫t0tG∗sds=0. Assume also that matrix A is diagonalizable and eigenvalue and eigenvector matrices of A are Λ=diagλ,λ∗, V=vkj,k,j=1,2, where λ=α+iβ,λ∗=α−iβ, f∗=0εx23T, F∗=0F∗2T.(13) can be applied to the stability analysis of such a planar system if (11) implies that ct<∞,∀t≥t0 which, in general, is difficult to warranty in advance. Thus, we outline a different technique recasting the auxiliary equation in the form, where ct≡1 and pt≡const.(1)
Firstly, let us write our planar model of system (1) in eigenbasis of A as follows:(18)y˙=Λy+Gty+V−10−ε∑j=12v2jyj3T+F2t,wherey∈ℂ2, ℂ2 is a two-dimensional space of complex numbers, x=Vy, Gt=V−1G∗tV, and Ft=V−1F∗t.(2)
Next, we rewrite the last equation as follows:(19)y˙=λI+iby+a−λIy+Gty+V−10−ε∑j=12v2jyj3T+Ft,whereI is an identity matrix, b=diagβ−β, a=diagαα, and λ∈ℝ is going to be determined subsequently.In order, let us select a linear subsystem of (19) with an underlined diagonal matrix as(20)y˙=λI+iby.For this subsystem, a diagonal fundamental solution matrix can be written as follows:wt=eλI+ibt−t0=eλt−t0eibt−t0, which implies that wt=eλt−t0eibt−t0=eλt−t0 since eibt−t0=1. Same reasoning yields that w−1t=e−λt−t0 which shows that ct=wtw−1t=1, and due to (10), pt=λ=const.Thus, if we use (18) as the underlying linear system, then the first addition in right side of the auxiliary equation (13), which is derived for (17), is ptz=λz.(3)
Further application of (13) to (17) brings the second term in the right side of our modified auxiliary equation in the following form:(21)a−λIy+Gty+V−10−ε∑j=12v2jyj3T+Ft≤a−λI+Gty+V−10−ε∑j=12v2jyj3T+Ft.Then, the utility of standard inequality yields(22)0−ε∑j=12v2jyj3T≤absε∑j=12v2jyj3≤εabs∑j=12v2jyj3≤ε∑j=12absv2jabsyj3≤ε∑j=12absv2jy3≤ε∑j=12absv2j3y3,whereabsq=a2+b2 if q=a+ib and absq=q if Imq=0.This lets us to write (13) as follows:(23)z˙=λ+a−λIz+Gtz+σz3+Ft,z0=V−1x0,whereσ=εV−1∑j=12absv2j3.(4)
Lastly, we selectλ to maximize the degree of stability of a linear equation z˙=λ+a−λIz which minimizes conservatism of our estimates and reduces to the following condition.(24)minλλ+a−λI=minλλ+α−λ,since a is a diagonal matrix with equal eigenvalues. Resolving (24) yields that λ=α which casts (23) in the following form:(25)z˙=α+Gtz+σz3+Ft,∀t≥t0,zt0,z0=z0=V−1x0.For further references, we write a homogeneous counterpart to (25) as follows:(26)z˙=α+Gtz+σz3,∀t≥t0,z0=V−1x0.Let us recall now that, due to (14), xt,x0≤zt,z0,z0=V−1x0,∀t≥t0, where xt,x0 and zt,z0 are solutions to planar versions of either equation (1) or (2) and equation (25) or (26), respectively.To further abridge the stability analysis of the trivial solution to (26), we are going to develop for this equation its linear, scalar, and, thus, integrable upper bound. Due to the comparison principle [4], solutions to such equations resolve the essential boundedness/stability properties of equation (26). In general, Lipschitz continuity condition (4) can be used to develop such a linear equation. Nonetheless, in our case, we apply a less conservative bound, z3≤z^2z,z∈0,z^,z^>0. Application of the last inequality brings the following linear equation:(27)Z˙=st,z^Z,∀t≥t0,Zt0,z0,z^=z0=V−1x0,where st,z^=α+Gt+σz^2. Note that solutions to (27) can be used to estimate solutions to (26) if z∈0,z^.In turn, to simplify the assessment of stability of (27), we set that supt0≤tst,z^=−α∗+z^2σ+Gs=−s¯z^, where we assumed that α=−α∗, α∗≥0 and Gs=sup∀t≥t0Gt<∞.Subsequently, the solutions to (26) can be bounded as follows: Z=z0exp∫t0tsτ,z^dτ≤z0exp−∫t0ts¯z^dτ≤z0exp−s¯z^t−t0 which implies that (26) is asymptotically stable if(28)s¯z^>0,z^∈0,z^max,where z^max is the maximal value of z^ for which (28) holds.The last assessment can be drawn through the utility of a more general reasoning which is repeatedly applied in this paper. In fact, determine a comparison equation to (27) as follows:(29)Z¯˙=−s¯z^Z¯,∀t≥t0,Z¯t0,z0,z^=z0=V−1x0.Clearly, the right side of (29) bounds from above the right side of (27), whereas the initial values for both equations are equal. Thus, Zt,z0,z^≤Z¯t,z0,z^=z0exp−s¯z^t−t0,∀t≥t0, due to comparison principle.Next, we infer that (28) implies the asymptotic stability of the trivial solution to (26). In fact, due to the comparison principle,(30)zt,z0≤Zt,z0,z^,∀t≥t0,if zt,z0∈0,z^,∀t≥t0, where zt,z0 and Zt,z0,z^ are solutions to (26) and (27).Let us show thatzt,z0∈0,z^,∀t≥t0 if z0<z^ and (28) holds. Clearly, zt0,z0=z0∈0,z^ under our assumption. Assume, on the contrary, that t=t∗>t0 is the smallest value of t for which zt∗,z0=z^. Then, (30) holds for ∀t∈t0,t∗, and thus zt∗,z0≤Zt∗,z0,z^≤Z¯t∗,z0,z^=z0exp−s¯z^t∗−t0<z0<z^. Hence, zt,z0<z^,∀t≥t0 if (28) holds and z0<z^. Consequently, the trivial solutions to (25) and, in turn, to our planar model of equation (2) are asymptotically stable.Furthermore, solving forz^ equation s¯z^=0 trumps estimation of the range of values of z^ which assure asymptotic stability of (29), i.e., values of z^∈0,z^max, where(31)z^max=α∗−Gsσ12.The last formula leads to the estimation of the stability regions for (26) and the corresponding model of equation (2). Really, (28) and, subsequently, (30) hold if z0<z^∈0,z^max. Hence, the last inequality implies that (29) and, in turn, (27) are asymptotically stable under that condition which, in turn, yields that the trivial solution to (26) is asymptotically stable as well, and a solution to (26) approaches zero as t⟶∞ if z0∈0,z^max. In order, this condition ensures asymptotic stability of the trivial solutions to our simplified model of (2) which relates to (26). Hence, due to the second relation in (26), the region of asymptotic stability of the trivial solution to our planar model of (2) includes interior points of the ellipsoid, which is defined as follows: V−1x0=z^max.Section3 generalizes this last formula and derives the boundedness/stability criteria and estimates the trapping/stability regions for multidimensional (1) and (2).


## 2. Modified Auxiliary Equation

This section derives a modified auxiliary equation withct=1 and pt=const for a broad class of nonlinear systems that are frequently found in applications. Solutions to the auxiliary equation bound from above the matching solutions (1) or (2) under some conditions which we specify below in a theorem which encapsulates our inferences and the underlined assumptions:(1)
At this point, we assume that the average ofBt exists and is defined as follows: At0=limt⟶∞sup∀t≥t0t−t0−1∫t0tBsds<∞,∀t0∈Τ. Additionally, we assume that At0∈ℝn×n is nonzero and diagonalizable matrix ∀t0∈Τ except possibly for some isolated values of t0.This last condition is met for any generic set of matricesAt0 which depends upon a scalar parameter, i.e., such a set of matrices that is structurally stable under small perturbations (see, e.g., [41], chapter 6, pp.235–256 for more details).In order, we set that complex conjugate eigenvalues ofAt0 , λkt0=αkt0±iβkt0,i=−1, 1≤k≤n1≤n,n1∈ℕ, and real eigenvalues of At0, λkt0=αkt0, n1<k≤n, where αk∈ℝ and βk∈ℝ>0, where ℝ>0 is a set of positive numbers. Additionally, we presume that αk≥αk+1,k∈1,n−1 and also define a square diagonal matrix, Λ=α+iβ with α=diagα1,α1,...,αn1,αn1,αn1+1,...,αn, β=diagβ1,−β1,...,βn1,−βn1,0,...,0, α,β∈ℝn×n.(2)
Next, we write (1) as follows:(32)x˙=Ax+G∗tx+f∗t,x+F∗t,∀t≥t0,xt0=x0,whereG∗t=Bt−At0 is a zero mean and continuous matrix.Note that the immediate application of (13) to (32), which is based on utility of the fundamental matrix of solutions to equation, x˙=Ax, i.e., wt=expAt, frequently fails since in this case ct=expα1−αnt and limt⟶∞ct=∞ if α1≠αn.To escape this shortcoming, we rewrite (32) into the eigenbasis of A as follows:(33)y˙=Λy+Gty+ft,y+Ft,∀t≥t0,y∈ℂn,yt0,y0=y0=V−1x0,wherey=V−1x, V∈ℂn×n is eigenmatrix of A, ℂn×n is n2-dimensional space of complex numbers, G=V−1G∗V, ft,y=V−1f∗t,Vy, and Ft=V−1F∗. To reduce the notation, we adopted in (33) that yt,t0,y0=yt,y0.(i)Subsequently, we rewrite (33) as(34)y˙=λI+iβy+α−λIy+Gty+ft,y+Ft,y0=V−1x0.whereλ∈ℝ is defined below. Then, we select y˙=λI+iβy as our underlying linear equation with diagonal matrix, λI+iβ=const, and write the fundamental matrix of solutions to this equation as follows:(35)wt=expλI+iβt−t0=eλt−t0expiβt−t0,whereI is identity matrix. Henceforth, wt=eλt−t0, w−1t=e−λt−t0 since expiβt−t0=1. In turn, application of (10) and (11) yields that pwt=λ and cwt=1 if wt is defined by (35).(3)
Next, adopting (35), we write equation (13) for (34) as follows:(36)z˙=λ+α−λIz+Gtz+Lt,z+Ft,∀λ∈ℝ,zt0,z0=z0=V−1x0,wherezt,z0≥yt,y0,∀t≥t0, y∈Ω2y⊂ℂn, Ω2y is a bounded neighborhood of y≡0, and L:=t0,∞×ℝ≥0⟶ℝ≥0 is continuous in t and z and is locally Lipschitz in z function with Lt,0≡0. Note that this function can be developed through application of inequality (12) to function, ft,y=V−1f∗t,Vy. Hence, (12) in this case takes the following form:(37)ft,y≤Lt,y,∀y∈Ω2y,∀t≥t0.For developingLt,y for polynomial vector fields, see examples in Section 1.1 and Section 1.2 and more inclusive examples in Section 5 as well as examples in [13, 14]. As mentioned earlier, inequality (37), for instance, holds if ft,y is a piece-wise polynomial function or can be approximated by such function with bounded Ω2y Lagrange-type error term. In the former case, Ω2y≡ℂn, and in the latter case, Ω2y≡ℂn if the error term is bounded in ℂn.To straightforward further references, we assume in remainder of this paper thatΩ2y≡ℂn.(4)
In the sequel, we defineλ using the following condition: minλλ+α−λI, which maximizes the degree of stability of a scalar linear equation, z˙=λ+α−λIz. Since α−λI is a diagonal matrix, α−λI=maxkαk−λ which yields that(38)minλλ+α−λI=minλλ+maxkαk−λ.Clearly, the last condition yields thatλ=αn and, consequently, λ+α−λI|λ=αn=α1 (see proof in appendix of this paper).This lets us to write (36) as follows:(39)z˙=α1+Gtz+Lt,z+Ft,zt0,z0=z0=V−1x0,(5)
To develop a less conservative version of (39), we set that matrix Dt=ImdiagGt, D∈ℝn×n, β+t=β+Dt, and G−=G−iD and rewrite (34) as follows:(40)y˙=λI+iβ+ty+α−λIy+G−ty+ft,y+Ft,yt0,y0=V−1x0.The fundamental matrix of solution for equation,y˙=λI+iβ+ty, again can be written as follows: w+t=expλI+iβ+tt−t0 since β+ is a diagonal matrix. As prior, this implies that w+t=eλt−t0, w+−1t=e−λt−t0, cw+t=1, and pw+t=α1. Lastly, a less conservative counterpart of (39) can be written as follows:(41)z˙=α1+G−tz+Lt,z+Ft,zt0,z0=z0=V−1x0.In order, we encapsulate our derivations and the underlying assumptions in the following.Theorem 1.
Assume thatBt is a continuous matrix, and the average of Bt, i.e., matrix At0≠0, exists and is diagonalizable for ∀t0∈Τ except possibly some isolated values of t0. Also, assume that function F∗t is continuous and function f∗t,x is continuous in both variables and inequality (37) holds with Ω2y≡ℂn. Lastly, we assume that equations (1), (39), and (41) possess unique solutions for ∀t≥t0, ∀x0∈H⊂ℝn, and ∀z0=V−1x0. Then, the norms of solutions to equation (1) are bounded by matching solutions to the scalar equations (39) or (41) as follows:(42)xt,x0≤Vzt,z0,z0=V−1x0,∀t≥t0,x0∈H.Proof.
The proof of this statement directly follows from the above assumptions and inferences made prior in this section. In fact, our first set of conditions allows to set up equations (32)–(35). The next set of conditions implies inequality (37) with Ω2y≡ℝn and, in turn, equations (36), (39), and (41). The last assumption assures that functions xt,x0 and zt,z0 exist for ∀t≥t0 and the specified initial values.
Note that the above inferences naturally filter out the effect of imaginary components in matrixΛ on evolution of zt,z0. This equates the degree of stability of equation x˙=Ax, which is measured by the maximal real part of eigenvalues of A, and equation z˙=α1z. Furthermore, (41) is currently defined explicitly, whereas the definition of its counterpart given in [13] involves numerical evaluation of wt, pwt, and cwt.
Albeit that the current version of the auxiliary equation is gaining computational efficiency and widening the scope of the relevant applications, its former counterpart can provide sharper estimates in the common application domain. In fact, our former technique [13] simulates wt for equation (3) which encapsulates the effects of matrices A and G∗t on the time histories of the corresponding solutions. In contrast, our current technique treats G∗t as a conservative perturbation.
To straight-up further referencing, we write a homogeneous counterpart of (41) as follows:(43)z˙=α1+G−tz+Lt,z,zt0,z0=z0=V−1x0.
Analysis of boundedness/stability of solutions to (41) or (43) is simplified under a condition that brings their linear upper bounds which are defined through utility of the following inequality:(44)Lt,z≤l2t,z^z,∀t≥t0,∀z∈0,z^,z^>0,where we assumed that function l2:t0,∞×ℝ≥0⟶ℝ≥0 is continuous in t and bounded in z^. For instance, (44) can be developed via the application of Lipschitz continuity condition to a scalar nonnegative function Lt,z. Note that definition of l2t,z^ is exceedingly simplified in this scalar case.
Furthermore, definition ofl2t,z^ is also abridged if Lt,z is a polynomial or piece-wise polynomial function in z which, in former case, can be written as follows: Lt,z=∑i=1Naitzi, where we assume that ait≥0. In this case, l2t,z^=∑i=1Naitz^i−1 since zi≤z^i−1z,∀z∈0,z^,z^>0.
Clearly, in some cases,l2 may merely depend upon t or be a scalar parameter which, however, complies with our further inferences that prompt closed-form boundedness/stability criteria and estimation of boundedness/stability regions for (1) and (2).


## 3. Boundedness/Stability of Nonautonomous Nonlinear Systems via Application of Modified Linear Auxiliary Equation

Application of (44) to (41) yields its linear counterpart that can be written as follows:(45)Z˙=μt,z^Z+Ft,∀t≥t0,Zt0,z0,z^=z0=V−1x0,where μt,z^=α1+G−t+l2t,z^ and Ft=F0V−1ηt with ηt=1.To shorten the notation, we adopted in (45) and frequently, but always, will set below that Zt,t0,z0,z^=Zt,z0,z^.Clearly, under the conditions of Theorem1 and inequality (44), G−t,Ft, and l2t,z^ are continuous functions which implies that μt,z^ is continuous in both variables as well.Thus, under the condition of Theorem1, (45) admits the following solution:(46)Zt,z0,z^=z0Zht,z^+F0ZFt,z^,where function Zt,z0,z^ is continuous in t and z0, Zht,t0,z^=exp∫t0tμτ,z^dτ, θt,τ,z^=exp∫τtμe,z^de, and ZFt,t0,z^=∫t0tθt,τ,z^V−1ηtdτ. As prior, we frequently, but always, are going to use reduce notations, i.e., Zht,z^=Zht,t0,z^ and ZFt,z^=ZFt,t0,z^. Clearly, Zh>0 and θ>0. Thus, the last inequality yields that ZF>0 as well. We are going to use these inequalities in the proof of Lemma 2.


### 3.1. Stability Criteria and Estimation of Stability Region

To simplify the subsequent references, we will write a homogeneous counterpart to (45) as follows:(47)Z˙=μt,z^Z,Zt=t0,z0,z^=z0=V−1x0.Necessary and sufficient conditions for the stability or asymptotic stability of a linear system directly follow from accessing the behavior of its fundamental solution matrix (see, e.g., a recent paper [42], Lemma 1, and additional references therein). For (47), with continuous in t and bounded in z^ function μt,z^, these conditions can be formulated as follows. Equation (47) is stable if and only if(48)∫t0tμτ,z^dτ<∞,∀z^∈0,z^b,∀t0∈Τ,and asymptotically stable if and only if(49)limt⟶∞sup∫t0tμτ,z^dτ=−∞,∀z^∈0,z^B,∀t0∈Τ,where z^b and z^B are some values of z^ for which either (48) or (49) holds.In turn, due to comparison principal, the stability of a linear scalar equation (47) implies stability of the trivial solutions to nonlinear equation (43) and, in turn, to (2) if zt,z0∈0,z^,∀t≥t0 since compliance with this condition enables the linearization of (43) via application of (44). Next lemma provides an upper bound of Zt,z0 which, subsequently, aids estimation of zt,z0 and, in turn, xt,x0.Lemma 1.
Assume that functionμt,z^ is continuous in t and bounded in z^, (47) is stable for some z^∈0,z^max, where z^max is the maximal value of z^ for which (47) is stable that can be infinity, and z0<z^/Zst0,z^, where Zst0,z^=supt≥t0Zht,t0,z^.Then,(50)supt≥t0Zt,z0,z^<z^,∀t≥t0,∀z0<z^Zst0,z^,z^∈0,z^max.Proof.
Clearly,supt≥t0Zht,t0,z^=Zst0,z^<∞ since Zht,t0,z^ is continuous in t and (47) is stable which yields that supt≥t0Zt,t0,z0,z^=Zst0,z^z0 and consequently (50) holds.Theorem 2.
Assume that assumptions of Theorem1 are met, l2t,z^ is continuous in t and bounded in z^ function, and equation (47) is stable for some z^∈0,z^max. Then, the trivial solution to (2) is stable and inequality (42) takes the form(51)xt,x0≤Vzt,z0≤VZt,z0,z^≤Vz0Zst0,z^,∀t≥t0∈Τ,∀z0<z^Zst0,z^,x0|V−1x0<z^Zst0,z^,z^∈0,z^max,where Zt,z0,z^, zt,z0 and xt,x0 are solutions to equations (47), (43), and (2), respectively, and the region of stability of the trivial solution to (2) includes all values of x0 containing the interior of ellipsoid which is defined as follows:(52)V−1x0<supz^∈0,z^maxz^Zst0,z^,t0∈Τ.Proof.
Let us show that under conditions of this theorem,zt,z0<z^,∀t≥t0,z^∈0,z^max if z0<z^/Zs. Pretend on the contrary that t=t∗>t0 is the smallest value of t such that zt∗,z0=z^ under prior conditions. Then, (44) and, thus, equation (47) hold for ∀t∈t0,t∗ which, due to comparison principle [4] and (50), implies that zt∗,z0≤Zt∗,z0,z^<z^. This contradiction shows that zt,z0<z^,∀t≥t0,∀z0<z^/Zs which, in turn, enables linearization of (43) by (44) and brings a linear equation (47) that enables our inferences. In order, the application of comparison principal prompts that zt,z0≤Zt,z0,z^,∀t≥t0,z0<z^/Zs,z^∈0,z^max which implies that the trivial solution to (43) is stable.
The last inequality lets us to write (42) as (51) which shows that the stability of (47) implies the stability of the trivial solutions to (43) and (2) under the conditions of this theorem.
Let us now estimate the regions of stability of the trivial solutions to (2) and (43). The solution to (43) is stable for all values of z0 which are contingent by the following inequality: z0<supz^∈0,z^maxz^/Zst0,z^. In sequel, since z0=V−1x0, region of stability of the trivial solution to (2) includes all values of x0 meeting condition (52).
Let us recall that (52) depends upon t0—a characteristic property of time-dependent dynamic systems.Example 1.
Let us apply inequality (52) for the estimation of the stability basin of the planar, nonlinear, and nonautonomous equation considered in Section 1.2. Homogeneous counterpart of this planar equation corresponds to nonlinear auxiliary equation (26), which, in turn, relates to its linear complement (27). For this last equation, we derived in Section 1.2 that Zh=exp∫t0tsτ,z^dτ≤exp−s¯z^t−t0 which implies that supt≥t0Zht,t0,z^=Zst0,z^≤supt≥t0exp−s¯z^t−t0=1.
Hence, setting in (52) that Zst0,z^=1 prompts more conservative but abridged estimation of the stability region for our planar system which can be written as follows: V−1x0<z^max. The last formula was derived also in Section 1.2 by application of a simplified inference.Theorem 3.
Assume that assumptions of Theorem1 are met, l2t,z^ is continuous in t and bounded in z^ function, and equation (47) is asymptotically stable for some z^∈0,z^max. Then, the trivial solution to (2) is asymptotically stable, inequality (51) holds where Zt,z0,z^, zt,z0, and xt,x0 are solutions to (47), (43), and (2), respectively, and the region of asymptotic stability of the trivial solution to (2) is defined by (52).Proof.
Literally, under these more conservative conditions, Theorem2 and, hence, (51) hold, which affirms the current statement since in this case limt⟶∞Zt,z0,z^=0 and, due to (51), limt⟶∞zt,z0=0 as well, and, in turn, limt⟶∞xt,x0=0 if z0,z^, and x0 are contingent by conditions of inequalities (51) and (52).
Clearly, in this case, (52) defines the region of asymptotic stability of the trivial solution to (2).
Next, we present some sufficient stability conditions appending the above statements which can be readily attested in a virtually closed-form format.Corollary 1.
Assume that conditions of Theorem1 hold, l2t,z^ is continuous in t and bounded in z^ function, and either (48) or (49) holds as well. Then, inequality (51) holds, where Zt,z0,z^, zt,z0, and xt,x0 are solutions to (47), (43), and (2), respectively, and the trivial solution to (2) is either stable or asymptotically stable, and the region of stability or asymptotic stability of (2) is estimated by (52), where in (52), z^max should be exchanged on either z^b or z^B, respectively.Proof.
Really, under the above conditions, (47) is either stable or asymptotically stable, which, due to Theorem 2 and Theorem 3, assures this statementCorollary 2.
Assume that conditions of Theorem1 hold, l2t,z^ is continuous in t and bounded in z^ function, and(53)μt,z^≤−νz^,∀t≥t0,νz^>0,z^∈0,z^v,where z^v is the maximal value of z^ for which the last inequality holds. Then, the trivial solution of (2) is asymptotically stable. Inequalities (51) and (52) hold, where in (52) z^max should be exchanged on z^v.Proof.
In fact, under the above conditions, (47) is asymptotically stable. Thus, this statement holds due to Theorem 3.
In sequel, let us set thatγt0,z^=limt⟶∞supt−t0−1∫t0tmτ,z^dτ, where mt,z^=G−t+l2t,z^ which steers us to the following.Corollary 3.
Assume that conditions of Theorem1 hold, l2t,z^ is continuous in t and bounded in z^ function, and ϕt0,z^=α1+γt0,z^<0,∀z^∈0,z^ϕ,t0∈Τ, where z^ϕ is the maximal value of z^ for which last inequality holds. Then, the trivial solution to (2) is asymptotically stable and (51) and (52) hold, where in (52), z^max should be exchanged on z^ϕ.Proof.
In fact, conditionϕt0,z^<0 yields that the Lyapunov exponent of solutions to (47) is negative, which assures that (47) is asymptotically stable [1]. Hence, due to Theorem 3, (51) and (52) hold, where Zt,z0,z^, zt,z0, and xt,x0 are solutions to (47), (43), and (2), respectively, which attests this statement.Remark.
Let us represent the last stability criterion in an abridged form. Pretend thatG−t=gavt0+g˜t0,t, where gav=limt⟶∞supt−t0−1∫t0tG−ττdτ, g˜=G−t−gav. Also suppose that l2t,z^=l2avt0,z^+l˜2t0,t,z^, l2avt0,z^=limt⟶∞supt−t0−1∫t0tl2t0,τ,z^dτ, l˜2=l2t,z^−l2avt0,z^. Then, it follows from definition of γ that γt0,z^=gavt0+l2avt0,z0 which implies that stability condition of Lemma 2 can be written as follows:(54)α1+gavt0+l2avt0,z^<0,t0∈Τ.
Note that the above statements grant stability criteria of the trivial solutions to (2) through the application of readily verifiable and frequently closed form formulas if ft,x is a piece-wise polynomial in x since in this case functions Lt,z and l2t,z^ can be defined in analytical form. In contrast, the utilities of either (6) or (8) require simulations of wt and subsequent estimation of the corresponding parameters or function from the simulated data. The solution to such inverse problems can be sensitive to the variation of some of their parameters, which, for instance, can include the length of time interval, etc. Consequently, to our knowledge, these problems have not been attended yet for practically feasible systems.


### 3.2. Boundedness Criteria and Estimation of Region of Boundedness

Next, let us develop some criteria of boundedness of solutions to (1) by assessing the relevant behavior of solutions to (45). These criteria frequently resemble the criteria of stability/asymptotic stability of forced solutions to nonhomogeneous systems. To our knowledge, such criteria and especially the estimation of regions of boundedness of solutions to multidimensional, time-varying nonlinear systems are rare in the current literature. Yet, there is some relation between our boundedness criteria and the concept of input-to-state stability developed by E. D. Sontag under more restrictive conditions (see [4, 43, 44]).A solution to (45) is Zt,t0,z0,z^<∞,∀t≥t0,z^∈0,z^1 if (47) is either stable or asymptotically stable and(55)ρt0,z^=supt≥t0ZFt,t0,z^<∞,z^∈0,z^ρ,∀t0∈Τ,where z^ρ is the maximal value of z^ for which (55) holds, and z^1=minz^max,z^ρ.Firstly, let us present an analog to Lemma1 for the nonhomogeneous equation (45) as follows.Lemma 2.
Assume that functionμt,z^ is continuous in t and bounded in z^, (47) is stable, and(56)z0<z^−F0ρt0,z^/Zst0,z^=zFt0,z^,z^∈0,z^1.Then,(57)supt≥t0Zt,z0,z^<Zst0,z^z0+F0ρt0,z^<z^,∀t≥t0,∀z0<zF,where, as prior, Zst0,z^=supt≥t0Zht,t0,z^.Proof.
As is mentioned in Lemma1, supt≥t0Zht,t0,z^=Zst0,z^<∞ under conditions of this statement. Then, supt≥t0Zt,t0,z0,z^≤z0supt≥t0Zht,t0,z^+supt≥t0ZFt,t0,z^=z0Zst0,z^+F0ρt0,z^ since Zht0,t,z^≥0 and ZFt0,t,z^≥0 which yields (56)
The next statement assesses the boundedness of solutions to the nonhomogeneous equation (41) and, in turn, (1).Theorem 4.
Assume that the assumptions of Theorem1 hold, equation (47) is stable for some z^∈0,z^max, inequalities (55) and (56) hold, and l2t,z^ is continuous in t and bounded in z^ function. Then,(i)
xt,x0≤∞,∀t≥t0 if these solutions are emanated from the interior of the ellipsoid that is defined as follows:(58)x0|V−1x0<supz^∈0,z^1zFt0,z^.(ii)
Inequality (42) takes the form(59)xt,x0≤Vzt,z0≤VZt,z0,z^=Vz0Zht,z0+F0ZFt,z0≤Vz0Zst0,z^+F0ρt0,z^;∀t≥t0,x0|V−1x0<zF,z0<zF,z^∈0,z^1,∀t0∈Τ,whereZt,t0,z0,z^, zt,z0, and xt,t0,x0 are solutions to (45), (41), and (1), respectively.(iii)
If (47) is asymptotically stable and other conditions of this statement hold, then(60)limt⟶∞supxt,x0≤F0Vρt0,z^,x0|V−1x0<zF,z^∈0,z^1,∀t0∈Τ.Proof.
Firstly, we show that under the conditions of this theorem,zt,z0<z^,∀t≥t0,z^∈0,z^1 if z0<zF. In fact, pretend, on the contrary, that t=t∗>t0 is the smallest value of t such that zt∗,z0=z^ under prior conditions. Then, (44) and, thus, (45) hold for ∀t∈t0,t∗ which, due to comparison principle and (56), implies that zt∗,z0≤Zt∗,z0,z^<z^. This contradiction prompts that zt,z0<z^,∀t≥t0,∀z0<zF, z^∈0,z^1 which, in turn, enables a linearization of (41) with the utility of (44) that, subsequently, yields equation (45). In order, the application of comparison principal returns that zt,z0≤Zt,z0,z^<∞,∀t≥t0,z0<zF,z^∈0,z^1 since (47) is stable and (56) holds.
Clearly, the last inequality lets us to write (42) as (59), which demonstrates that the solutions of (1) are bounded within the region specified by (58).
Lastly, (60) directly follows from (59) if (47) is asymptotically stable since in this case limt⟶∞Zht,x0=0.
Next, two corollaries present simplified sufficient conditions under which the presumptions of Theorem4 are met.Corollary 5.
Assume that the conditions of Theorem1 and inequalities (53), (55), and (56) hold, and l2t,z^ is continuous in t and bounded in z^ function. Then, the solutions to (1) are bounded in the region that is defined by (58), inequality (59) holds, and(61)limt→∞supxt,x0≤F0VV−1v,∀x0|V−1x0=zF.Note that this statement assumes that in (58) and (59), z^1=minz^v,z^ρ.Proof.
In fact, the condition of this statement implies the asymptotic stability of (47) for ∀z^∈0,z^v which, together with other assumptions, yields that the conditions of Theorem 4 hold. In turn, this implies that (59) holds which yields that limt⟶∞supxt,x0≤F0Vsupt≥t0ZFt,x0. Next, ZFt,z^=∫t0tθt,τ,z^V−1ηtdτ≤∫t0tθt,τ,z^V−1ηtdτ≤V−1∫t0te−vt−τdτ≤V−1/v.
In the sequel, let us consider the Lyapunov exponent of solutions to (47), i.e., ϕz^=α1+γz^<0,z^∈0,z^ϕ,z^ϕ>0 (see Corollary 3 for additional details). This implies that Zht,z0,z^≤D`ϑ,z^e−χt−t0 and, due to (15), that θt,τ,z^≤Dϑ,z^e−χt−τ, where −χ=ϕ+ϑ , χ>0, ϑ>0 is a small value, and constants Dϑ,z^>0 [6]. This leads to the following.Corollary 6.
Assume that the conditions of Theorem1 and inequalities (55) and (47) hold, ϕt0,z^<0,∀z^∈0,z^ϕ, and function l2t,z^ is continuous in t and bounded in z^ function. Then, (59) and (60) hold and limt⟶∞xt,x0≤F0VV−1Dϑ,z^/χ,∀x0|V−1x0=zF.Note that in this statement, it is assumed that in (58) and 59, z^1=minz^ϕ,z^ρ.Proof.
In fact, in this case, (47) is asymptotically stable [1] and the conditions of Theorem 4 hold. Thus, (59) yields that limt⟶∞supxt,x0≤F0Vsupt≥t0ZFt,x0. In turn, ZFt,t0,z^=∫t0tθt,τ,z^V−1ηtdτ≤∫t0tθt,τ,z^V−1ηtdτ≤V−1D`ϑ,z^∫t0te−χt−τdτ≤V−1Dϑ,z^/χ.Example 2.
To abridge our derivation, we apply this methodology to estimation of the region of boundedness of solutions to (25), which, in turn, prompts estimation of the corresponding region for the planar system considered in Section 1.2. Yet, (25) is a scalar but nonintegrable equation and direct estimation of its boundedness region requires repeated simulations of this equation. Thus, to simplify our further inferences, we bring an autonomous counterpart to (25) which can be written as follows:(62)z˙1=−α∗+Gsz1+σz13+F˜,z1t0,z0=z0=V−1x0,where, as prior, we assume that Gs=sup∀t≥t0G−t<∞ and Ft≤F0V−1=F˜=const. Thus, due to comparison principle, zt,z0≤z1t,z0,∀t≥t0, where zt,z0 is a solution to (25). In turn, the boundedness of solutions to the autonomous and integrable equation (62) is determined by the location of its readily accessible fixed solutions.
As in Section1.2, to simplify the assessment of the behavior of solutions to (62), we consider its linear counterpart(63)Z˙1=−s¯z^Z1+F˜,∀t≥t0,Z1t0,z0,z^=z0=V−1x0,where −s¯z^=−α∗+Gs+z^2σ. Thus, (63) can be viewed as a nonhomogeneous analog to (29).
Furthermore, let us write a nonhomogeneous counterpart to (27) as follows:(64)Z˙2=st,z^Z2+F˜,∀t≥t0,z0<z^,Z2t0,z0,z^=z0=V−1x0.
Thus, due to comparison principle,z1t,z0≤Z1t,z0,z^,∀t≥t0, and Z2t,z0,z^≤Z1t,z0,z^, ∀t≥t0.
Clearly, the general solutions to the linear scalar (63) and (64) take the form Zit,z0,z^=z0Zi,ht,z0,z^+F˜Zi,Ft,z0,z^,i=1,2,∀t≥t0, where Zi,ht,z0,z^,i=1,2, are solutions to (29) and (27) with z0=1 and Zi,Ft,z0,z^ are components of forced solution to the corresponding equations. Let us recall that sup∀t≥t0Z1,ht,z0,z^=Z1,sz0,z^=sup∀t≥t0exp−s¯z^t−t0=1 (see Section 1.2).
To estimate the right side of (56), we recall that Zsz0,z^=sup∀t≥t0Z2,ht,z0,z^≤sup∀t≥t0Z1,ht,z0,z^=1.
In sequence,θt,τ,z^=exp∫τtse,z^de≤exp∫τt−s¯z^de≤exp−s¯z^t−τ which, with use of (55), lets us bound ρt0,z^ as follows:(65)ρt0,z^=sup∀t≥t0∫t0tθt,τ,z^V−1ηtdτ≤V−1sup∀t≥t0∫t0tθt,τ,z^dτ≤V−1sup∀t≥t0∫t0texp−s¯z^τ−t0dτ≤V−1sz^sup∀t≥t01−exp−sz^t−t0=V−1s¯z^.
Then, using the developed inequalities, we bound the right side of (58) as follows:(66)supz^∈0,z^1zFz^=supz^∈0,z^1z^−F0ρt0,z^Zst0,z^≤supz^∈0,z^1z^−F˜s¯z^,where in this case, z^1=z^∗ and z^∗ is the smallest positive root of equation s¯z^=0.
Next, we infer thatsupz^∈0,z^1z^−F˜/s¯z^ is reached for z^=z^s<z∗ since limz^⟶z∗1/s¯z^=∞. Note that zs can be readily found numerically. Hence, the region of boundedness of solutions to (63) and, thus, to (62) and, subsequently, (25) can be estimated as follows:(67)z0<z^s−F˜s¯z^s=R.
At this point, we consider a numerical example withα∗=0.3 and Gs=σ=F˜=0.1. For this set of parameters, (62) has three real fixed solutions in which two solutions are positive and equal to 0.618 and one. Furthermore, in this case, the smallest positive fixed solution is stable, whereas the largest is unstable. Thus, the solutions to (62) are bounded and approach the stable fixed solution as t⟶∞ if z0∈0,1 and approach infinity as t⟶∞ if z0>1. Yet, for this set of parameters, z∗=2,zs=0.848, and R=0.07. The last value provides a rather conservative estimate of the boundedness interval of solutions to (62) since our current derivation rests on a sequence of upper bounds which simplifies our inference. More accurate estimates can be obtained through direct simulation of (58). Nonetheless, the above reasoning authenticates our theory.


## 3.1. Stability Criteria and Estimation of Stability Region

To simplify the subsequent references, we will write a homogeneous counterpart to (45) as follows:(47)Z˙=μt,z^Z,Zt=t0,z0,z^=z0=V−1x0.Necessary and sufficient conditions for the stability or asymptotic stability of a linear system directly follow from accessing the behavior of its fundamental solution matrix (see, e.g., a recent paper [42], Lemma 1, and additional references therein). For (47), with continuous in t and bounded in z^ function μt,z^, these conditions can be formulated as follows. Equation (47) is stable if and only if(48)∫t0tμτ,z^dτ<∞,∀z^∈0,z^b,∀t0∈Τ,and asymptotically stable if and only if(49)limt⟶∞sup∫t0tμτ,z^dτ=−∞,∀z^∈0,z^B,∀t0∈Τ,where z^b and z^B are some values of z^ for which either (48) or (49) holds.In turn, due to comparison principal, the stability of a linear scalar equation (47) implies stability of the trivial solutions to nonlinear equation (43) and, in turn, to (2) if zt,z0∈0,z^,∀t≥t0 since compliance with this condition enables the linearization of (43) via application of (44). Next lemma provides an upper bound of Zt,z0 which, subsequently, aids estimation of zt,z0 and, in turn, xt,x0.Lemma 1.
Assume that functionμt,z^ is continuous in t and bounded in z^, (47) is stable for some z^∈0,z^max, where z^max is the maximal value of z^ for which (47) is stable that can be infinity, and z0<z^/Zst0,z^, where Zst0,z^=supt≥t0Zht,t0,z^.Then,(50)supt≥t0Zt,z0,z^<z^,∀t≥t0,∀z0<z^Zst0,z^,z^∈0,z^max.Proof.
Clearly,supt≥t0Zht,t0,z^=Zst0,z^<∞ since Zht,t0,z^ is continuous in t and (47) is stable which yields that supt≥t0Zt,t0,z0,z^=Zst0,z^z0 and consequently (50) holds.Theorem 2.
Assume that assumptions of Theorem1 are met, l2t,z^ is continuous in t and bounded in z^ function, and equation (47) is stable for some z^∈0,z^max. Then, the trivial solution to (2) is stable and inequality (42) takes the form(51)xt,x0≤Vzt,z0≤VZt,z0,z^≤Vz0Zst0,z^,∀t≥t0∈Τ,∀z0<z^Zst0,z^,x0|V−1x0<z^Zst0,z^,z^∈0,z^max,where Zt,z0,z^, zt,z0 and xt,x0 are solutions to equations (47), (43), and (2), respectively, and the region of stability of the trivial solution to (2) includes all values of x0 containing the interior of ellipsoid which is defined as follows:(52)V−1x0<supz^∈0,z^maxz^Zst0,z^,t0∈Τ.Proof.
Let us show that under conditions of this theorem,zt,z0<z^,∀t≥t0,z^∈0,z^max if z0<z^/Zs. Pretend on the contrary that t=t∗>t0 is the smallest value of t such that zt∗,z0=z^ under prior conditions. Then, (44) and, thus, equation (47) hold for ∀t∈t0,t∗ which, due to comparison principle [4] and (50), implies that zt∗,z0≤Zt∗,z0,z^<z^. This contradiction shows that zt,z0<z^,∀t≥t0,∀z0<z^/Zs which, in turn, enables linearization of (43) by (44) and brings a linear equation (47) that enables our inferences. In order, the application of comparison principal prompts that zt,z0≤Zt,z0,z^,∀t≥t0,z0<z^/Zs,z^∈0,z^max which implies that the trivial solution to (43) is stable.
The last inequality lets us to write (42) as (51) which shows that the stability of (47) implies the stability of the trivial solutions to (43) and (2) under the conditions of this theorem.
Let us now estimate the regions of stability of the trivial solutions to (2) and (43). The solution to (43) is stable for all values of z0 which are contingent by the following inequality: z0<supz^∈0,z^maxz^/Zst0,z^. In sequel, since z0=V−1x0, region of stability of the trivial solution to (2) includes all values of x0 meeting condition (52).
Let us recall that (52) depends upon t0—a characteristic property of time-dependent dynamic systems.Example 1.
Let us apply inequality (52) for the estimation of the stability basin of the planar, nonlinear, and nonautonomous equation considered in Section 1.2. Homogeneous counterpart of this planar equation corresponds to nonlinear auxiliary equation (26), which, in turn, relates to its linear complement (27). For this last equation, we derived in Section 1.2 that Zh=exp∫t0tsτ,z^dτ≤exp−s¯z^t−t0 which implies that supt≥t0Zht,t0,z^=Zst0,z^≤supt≥t0exp−s¯z^t−t0=1.
Hence, setting in (52) that Zst0,z^=1 prompts more conservative but abridged estimation of the stability region for our planar system which can be written as follows: V−1x0<z^max. The last formula was derived also in Section 1.2 by application of a simplified inference.Theorem 3.
Assume that assumptions of Theorem1 are met, l2t,z^ is continuous in t and bounded in z^ function, and equation (47) is asymptotically stable for some z^∈0,z^max. Then, the trivial solution to (2) is asymptotically stable, inequality (51) holds where Zt,z0,z^, zt,z0, and xt,x0 are solutions to (47), (43), and (2), respectively, and the region of asymptotic stability of the trivial solution to (2) is defined by (52).Proof.
Literally, under these more conservative conditions, Theorem2 and, hence, (51) hold, which affirms the current statement since in this case limt⟶∞Zt,z0,z^=0 and, due to (51), limt⟶∞zt,z0=0 as well, and, in turn, limt⟶∞xt,x0=0 if z0,z^, and x0 are contingent by conditions of inequalities (51) and (52).
Clearly, in this case, (52) defines the region of asymptotic stability of the trivial solution to (2).
Next, we present some sufficient stability conditions appending the above statements which can be readily attested in a virtually closed-form format.Corollary 1.
Assume that conditions of Theorem1 hold, l2t,z^ is continuous in t and bounded in z^ function, and either (48) or (49) holds as well. Then, inequality (51) holds, where Zt,z0,z^, zt,z0, and xt,x0 are solutions to (47), (43), and (2), respectively, and the trivial solution to (2) is either stable or asymptotically stable, and the region of stability or asymptotic stability of (2) is estimated by (52), where in (52), z^max should be exchanged on either z^b or z^B, respectively.Proof.
Really, under the above conditions, (47) is either stable or asymptotically stable, which, due to Theorem 2 and Theorem 3, assures this statementCorollary 2.
Assume that conditions of Theorem1 hold, l2t,z^ is continuous in t and bounded in z^ function, and(53)μt,z^≤−νz^,∀t≥t0,νz^>0,z^∈0,z^v,where z^v is the maximal value of z^ for which the last inequality holds. Then, the trivial solution of (2) is asymptotically stable. Inequalities (51) and (52) hold, where in (52) z^max should be exchanged on z^v.Proof.
In fact, under the above conditions, (47) is asymptotically stable. Thus, this statement holds due to Theorem 3.
In sequel, let us set thatγt0,z^=limt⟶∞supt−t0−1∫t0tmτ,z^dτ, where mt,z^=G−t+l2t,z^ which steers us to the following.Corollary 3.
Assume that conditions of Theorem1 hold, l2t,z^ is continuous in t and bounded in z^ function, and ϕt0,z^=α1+γt0,z^<0,∀z^∈0,z^ϕ,t0∈Τ, where z^ϕ is the maximal value of z^ for which last inequality holds. Then, the trivial solution to (2) is asymptotically stable and (51) and (52) hold, where in (52), z^max should be exchanged on z^ϕ.Proof.
In fact, conditionϕt0,z^<0 yields that the Lyapunov exponent of solutions to (47) is negative, which assures that (47) is asymptotically stable [1]. Hence, due to Theorem 3, (51) and (52) hold, where Zt,z0,z^, zt,z0, and xt,x0 are solutions to (47), (43), and (2), respectively, which attests this statement.Remark.
Let us represent the last stability criterion in an abridged form. Pretend thatG−t=gavt0+g˜t0,t, where gav=limt⟶∞supt−t0−1∫t0tG−ττdτ, g˜=G−t−gav. Also suppose that l2t,z^=l2avt0,z^+l˜2t0,t,z^, l2avt0,z^=limt⟶∞supt−t0−1∫t0tl2t0,τ,z^dτ, l˜2=l2t,z^−l2avt0,z^. Then, it follows from definition of γ that γt0,z^=gavt0+l2avt0,z0 which implies that stability condition of Lemma 2 can be written as follows:(54)α1+gavt0+l2avt0,z^<0,t0∈Τ.
Note that the above statements grant stability criteria of the trivial solutions to (2) through the application of readily verifiable and frequently closed form formulas if ft,x is a piece-wise polynomial in x since in this case functions Lt,z and l2t,z^ can be defined in analytical form. In contrast, the utilities of either (6) or (8) require simulations of wt and subsequent estimation of the corresponding parameters or function from the simulated data. The solution to such inverse problems can be sensitive to the variation of some of their parameters, which, for instance, can include the length of time interval, etc. Consequently, to our knowledge, these problems have not been attended yet for practically feasible systems.


## 3.2. Boundedness Criteria and Estimation of Region of Boundedness

Next, let us develop some criteria of boundedness of solutions to (1) by assessing the relevant behavior of solutions to (45). These criteria frequently resemble the criteria of stability/asymptotic stability of forced solutions to nonhomogeneous systems. To our knowledge, such criteria and especially the estimation of regions of boundedness of solutions to multidimensional, time-varying nonlinear systems are rare in the current literature. Yet, there is some relation between our boundedness criteria and the concept of input-to-state stability developed by E. D. Sontag under more restrictive conditions (see [4, 43, 44]).A solution to (45) is Zt,t0,z0,z^<∞,∀t≥t0,z^∈0,z^1 if (47) is either stable or asymptotically stable and(55)ρt0,z^=supt≥t0ZFt,t0,z^<∞,z^∈0,z^ρ,∀t0∈Τ,where z^ρ is the maximal value of z^ for which (55) holds, and z^1=minz^max,z^ρ.Firstly, let us present an analog to Lemma1 for the nonhomogeneous equation (45) as follows.Lemma 2.
Assume that functionμt,z^ is continuous in t and bounded in z^, (47) is stable, and(56)z0<z^−F0ρt0,z^/Zst0,z^=zFt0,z^,z^∈0,z^1.Then,(57)supt≥t0Zt,z0,z^<Zst0,z^z0+F0ρt0,z^<z^,∀t≥t0,∀z0<zF,where, as prior, Zst0,z^=supt≥t0Zht,t0,z^.Proof.
As is mentioned in Lemma1, supt≥t0Zht,t0,z^=Zst0,z^<∞ under conditions of this statement. Then, supt≥t0Zt,t0,z0,z^≤z0supt≥t0Zht,t0,z^+supt≥t0ZFt,t0,z^=z0Zst0,z^+F0ρt0,z^ since Zht0,t,z^≥0 and ZFt0,t,z^≥0 which yields (56)
The next statement assesses the boundedness of solutions to the nonhomogeneous equation (41) and, in turn, (1).Theorem 4.
Assume that the assumptions of Theorem1 hold, equation (47) is stable for some z^∈0,z^max, inequalities (55) and (56) hold, and l2t,z^ is continuous in t and bounded in z^ function. Then,(i)
xt,x0≤∞,∀t≥t0 if these solutions are emanated from the interior of the ellipsoid that is defined as follows:(58)x0|V−1x0<supz^∈0,z^1zFt0,z^.(ii)
Inequality (42) takes the form(59)xt,x0≤Vzt,z0≤VZt,z0,z^=Vz0Zht,z0+F0ZFt,z0≤Vz0Zst0,z^+F0ρt0,z^;∀t≥t0,x0|V−1x0<zF,z0<zF,z^∈0,z^1,∀t0∈Τ,whereZt,t0,z0,z^, zt,z0, and xt,t0,x0 are solutions to (45), (41), and (1), respectively.(iii)
If (47) is asymptotically stable and other conditions of this statement hold, then(60)limt⟶∞supxt,x0≤F0Vρt0,z^,x0|V−1x0<zF,z^∈0,z^1,∀t0∈Τ.Proof.
Firstly, we show that under the conditions of this theorem,zt,z0<z^,∀t≥t0,z^∈0,z^1 if z0<zF. In fact, pretend, on the contrary, that t=t∗>t0 is the smallest value of t such that zt∗,z0=z^ under prior conditions. Then, (44) and, thus, (45) hold for ∀t∈t0,t∗ which, due to comparison principle and (56), implies that zt∗,z0≤Zt∗,z0,z^<z^. This contradiction prompts that zt,z0<z^,∀t≥t0,∀z0<zF, z^∈0,z^1 which, in turn, enables a linearization of (41) with the utility of (44) that, subsequently, yields equation (45). In order, the application of comparison principal returns that zt,z0≤Zt,z0,z^<∞,∀t≥t0,z0<zF,z^∈0,z^1 since (47) is stable and (56) holds.
Clearly, the last inequality lets us to write (42) as (59), which demonstrates that the solutions of (1) are bounded within the region specified by (58).
Lastly, (60) directly follows from (59) if (47) is asymptotically stable since in this case limt⟶∞Zht,x0=0.
Next, two corollaries present simplified sufficient conditions under which the presumptions of Theorem4 are met.Corollary 5.
Assume that the conditions of Theorem1 and inequalities (53), (55), and (56) hold, and l2t,z^ is continuous in t and bounded in z^ function. Then, the solutions to (1) are bounded in the region that is defined by (58), inequality (59) holds, and(61)limt→∞supxt,x0≤F0VV−1v,∀x0|V−1x0=zF.Note that this statement assumes that in (58) and (59), z^1=minz^v,z^ρ.Proof.
In fact, the condition of this statement implies the asymptotic stability of (47) for ∀z^∈0,z^v which, together with other assumptions, yields that the conditions of Theorem 4 hold. In turn, this implies that (59) holds which yields that limt⟶∞supxt,x0≤F0Vsupt≥t0ZFt,x0. Next, ZFt,z^=∫t0tθt,τ,z^V−1ηtdτ≤∫t0tθt,τ,z^V−1ηtdτ≤V−1∫t0te−vt−τdτ≤V−1/v.
In the sequel, let us consider the Lyapunov exponent of solutions to (47), i.e., ϕz^=α1+γz^<0,z^∈0,z^ϕ,z^ϕ>0 (see Corollary 3 for additional details). This implies that Zht,z0,z^≤D`ϑ,z^e−χt−t0 and, due to (15), that θt,τ,z^≤Dϑ,z^e−χt−τ, where −χ=ϕ+ϑ , χ>0, ϑ>0 is a small value, and constants Dϑ,z^>0 [6]. This leads to the following.Corollary 6.
Assume that the conditions of Theorem1 and inequalities (55) and (47) hold, ϕt0,z^<0,∀z^∈0,z^ϕ, and function l2t,z^ is continuous in t and bounded in z^ function. Then, (59) and (60) hold and limt⟶∞xt,x0≤F0VV−1Dϑ,z^/χ,∀x0|V−1x0=zF.Note that in this statement, it is assumed that in (58) and 59, z^1=minz^ϕ,z^ρ.Proof.
In fact, in this case, (47) is asymptotically stable [1] and the conditions of Theorem 4 hold. Thus, (59) yields that limt⟶∞supxt,x0≤F0Vsupt≥t0ZFt,x0. In turn, ZFt,t0,z^=∫t0tθt,τ,z^V−1ηtdτ≤∫t0tθt,τ,z^V−1ηtdτ≤V−1D`ϑ,z^∫t0te−χt−τdτ≤V−1Dϑ,z^/χ.Example 2.
To abridge our derivation, we apply this methodology to estimation of the region of boundedness of solutions to (25), which, in turn, prompts estimation of the corresponding region for the planar system considered in Section 1.2. Yet, (25) is a scalar but nonintegrable equation and direct estimation of its boundedness region requires repeated simulations of this equation. Thus, to simplify our further inferences, we bring an autonomous counterpart to (25) which can be written as follows:(62)z˙1=−α∗+Gsz1+σz13+F˜,z1t0,z0=z0=V−1x0,where, as prior, we assume that Gs=sup∀t≥t0G−t<∞ and Ft≤F0V−1=F˜=const. Thus, due to comparison principle, zt,z0≤z1t,z0,∀t≥t0, where zt,z0 is a solution to (25). In turn, the boundedness of solutions to the autonomous and integrable equation (62) is determined by the location of its readily accessible fixed solutions.
As in Section1.2, to simplify the assessment of the behavior of solutions to (62), we consider its linear counterpart(63)Z˙1=−s¯z^Z1+F˜,∀t≥t0,Z1t0,z0,z^=z0=V−1x0,where −s¯z^=−α∗+Gs+z^2σ. Thus, (63) can be viewed as a nonhomogeneous analog to (29).
Furthermore, let us write a nonhomogeneous counterpart to (27) as follows:(64)Z˙2=st,z^Z2+F˜,∀t≥t0,z0<z^,Z2t0,z0,z^=z0=V−1x0.
Thus, due to comparison principle,z1t,z0≤Z1t,z0,z^,∀t≥t0, and Z2t,z0,z^≤Z1t,z0,z^, ∀t≥t0.
Clearly, the general solutions to the linear scalar (63) and (64) take the form Zit,z0,z^=z0Zi,ht,z0,z^+F˜Zi,Ft,z0,z^,i=1,2,∀t≥t0, where Zi,ht,z0,z^,i=1,2, are solutions to (29) and (27) with z0=1 and Zi,Ft,z0,z^ are components of forced solution to the corresponding equations. Let us recall that sup∀t≥t0Z1,ht,z0,z^=Z1,sz0,z^=sup∀t≥t0exp−s¯z^t−t0=1 (see Section 1.2).
To estimate the right side of (56), we recall that Zsz0,z^=sup∀t≥t0Z2,ht,z0,z^≤sup∀t≥t0Z1,ht,z0,z^=1.
In sequence,θt,τ,z^=exp∫τtse,z^de≤exp∫τt−s¯z^de≤exp−s¯z^t−τ which, with use of (55), lets us bound ρt0,z^ as follows:(65)ρt0,z^=sup∀t≥t0∫t0tθt,τ,z^V−1ηtdτ≤V−1sup∀t≥t0∫t0tθt,τ,z^dτ≤V−1sup∀t≥t0∫t0texp−s¯z^τ−t0dτ≤V−1sz^sup∀t≥t01−exp−sz^t−t0=V−1s¯z^.
Then, using the developed inequalities, we bound the right side of (58) as follows:(66)supz^∈0,z^1zFz^=supz^∈0,z^1z^−F0ρt0,z^Zst0,z^≤supz^∈0,z^1z^−F˜s¯z^,where in this case, z^1=z^∗ and z^∗ is the smallest positive root of equation s¯z^=0.
Next, we infer thatsupz^∈0,z^1z^−F˜/s¯z^ is reached for z^=z^s<z∗ since limz^⟶z∗1/s¯z^=∞. Note that zs can be readily found numerically. Hence, the region of boundedness of solutions to (63) and, thus, to (62) and, subsequently, (25) can be estimated as follows:(67)z0<z^s−F˜s¯z^s=R.
At this point, we consider a numerical example withα∗=0.3 and Gs=σ=F˜=0.1. For this set of parameters, (62) has three real fixed solutions in which two solutions are positive and equal to 0.618 and one. Furthermore, in this case, the smallest positive fixed solution is stable, whereas the largest is unstable. Thus, the solutions to (62) are bounded and approach the stable fixed solution as t⟶∞ if z0∈0,1 and approach infinity as t⟶∞ if z0>1. Yet, for this set of parameters, z∗=2,zs=0.848, and R=0.07. The last value provides a rather conservative estimate of the boundedness interval of solutions to (62) since our current derivation rests on a sequence of upper bounds which simplifies our inference. More accurate estimates can be obtained through direct simulation of (58). Nonetheless, the above reasoning authenticates our theory.


## 4. Modified Nonlinear Auxiliary Equation and Boundedness/Stability of Nonautonomous and Nonlinear Systems

Analysis of solutions to the scalar nonlinear auxiliary (41) and (43) conveys the less conservative criteria of boundedness/stability of the initial nonlinear and multidimensional systems and naturally enfolds the estimation of their trapping/stability regions.In general, (41) is a nonintegrable, scalar, and nonlinear equation with time-varying coefficients and the behavior of its solutions can be readily assessed either in simulations or by laying out some integrable nonlinear and autonomous equations bounding (41) from above and below. The former approach is substantially abridged due to the following.Theorem 5.
Assume that equation (41) possesses a unique solution zt,z0,∀t≥t0,z0≥0, and zt,z0′ and zt,z0′′ are solutions to (41) with 0≤z0′≤z0′′. Then, zt,z0′≤zt,z0′′,∀t≥t0.Proof.
In fact, the solutions to (41) do not intersect in t×z-plane due to uniqueness property of this equation.
Application of the above statement grants some boundedness/stability criteria for (1) and (2), respectively. Firstly, let us define a set centered at zero concentric ellipsoids, Ez⊂ℝn, as follows:(68)Ez:V−1x=z≥0.
Also, we assume that∂Ez⊂ℝn−1 defines the boundaries of these ellipsoids and E−z=Ez−∂Ez. This prompts the following.Theorem 6.
Assume that equations (2) and (43) possess unique solutions for ∀x0∈H⊂ℝn and ∀z0=V−1x0, and that the trivial solution to (43) is either stable or asymptotically stable. Then, the trivial solution to (2) is stable/asymptotically stable as well, respectively. Furthermore, assume that the interval 0,z¯ defines the stability region of the stable/asymptotically stable trivial solution to (24). Then, set E−z¯ is contained in stability region of (2).Proof.
Really, the proof of this statement follows from inequality (42), where it is presumed that xt,x0 and zt,z0 are solutions to equations (2) and (43), respectively.Theorem 7.
Assume that equations (2) and (41) possess unique solutions for ∀x0∈H⊂ℝn and ∀z0=V−1x0, and interval 0,z¯ defines a trapping region of (41) about z≡0, i.e., zt,z0≤z¯, ∀z0∈0,z¯,∀t≥t0. Then, the ellipsoid Ez¯ is included into the trapping region of equation (1) about zero, i.e., xt,x0<∞,∀x0∈Ez¯,∀t≥t0.Proof.
The proof of this statement follows from inequality (42), where xt,x0 and zt,z0 are solutions to equations (1) and (41), respectively.
Thus, the problem of estimating the trapping/stability regions about pointx≡0 of multidimensional equation (1) or (2) is reduced to reckoning the threshold value z¯ splitting the interval of initial values of solutions to (41) or (43) in two parts associated with qualitatively distinct behavior of such solutions on long time intervals. Then, the boundary of the trapping/stability region for either (1) or (2) is given by the following formula: z¯=V−1x0. In turn, the task of simulating the threshold value is markedly simplified since zt,z0 monotonically increases in z0,∀t≥t0 due to Theorem 5.
Yet, the structure of solutions to (41) or (43) and, in turn, to their multidimensional counterparts (1) and (2) can be further divulged through the analysis of solutions to some simplified equations which bound (41) or (43). For instance, we acknowledge again that under (44) the solutions to linear equations (45) and (47) bound from above the corresponding solutions to equations (41) and (43) that are stemmed from the same initial values.
A distinct stability criterion can be recouped using the Massera, Chetaev, and Malkin theorem, which replaces (44) by a complementary condition (see [1] for contemporary references and citations). In this regard, we assume that Lt,z=L1tz+Lrt,z, where L1:t0,∞⟶ℝ≥0 is continuous in t and Lr:t0,∞×ℝ≥0⟶ℝ≥0 is continuous in both variables and also Lipschitz in z, and(69)Lrt,z≤l3zr,∀t≥t0,r>1,∀z∈0,z˜,l3,z˜>0.
Then, (43) can be written as(70)z˙=ψtz+Lrt,z,zt0=z0=V−1x0,where ψt=α1+G−t+L1t.
Next, letςt,τ=∫τtψsds, where(71)ςt,τ≤ς0−ς1t−τ+ς2τ,∀t≥τ≥0,ς0,ς1>0,ς2≥0,and(72)r−1ς1≥ς2.
This comprises the following.Theorem 8.
Assume that equations (2) and (43) possess unique solutions for ∀x0∈H⊂ℝn and ∀z0≥0, and that (69) is contended with sufficiently small l3 and (71) and (72) hold as well. Then, the trivial solution to (70) is asymptotically stable, which, in turn, due to the comparison principle, implies that the trivial solution to the corresponding equation (2) is asymptotically stable as well.Proof.
Indeed, the proof of this statement directly follows from application of the above-mentioned theorem to equation (70) which represents equation (43) in this case. Then, this statement is assured due to application of the comparison principle
Some additional approaches to analysis of the structure of solutions to a scalar auxiliary equation were outlined in [13]. Below we complement and apply some of these techniques to a modified auxiliary equation (41) which is developed in this paper. To this end, let us write two scalar and autonomous equations as follows:(73)z˙2=κ+z2+L+z2+F0,z2t0=z0=V−1x0,(74)z˙3=κ−z3+L−z3+F−,z3t0=z0=V−1x0,where κ+=α1+supt≥t0G−t,κ−=α1+inft≥t0G−t,L+z=supt≥t0Lt,z,L−z=inft≥t0Lt,z,z≥0, F−=F0inft≥t0ηt.
In turn, we assume in the remainder of this section that (73) and (74) possess unique solutions, zit,z0,∀t≥t0,z0≥0,i=2,3.
Clearly, the right sides of the last two scalar, autonomous, and integrable equations are bound from above and below the right side in (41), which, due to the comparison principle, implies that z3t,z0≤zt,z0≤z2t,z0, ∀t≥t0, where zt,z0 is a solution to (41). In turn, the structures of solutions to the autonomous and scalar equations (73) and (75) are determined by the location and stability of their fixed solutions. Application of such reasoning to (73) yields sufficient conditions for the boundedness and stability of solutions to (41) and (43), respectively, which, in turn, embraces the corresponding statements for solutions to (1) and (2). Conversely, resolving the behavior of solutions to (74) brings the necessary conditions for the boundedness and stability of solutions to (41) and (43), respectively.
First, we presume below thatκ+<0 and κ−<0 since, otherwise, the right sides of (73) and (74) become positive which entails that limt⟶∞zit,z0=∞,∀z0≥0i=2,3. The last conditions imply that in these cases, our bounds, i e., zit,z0, become overconservative if values of t become relatively large.
Next, we recap the akin sufficient conditions in the following.Theorem 9.
Presume that equations (2), (41), (43), and (73) possess unique solutions for ∀x0∈H⊂ℝn and ∀z0=V−1x0. Additionally, assume (I) that F0=0, κ+<0, and zc≠0 is a unique unstable fixed solution to (73). Then, the trivial solution to (43) is asymptotically stable, i.e., limt⟶∞zt,z0=0,∀z0∈0,zc, and, consequently, zc bounds from above the stability region of the trivial solution to (43). This, in turn, implies that the trivial solution to (2) is asymptotically stable and limt⟶∞xt,x0=0,∀x0∈E−zc, i.e., E−zc is enclosed into stability region of the trivial solution to (2).Next, let us assume thatF0>0, κ+<0, and (II) (73) has two fixed solutions 0<zc2<zc1 corresponding to simple roots of its right side. Then, zc1 and zc2 are unstable and stable fixed solutions to (73), respectively. Furthermore, a solution to (41) is zt,z0≤zci,∀z0≤zci,i=1,2, and limt⟶∞z2t,z0=zc2,∀z0<zc1. This, consequently, yields that(75)xt,x0≤Vzc1,∀x0∈Ezc1Vzc2,∀x0∈Ezc2,∀t≥t0,limt⟶∞supxt,x0≤Vzc2,∀x0∈E−zc1.Lastly, let us assume thatF0>0, κ+<0, and (III) zc1=zc2=zC. Then, a solution to (41), zt,z0≤zC,∀z0≤zC, in turn, implies that xt,x0≤VzC,∀x0∈E−zC,∀t≥t0.Proof.
The proof of this statement immediately follows from the qualitative analysis of the behavior of solutions to (73) which bound from above solutions to either (41) or (43) and the subsequent application of inequality (42) that, in turn, bounds time histories of xt,x0.
Let us, for instance, show thatzc1 and zc2 are unstable and stable fixed solutions to (73). In fact, we get from (73) that z˙2|z2=0=F0>0 which implies that continuous right side of (73) is a positive function for ∀z<zc2 since zc2 is a simple and smallest root of equation, κ+z2+L+z2+F0=0. Hence, zc2 is a stable and, in turn, zc1 is an unstable fixed solution to (73).
Application of similar reasoning to (74) leads to necessary conditions for the boundedness and stability of solutions to (41) and (43), respectively, which can aid simulations of these equations. Nonetheless, these conditions do not directly endorse the consistent properties of solutions to (1) or (2). Yet, the utility of the lower bounds of solutions to (41) can simplify numerical resembling of the threshold initial value of solutions to this equation. In fact, assume that z⌢0 is a threshold initial value for (74) such that, for instance, limt⟶∞z3t,z0<∞,∀z0≤z⌢0 and limt⟶∞z3t,z0=∞,∀z0>z⌢0. Then, limt⟶∞zt,z0=∞,∀z0>z⌢0 as well. Thus, z⌢0 yields an upper estimate of the actual threshold initial value of solutions to (41) which can be accessed through analytical reasoning and aid the approximate resembling of the boundaries of trapping/stability regions for this equation.
Next, we endorse an additional lower bound for solutions to (43) through the utility of the integrable Bernoulli equation which is bounded from below the right side of (43). In fact, assume that function Lt,z in this equation can be written as follows: Lt,z=L1tz+Lntzn+Pt,z, where functions, L1t,Lnt>0,∀t≥t0, 1<n∈ℕ, are continous and function Pt,z≥0, z≥0, ∀t≥t0 is continuous and Lipschitz in z. This enfolds the following Bernoulli equation which assumes analytical solutions:(76)z˙4=α1+G−t+L1tz4+L2tz4n,z4t0,z0=z0=V−1x0.
Clearly,z4t,z0≤zt,z0,∀t≥t0, and z4t,z0 should be close to zt,z0 if Pt,z is sufficiently small. Thus, the threshold initial value procured in the testing of solutions to the integrable equation (76) bounds from above actual threshold value which can be retrieved in simulations of (43).


## 5. Simulations

This section applies the developed methodology to decipher the evolutions of solutions’ norms and estimate the regions of stability/boundedness of two nonlinear systems with time-dependent nonperiodic coefficients which are common in various applications [40]. These systems comprise of two coupled Van der Pol-like or Duffing-like oscillators with variable coefficients. In both cases, the direct utility of our prior technique [13] is compromised since in these cases, limt⟶∞ct=∞.


### 5.1. Coupled System of Van der Pol-Like Oscillators

(32) can be turned into such a system if, for instance, we assume that(77)A=0100−ω12+d−α1d00001d0−ω22+d−α2,G∗=0000−g21t000000000−g43t0,F∗t=0F∗1t0F∗2tT,(78)f∗=0−μ1x230−μ2x43T,where α1=0.4,α2=0.2,ω12=1,ω22=4, g21=a1sinr1t+a2sinr2t,g43=b1sins1t+b2sins2t, a1=a2=b1=b2=0.1,r1=3.14,r2=6.15,s1=3.1,s2=6.28, F∗i=F0isinqit,i=1,2, q1=5.43,q2=10. Note that most of our simulations assume that μ1=μ2=18.2.Next, letV=vij, i,j=1,...,4, be the eigenvector matrix for matrix A that is defined above, and in (33), Λ=diagλ1,λ1∗,λ2,λ2∗,λk=αk+iβk,λk∗=αk−iβk,k=1,2,α1≥α2:(79)f=V−10−μ1∑k=14v2kyk30−μ2∑k=14v4kyk3T.As in Section2, Gt=V−1G∗V , Ft=V−1F∗t.Finally, we set in (41), G−=G−ImdiagG and Lt,z=V−1μ1δ1+μ2δ2z3, where(80)δ1=μ1∑k=14absv2k3,δ2=μ2∑k=14absv4k3.As in Section1.2, absq=a2+b2 if q=a+ib, and absq=q if Imq=0.Let us show briefly how to derive listed above formula forLt,z. Assume that Yi=∑k=14vikyk,i=2,4. Then, (79) can be written as follows: f=V−1Y, where Y=0−μ1Y230−μ2Y43T. Subsequently, application of standard inequalities, V−1Y≤V−1Y, Y2≤Y1=absμ1Y23+absμ1Y43 and ykm≤ym, m∈ℕ, returns that Y1=absμ1∑k=14v2kyk3+absμ2∑k=14v4kyk3=μ1δ2+μ2δ4y3. Consequently, V−1Y≤V−1μ1δ2+μ2δ4z3.The fourth-order Runge–Kutta method with variable step size was used in our simulations. To estimate the boundaries of the trapping/stability regions, the initial fourth-order systems were written in double polar-like coordinates as follows:(81)x1t,x0=r1t,r0,φ0cosφ1t,r0,φ0,x2t,x0=r1t,r0,φ0sinφ1t,r0,φ0,x3t,x0=r2t,r0,φ0cosφ2t,r0,φ0,x4t,x0=r2t,r0,φ0cosφ2t,r0,φ0,r0=r1,0r2,0T,φ0=φ1,0φ2,0T.To reduce the scope of simulations, we present below only computations of projections of the boundedness/stability regions on some coordinate planes. In our simulations, both angle coordinates were discretized with step equal toπ/60. For each fixed set of angle coordinates, the radial ones were adjusted sequentially to approximate the values located on the region’s boundary. We started from small radial values and advanced them until rapid increases in xt,x0 in two consecutive time steps were numerically detected. Similar approach was used to approximate the threshold values in simulations of our scalar auxiliary equations. The running time in the latter computations practically does not depend upon dimension of the original system and is reduced further since zt,z0 monotonically increases in z0 for ∀t≥t0. In contrast, the running time required to estimate the boundary of trapping/stability regions increases with dimensions as mn, where m is the number of points that were taken to discretize a phase-space variable and n is the number of these coordinates.The results of simulation of the Van der Pol-like system with the indicated parameters are partly shown in Figure 1. Figure1(a) contrasts in solid, dashed, and dotted line time-histories of the norms of solutions to equations (1.2), (2.10), and (2.7), respectively, which correspond to a coupled Van der Pol-like system. Note that in this case, x0 is located in the main part of the actual stability region and afar from its boundary. Clearly, the utility of (43) notably enhances the estimates delivered through the application of a more conservative equation (39) with F0=0 and adequately estimates the norm of actual solutions to (2). In the subsequent simulations presented in this section, we set that μ1=μ2=18.2. Figures 1(b)–1(d) display in dashed and solid lines the projections on some coordinate planes of the boundaries of stability regions of trivial solutions to equations (2) and (43) corresponding to the Van der Pol-like system. In fact, simulations of (43) estimate the threshold value of z0, i.e., z¯, which, in turn, defines the ellipsoid imbedded in the stability region through application of the formula, z¯=V−1x0. Clearly, the attained estimates of the boundary of stability region turn out to be fairly conservative if the structure and parameters of the model are defined precisely. Yet, the practical merit of these estimates becomes more apparent for systems under uncertainty, where, for instance, it is assumed that at most the norms of some model components or parameters are defined precisely.Figure 1
(a) contrasts in solid, dashed, and dotted lines time histories of the norms of solutions to equation (2) corresponding to a coupled system of Van der Pol-like oscillators and the relevant equations (43) and (39) with F0=0, respectively. (b) displays in dashed and solid lines projections of the boundary of stability region on x1×x2 -plane obtained in simulations of equations (2) corresponding to a coupled system of Van der Pol-like oscillators and the relevant equation (43), respectively. (c) displays in dashed and solid lines projections of the boundary of stability region on x1×x3 -plane obtained in simulations of equations (2) corresponding to a coupled system of Van der Pol-like oscillators and the relevant equation (43), respectively. (d) displays in dashed and solid lines projections of the boundary of stability region on x1×x4-plane obtained in simulations of equations (2) corresponding to a coupled system of Van der Pol-like oscillators and the relevant equation (43), respectively.
(a)(b)(c)(d)


### 5.2. Coupled System of Duffing-Like Oscillators

For this system, (77) remains intact, but (78) should be altered as follows:(82)f×=0−μ1x130−μ2x33T,which, in turn, requires some alteration in (80), i.e., v2k⟶v1k and v4k⟶v3k, respectively. Figure 2 displays in dashed and solid lines the projections on some coordinate planes of the boundaries of stability regions of the trivial solutions to (2) and (41) corresponding to a coupled system of Duffing-like oscillators, respectively. As in Section 5.1, the current estimates turn out to be rather conservative for precisely defined systems but appear to be more compelling for systems under uncertainties. Note that the simulations on Figure 2(b) are started at t0=1.43, whereas other plots in this figure are simulated with t0=0. Clearly, comparison of Figures 2(a) and 2(b) shows that the shapes of the attractors’ boundaries in these figures are somewhat affected by the initial time moment. However, such dependence might become less pronounced for other values of t0 since time-dependent coefficients are almost periodic in this model.Figure 2
(a) Displays in dashed and solid lines projections of the boundary of stability region onx1×x2 -plane obtained in simulations of equations (2) corresponding to a coupled system of Duffing-like oscillators and the relevant equation (43), respectively (b) Displays in dashed and solid lines projections of the boundary of stability region on x1×x2 -plane obtained in simulations of equations (2) corresponding to a coupled system of Duffing-like oscillators and relevant to equation (43), respectively. Note that in this case the simulations start at t0=1.43.(c) Displays in dashed and solid lines projections of the boundary of stability region on x1×x3 -plane obtained in simulations of equations (2) corresponding to a coupled system of Duffing-like oscillators and the relevant equation (43), respectively.(d) Displays in dashed and solid lines projections of the boundary of stability region on x1×x4 -plane obtained in simulations of equations (2) corresponding to a coupled system of Duffing-like oscillators and the relevant equation (43), respectively.
(a)(b)(c)(d)Next, we review some results of simulations on estimating the trapping regions for nonhomogeneous equation (1). Figure 3 displays in dashed and solid lines projections on the coordinate planes of the boundaries of trapping regions about zero for (1) and (41) that correspond to a coupled Duffing-like system with variable coefficients. In these simulations. In these simulations, μ1=μ2=0.3 and F01=0.01,F02=0. Our outcomes in these simulations are comparable to the ones that were presented prior and show that our fairly conservative estimates of the boundary of trapping region of (1) may become more practically appealing for a system under uncertainty. Yet, the utility of these estimates for precisely defined systems provides a numerical affirmation of our boundedness and stability criteria. Nonetheless, combining the techniques outlined in this paper with the pertinent successive approximations should considerably refine the accuracy of the corresponding estimates (see [14], where such a strategy was developed and implemented for suitable systems).Figure 3
(a) displays in dashed and solid lines projections onx1×x2-plane of the boundary of trapping region obtained from simulations of the nonhomogeneous equation (1) describing the coupled model of Duffing-like oscillators and the corresponding equation (43), respectively (b) Displays in dashed and solid lines projections on x1×x3 -plane of the boundary of trapping region obtained from simulations of the nonhomogeneous equation (1) describing the coupled model of Duffing-like oscillators and the corresponding equation (43), respectively. (c) Displays in dashed and solid lines projections on x1×x4-plane of the boundary of trapping region obtained from simulations of the nonhomogeneous equation (1) describing the coupled model of Duffing-like oscillators and the corresponding equation (43), respectively (d) Displays in dashed and solid lines projections on x2×x4-plane of the boundary of trapping region obtained from simulations of the nonhomogeneous equation (1) describing the coupled model of Duffing-like oscillators and the corresponding equation (43), respectively.
(a)(b)(c)(d)Comparison of our latter [13] and current techniques shows that the current one gains computational efficiency and is applicable to an essentially wider class of systems. Nonetheless, the accuracy of this technique decreases if α1 is a small value, sup∀t≥t0G−t is a large value, and sup∀t≥t0Lt,z becomes sufficiently large for relevant values of z. The accuracy of our prior technique is practically nonsensitive to variations of sup∀t≥t0G−t but decreases with growth of sup∀t≥t0Lt,z as well. However, this prior technique also requires numerical simulation of wt, pt, and ct which is avoided in the current approach.


## 5.1. Coupled System of Van der Pol-Like Oscillators

(32) can be turned into such a system if, for instance, we assume that(77)A=0100−ω12+d−α1d00001d0−ω22+d−α2,G∗=0000−g21t000000000−g43t0,F∗t=0F∗1t0F∗2tT,(78)f∗=0−μ1x230−μ2x43T,where α1=0.4,α2=0.2,ω12=1,ω22=4, g21=a1sinr1t+a2sinr2t,g43=b1sins1t+b2sins2t, a1=a2=b1=b2=0.1,r1=3.14,r2=6.15,s1=3.1,s2=6.28, F∗i=F0isinqit,i=1,2, q1=5.43,q2=10. Note that most of our simulations assume that μ1=μ2=18.2.Next, letV=vij, i,j=1,...,4, be the eigenvector matrix for matrix A that is defined above, and in (33), Λ=diagλ1,λ1∗,λ2,λ2∗,λk=αk+iβk,λk∗=αk−iβk,k=1,2,α1≥α2:(79)f=V−10−μ1∑k=14v2kyk30−μ2∑k=14v4kyk3T.As in Section2, Gt=V−1G∗V , Ft=V−1F∗t.Finally, we set in (41), G−=G−ImdiagG and Lt,z=V−1μ1δ1+μ2δ2z3, where(80)δ1=μ1∑k=14absv2k3,δ2=μ2∑k=14absv4k3.As in Section1.2, absq=a2+b2 if q=a+ib, and absq=q if Imq=0.Let us show briefly how to derive listed above formula forLt,z. Assume that Yi=∑k=14vikyk,i=2,4. Then, (79) can be written as follows: f=V−1Y, where Y=0−μ1Y230−μ2Y43T. Subsequently, application of standard inequalities, V−1Y≤V−1Y, Y2≤Y1=absμ1Y23+absμ1Y43 and ykm≤ym, m∈ℕ, returns that Y1=absμ1∑k=14v2kyk3+absμ2∑k=14v4kyk3=μ1δ2+μ2δ4y3. Consequently, V−1Y≤V−1μ1δ2+μ2δ4z3.The fourth-order Runge–Kutta method with variable step size was used in our simulations. To estimate the boundaries of the trapping/stability regions, the initial fourth-order systems were written in double polar-like coordinates as follows:(81)x1t,x0=r1t,r0,φ0cosφ1t,r0,φ0,x2t,x0=r1t,r0,φ0sinφ1t,r0,φ0,x3t,x0=r2t,r0,φ0cosφ2t,r0,φ0,x4t,x0=r2t,r0,φ0cosφ2t,r0,φ0,r0=r1,0r2,0T,φ0=φ1,0φ2,0T.To reduce the scope of simulations, we present below only computations of projections of the boundedness/stability regions on some coordinate planes. In our simulations, both angle coordinates were discretized with step equal toπ/60. For each fixed set of angle coordinates, the radial ones were adjusted sequentially to approximate the values located on the region’s boundary. We started from small radial values and advanced them until rapid increases in xt,x0 in two consecutive time steps were numerically detected. Similar approach was used to approximate the threshold values in simulations of our scalar auxiliary equations. The running time in the latter computations practically does not depend upon dimension of the original system and is reduced further since zt,z0 monotonically increases in z0 for ∀t≥t0. In contrast, the running time required to estimate the boundary of trapping/stability regions increases with dimensions as mn, where m is the number of points that were taken to discretize a phase-space variable and n is the number of these coordinates.The results of simulation of the Van der Pol-like system with the indicated parameters are partly shown in Figure 1. Figure1(a) contrasts in solid, dashed, and dotted line time-histories of the norms of solutions to equations (1.2), (2.10), and (2.7), respectively, which correspond to a coupled Van der Pol-like system. Note that in this case, x0 is located in the main part of the actual stability region and afar from its boundary. Clearly, the utility of (43) notably enhances the estimates delivered through the application of a more conservative equation (39) with F0=0 and adequately estimates the norm of actual solutions to (2). In the subsequent simulations presented in this section, we set that μ1=μ2=18.2. Figures 1(b)–1(d) display in dashed and solid lines the projections on some coordinate planes of the boundaries of stability regions of trivial solutions to equations (2) and (43) corresponding to the Van der Pol-like system. In fact, simulations of (43) estimate the threshold value of z0, i.e., z¯, which, in turn, defines the ellipsoid imbedded in the stability region through application of the formula, z¯=V−1x0. Clearly, the attained estimates of the boundary of stability region turn out to be fairly conservative if the structure and parameters of the model are defined precisely. Yet, the practical merit of these estimates becomes more apparent for systems under uncertainty, where, for instance, it is assumed that at most the norms of some model components or parameters are defined precisely.Figure 1
(a) contrasts in solid, dashed, and dotted lines time histories of the norms of solutions to equation (2) corresponding to a coupled system of Van der Pol-like oscillators and the relevant equations (43) and (39) with F0=0, respectively. (b) displays in dashed and solid lines projections of the boundary of stability region on x1×x2 -plane obtained in simulations of equations (2) corresponding to a coupled system of Van der Pol-like oscillators and the relevant equation (43), respectively. (c) displays in dashed and solid lines projections of the boundary of stability region on x1×x3 -plane obtained in simulations of equations (2) corresponding to a coupled system of Van der Pol-like oscillators and the relevant equation (43), respectively. (d) displays in dashed and solid lines projections of the boundary of stability region on x1×x4-plane obtained in simulations of equations (2) corresponding to a coupled system of Van der Pol-like oscillators and the relevant equation (43), respectively.
(a)(b)(c)(d)


## 5.2. Coupled System of Duffing-Like Oscillators

For this system, (77) remains intact, but (78) should be altered as follows:(82)f×=0−μ1x130−μ2x33T,which, in turn, requires some alteration in (80), i.e., v2k⟶v1k and v4k⟶v3k, respectively. Figure 2 displays in dashed and solid lines the projections on some coordinate planes of the boundaries of stability regions of the trivial solutions to (2) and (41) corresponding to a coupled system of Duffing-like oscillators, respectively. As in Section 5.1, the current estimates turn out to be rather conservative for precisely defined systems but appear to be more compelling for systems under uncertainties. Note that the simulations on Figure 2(b) are started at t0=1.43, whereas other plots in this figure are simulated with t0=0. Clearly, comparison of Figures 2(a) and 2(b) shows that the shapes of the attractors’ boundaries in these figures are somewhat affected by the initial time moment. However, such dependence might become less pronounced for other values of t0 since time-dependent coefficients are almost periodic in this model.Figure 2
(a) Displays in dashed and solid lines projections of the boundary of stability region onx1×x2 -plane obtained in simulations of equations (2) corresponding to a coupled system of Duffing-like oscillators and the relevant equation (43), respectively (b) Displays in dashed and solid lines projections of the boundary of stability region on x1×x2 -plane obtained in simulations of equations (2) corresponding to a coupled system of Duffing-like oscillators and relevant to equation (43), respectively. Note that in this case the simulations start at t0=1.43.(c) Displays in dashed and solid lines projections of the boundary of stability region on x1×x3 -plane obtained in simulations of equations (2) corresponding to a coupled system of Duffing-like oscillators and the relevant equation (43), respectively.(d) Displays in dashed and solid lines projections of the boundary of stability region on x1×x4 -plane obtained in simulations of equations (2) corresponding to a coupled system of Duffing-like oscillators and the relevant equation (43), respectively.
(a)(b)(c)(d)Next, we review some results of simulations on estimating the trapping regions for nonhomogeneous equation (1). Figure 3 displays in dashed and solid lines projections on the coordinate planes of the boundaries of trapping regions about zero for (1) and (41) that correspond to a coupled Duffing-like system with variable coefficients. In these simulations. In these simulations, μ1=μ2=0.3 and F01=0.01,F02=0. Our outcomes in these simulations are comparable to the ones that were presented prior and show that our fairly conservative estimates of the boundary of trapping region of (1) may become more practically appealing for a system under uncertainty. Yet, the utility of these estimates for precisely defined systems provides a numerical affirmation of our boundedness and stability criteria. Nonetheless, combining the techniques outlined in this paper with the pertinent successive approximations should considerably refine the accuracy of the corresponding estimates (see [14], where such a strategy was developed and implemented for suitable systems).Figure 3
(a) displays in dashed and solid lines projections onx1×x2-plane of the boundary of trapping region obtained from simulations of the nonhomogeneous equation (1) describing the coupled model of Duffing-like oscillators and the corresponding equation (43), respectively (b) Displays in dashed and solid lines projections on x1×x3 -plane of the boundary of trapping region obtained from simulations of the nonhomogeneous equation (1) describing the coupled model of Duffing-like oscillators and the corresponding equation (43), respectively. (c) Displays in dashed and solid lines projections on x1×x4-plane of the boundary of trapping region obtained from simulations of the nonhomogeneous equation (1) describing the coupled model of Duffing-like oscillators and the corresponding equation (43), respectively (d) Displays in dashed and solid lines projections on x2×x4-plane of the boundary of trapping region obtained from simulations of the nonhomogeneous equation (1) describing the coupled model of Duffing-like oscillators and the corresponding equation (43), respectively.
(a)(b)(c)(d)Comparison of our latter [13] and current techniques shows that the current one gains computational efficiency and is applicable to an essentially wider class of systems. Nonetheless, the accuracy of this technique decreases if α1 is a small value, sup∀t≥t0G−t is a large value, and sup∀t≥t0Lt,z becomes sufficiently large for relevant values of z. The accuracy of our prior technique is practically nonsensitive to variations of sup∀t≥t0G−t but decreases with growth of sup∀t≥t0Lt,z as well. However, this prior technique also requires numerical simulation of wt, pt, and ct which is avoided in the current approach.


## 6. Conclusion and Forthcoming Studies

In [13], we derived a scalar nonlinear auxiliary equation with solutions bounded from above the norm of solutions to some nonlinear and nonautonomous systems and applied this technique to the assessment of the boundedness and stability of solutions to such systems. Yet, the application domain of this approach is constrained since the underlying auxiliary equation contains a function that in some cases can approach infinity with t⟶∞. The current paper presents a novel technique casting the auxiliary equation in a modified form escaping this limitation for a wide class of systems arising in applications. Furthermore, the modified auxiliary equation is simpler than its former counterpart and more efficient in computation. Still, in the common application domain, the current estimates frequently turn out to be more conservative than their prior analogs.Next, we present various novel boundedness and stability criteria which are stemmed from analysis of the modified nonlinear auxiliary equation and its linear counterpart. The latter approach leads to close-form and simplified boundedness and stability criteria that should be valuable in resolving relevant control problems.Lastly, we authenticate our study using inclusive simulations of systems with typical dissipative and conservative nonlinear components that were intractable to our prior technique. The simulations show that our upper estimates of the norms of solutions to nonlinear and nonautonomous systems turn out to be adequate if they stem from the central parts of the corresponding trapping/stability regions. Nonetheless, the boundaries of trapping/stability regions are estimated rather conservatively by our current technique if the corresponding systems are defined precisely. Yet, such estimates turn out to be more appealing in applications for systems under uncertainty.Note also that the precision of our estimates can be substantially improved through the integration of our current technique with a methodology of successive approximations that was outlined in [14]. This should provide bilateral bounds for the norms of solutions to a broad class of time-varying and nonlinear systems and enhance recursively the accuracy of estimation of their trapping/stability regions.



---

*Source: 5098677-2022-06-29.xml*