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

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.


Introduction and Motivation Example
1.1. Background. Analysis of boundedness and stability of nonlinear systems with variable and nonperiodic coe cients 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 con rmable necessary and su cient conditions of local stability of the trivial solution to homogeneous systems of such kind (see, e.g., [1][2][3][4][5][6][7]). It was shown by Lyapunov [8] that under some additional conditions, the trivial solution to a homogeneous system with time-varying coe cients 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 con rmation 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 coe cients and alter its stability. is 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 su cient 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 su cient 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.
Due to continuity of the right side of (1), the last assumption implies that x(t, x 0 ) is a continuous and continuously differential function that is bounded for any t ∈ [t 0 , t * ], t * < ∞. Consequently, the remainder of this paper focuses on behavior of x(t, x 0 ) for t ⟶ ∞.
We also will assess stability of the trivial solution to the following homogeneous equation: and its linear counterpart: Clearly, the solution to (3) can be written as x(t, x 0 ) � W(t, t 0 )x 0 , where W(t, t 0 ) � w(t)w − 1 (t 0 ) and w(t) are transition and fundamental solution matrices for (3), respectively [1,4].
A more general but less tractable condition on the norm of the transition matrices was presented in [1,2] in the form where it was shown that (4) and (7) and the following condition: 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][16][17][18][19][20][21][22][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. e problem of estimating the stability regions of autonomous nonlinear systems has attracted numerous publications in the last few decades [24][25][26][27][28][29][30][31][32][33][34][35][36][37][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, ‖w(t 0 )‖ � 1, we derive the following inequality: ‖x(t, x 0 )‖ ≤ X(t, X 0 ), ∀t ≥ t 0 , where x(t, x 0 ) is a solution to (1) and X(t, X 0 ): [t 0 , ∞) × R ≥0 ⟶ R ≥0 (R ≥0 is a set of nonnegative real numbers) is a solution to a scalar equation we derived in [13]: Note that and 2 Mathematical Problems in Engineering In turn, Pinsky and Koblik [13] also introduced a nonlinear extension of the Lipschitz continuity condition which was presented as follows: where L: is continuous function in t and ‖x‖ and locally Lipschitz in ‖x‖, L(t, 0) � 0, and Ω x 2 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 Ω x 2 error term, which, for instance, can be written in Lagrange form. In the former case, Ω x 2 ≡ R n and in the latter case, Ω x 2 ≡ R n if the error term is also globally bounded for ∀x ∈ R n . us, condition Ω x 2 ≡ R 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 define L(t, ‖x‖) for a simple but representative vector function. Assume that x � [x 1 x 2 ] T and a vector function f ∈ R 2 is defined, e g., as follows: where we use that |x n i | ≤ ‖x‖ n , i � 1, 2, n ∈ N, and N 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: To reduce the notation, we adopt in (13) and throughout this paper that z(t, t 0 , z 0 ) � z(t, z 0 ).
In turn, it was shown in [13] that where x(t, x 0 ), X(t, X 0 ), and z(t, z 0 ) are solutions of (1), (9), and (13), respectively. us, (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][4][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.
is 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 c(t) < ∞, ∀t ≥ t 0 , which considerably limits the scope of its applications. Nonetheless, it follows from (11) that frequently lim t⟶∞ c(t) � ∞ even if A is a stable and time-invariant matrix. e 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 with c(t) ≡ 1 and p(t) ≡ const under some conditions that are frequently met in various applications. is 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 w(t), p(t), and c(t) 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][2][3][4][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].
x 0 )) , 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 � lim where σ i are the singular values of the fundamental solution matrix of (3), and for linear systems, max i ϕ l,i � ϕ [1]. is lets us to bound the transition matrix of (3) as follows: (15) 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 where function g(t, u 1 ) is continuous in t and locally Lipschitz in u 1 for ∀u 1 ∈ ℘ ⊂ R. Suppose that the solution to this equation is u 1 (t, u 10 ) ∈ ℘, ∀t ≥ t 0 . Next, consider a differential inequality where D + u 2 denotes the upper right-hand derivative in t of u 2 (t, u 20 ) [4] and u 2 (t, u 20 ) ∈ ℘, ∀t ≥ t 0 . en, u 1 (t,

Mathematical Problems in Engineering
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, I 1 (t 0 ), is called a trapping region of equation (1) if Clearly, this definition acknowledges that I 1 is the invariant set of (1) containing zero.

Definition 2.
A connected and open set of all initial vectors, I 2 (t 0 ), that includes zero vector, is called a region of stability of the trivial solution to (2) if condition x 0 ∈ I 2 (t 0 ) implies that x(t, x 0 ) is stable.

Definition 3.
A connected and open set of all initial vectors, I 3 (t 0 ), that includes zero vector, is called a region of asymptotic stability of the trivial solution to (2) if condition

Motivation Example.
e first part of this section derives a novel auxiliary equation with c(t) � 1 and p(t) � const for a simple planar system. e 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 break out our current approach into a sequence of straightforward steps for a simple version of system (1) where Assume also that matrix A is diagonalizable and eigenvalue and eigenvector matrices of can be applied to the stability analysis of such a planar system if (11) implies that c(t) < ∞, ∀t ≥ t 0 which, in general, is difficult to warranty in advance. us, we outline a different technique recasting the auxiliary equation in the form, where c(t) ≡ 1 and p(t) ≡ const.
(1) Firstly, let us write our planar model of system (1) in eigenbasis of A as follows: where y ∈ C 2 , C 2 is a two-dimensional space of complex numbers, x � Vy, (2) Next, we rewrite the last equation as follows: where I is an identity matrix, b � diag β − β , a � diag α α , and λ ∈ R is going to be determined subsequently.
In order, let us select a linear subsystem of (19) with an underlined diagonal matrix as For this subsystem, a diagonal fundamental solution matrix can be written as follows: Same reasoning yields that ‖w − 1 (t)‖ � e − λ(t− t 0 ) which shows that c(t) � ‖w(t)‖‖w − 1 (t)‖ � 1, and due to (10), us, 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 p(t)z � λ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: en, the utility of standard inequality yields where abs(q) � is lets us to write (13) as follows: where σ � ε‖V − 1 ‖( 2 j�1 abs(v 2j )) 3 . (4) Lastly, we select λ to maximize the degree of stability of a linear equation _ z � (λ + ‖a − λI‖)z which minimizes conservatism of our estimates and reduces to the following condition.
since a is a diagonal matrix with equal eigenvalues. Resolving (24) yields that λ � α which casts (23) in the following form: For further references, we write a homogeneous counterpart to (25) as follows: Let us recall now that, due to (14), x 0 ) and z(t, z 0 ) 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, z 3 ≤ (z) 2 z, z ∈ [0, z], z > 0. Application of the last inequality brings the following linear equation: where s(t, z) � α + ‖G(t)‖ + σ(z) 2 . Note that solutions to (27) can be used to estimate solutions to (26) In turn, to simplify the assessment of stability of (27), we set that sup Subsequently, the solutions to (26) can be bounded as follows: where z max is the maximal value of z for which (28) holds. e 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: Clearly, the right side of (29) bounds from above the right side of (27), whereas the initial values for both equations are equal. us, Next, we infer that (28) implies the asymptotic stability of the trivial solution to (26). In fact, due to the comparison principle, (26) and (27). (28) holds and z 0 < z. Consequently, the trivial solutions to (25) and, in turn, to our planar model of equation (2) are asymptotically stable.

Mathematical Problems in Engineering
Furthermore, solving for z 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 e 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 z 0 < 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 z 0 ∈ (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: Section 3 generalizes this last formula and derives the boundedness/stability criteria and estimates the trapping/ stability regions for multidimensional (1) and (2).

Modified Auxiliary Equation
is section derives a modified auxiliary equation with c(t) � 1 and p(t) � 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 of B(t) exists and is defined as follows: Additionally, we assume that A(t 0 ) ∈ R n×n is nonzero and diagonalizable matrix ∀t 0 ∈ Τ except possibly for some isolated values of t 0 . is last condition is met for any generic set of matrices A(t 0 ) 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 of where R >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, (2) Next, we write (1) as follows: where 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, To escape this shortcoming, we rewrite (32) into the eigenbasis of A as follows: where λ ∈ R is defined below. en, 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: where I is identity matrix. Henceforth, In turn, application of (10) and (11) (35), we write equation (13) for (34) as follows: where is a bounded neighborhood of y ≡ 0, and L: � [t 0 , t∞) × R ≥0 ⟶ R ≥0 is continuous in t and z and is locally Lipschitz in z function with L(t, 0) ≡ 0. Note that this function can be developed through application of inequality (12) to function, 6 Mathematical Problems in Engineering Vy). Hence, (12) in this case takes the following form: For developing L(t, ‖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 f(t, y) is a piece-wise polynomial function or can be approximated by such function with bounded Ω y 2 Lagrange-type error term. In the former case, Ω y 2 ≡ C n , and in the latter case, Ω y 2 ≡ C n if the error term is bounded in C n . To straightforward further references, we assume in remainder of this paper that Ω y 2 ≡ C n . (4) In the sequel, we define λ using the following condition: min λ (λ + ‖α − λI‖), which maximizes the degree of stability of a scalar linear equation, Clearly, the last condition yields that λ � α n and, consequently, (λ + ‖α − λI‖) |λ�α n � α 1 (see proof in appendix of this paper). is lets us to write (36) as follows: (5) To develop a less conservative version of (39), we set that matrix D(t) � Im(diagG(t)), D ∈ R n×n , β + (t) � β + D(t), and G − � G − i D and rewrite (34) as follows: e fundamental matrix of solution for equation, _ y � (λI + iβ + (t))y, again can be written as follows: w + (t) � exp(λI + iβ + (t))(t − t 0 ) since β + is a diagonal matrix. As prior, this implies that ‖w , c(w + (t)) � 1, and p(w + (t)) � α 1 . Lastly, a less conservative counterpart of (39) can be written as follows: In order, we encapsulate our derivations and the underlying assumptions in the following.

Theorem 1. Assume that B(t) is a continuous matrix, and the average of B(t),
i.e., matrix A(t 0 ) ≠ 0, exists and is diagonalizable for ∀t 0 ∈ Τ except possibly some isolated values of t 0 . Also, assume that function F * (t) is continuous and function f * (t, x) is continuous in both variables and inequality (37) holds with Ω y 2 ≡ C n . Lastly, we assume that equations (1), (39), and (41) possess unique solutions for ∀t ≥ t 0 , ∀x 0 ∈ H ⊂ R n , and ∀z 0 � ‖V − 1 x 0 ‖. en, the norms of solutions to equation (1) are bounded by matching solutions to the scalar equations (39) or (41) as follows: Proof. e 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). e next set of conditions implies inequality (37) with Ω y 2 ≡ R n and, in turn, equations (36), (39), and (41). e last assumption assures that functions x(t, x 0 ) and z(t, z 0 ) exist for ∀t ≥ t 0 and the specified initial values.
Note that the above inferences naturally filter out the effect of imaginary components in matrix Λ on evolution of z(t, z 0 ).
is equates the degree of stability of equation _ x � Ax, which is measured by the maximal real part of eigenvalues of A, and equation _ z � α 1 z. Furthermore, (41) is currently defined explicitly, whereas the definition of its counterpart given in [13] involves numerical evaluation of w(t), p(w(t)), and c(w(t)).
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 w(t) 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: 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: where we assumed that function l 2 : 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 L(t, z). Note that definition of l 2 (t, z) is exceedingly simplified in this scalar case. Furthermore, definition of l 2 (t, z) is also abridged if L(t, z) is a polynomial or piece-wise polynomial function in z which, in former case, can be written as follows: Mathematical Problems in Engineering Clearly, in some cases, l 2 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).

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: where To shorten the notation, we adopted in (45) and frequently, but always, will set below that Clearly, under the conditions of eorem 1 and inequality (44), ‖G − (t)‖, ‖F(t)‖, and l 2 (t, z) are continuous functions which implies that μ(t, z) is continuous in both variables as well.
us, under the condition of eorem 1, (45) admits the following solution: where function Z(t, z 0 , z) is continuous in t and z 0 , As prior, we frequently, but always, are going to use reduce notations, i.e., Z h (t, z) � Z h (t, t 0 , z) and Z F (t, z) � Z F (t, t 0 , z). Clearly, Z h > 0 and θ > 0. us, the last inequality yields that Z F > 0 as well. We are going to use these inequalities in the proof of Lemma 2.

Stability Criteria and Estimation of Stability Region.
To simplify the subsequent references, we will write a homogeneous counterpart to (45) as follows: 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 and asymptotically stable if and only if 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 z(t, z 0 ) ∈ [0, z], ∀t ≥ t 0 since compliance with this condition enables the linearization of (43) via application of (44). Next lemma provides an upper bound of Z(t, z 0 ) which, subsequently, aids estimation of z(t, z 0 ) and, in turn, ‖x(t, x 0 )‖. (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 (47) is stable which yields that  (47) is stable for some z ∈ (0, z max ). en, the trivial solution to (2) is stable and inequality (42) takes the form

Lemma 1. Assume that function μ(t, z) is continuous in t and bounded in z,
where Z(t, z 0 , z), z(t, z 0 ) and x(t, x 0 ) are solutions to equations (47), (43), and (2), respectively, and the region of stability of the trivial solution to (2) includes all values of x 0 containing the interior of ellipsoid which is defined as follows: Proof. Let us show that under conditions of this theorem, z(t, z 0 ) < z, ∀t ≥ t 0 , z ∈ (0, z max ) if z 0 < z/Z s . Pretend on the contrary that t � t * > t 0 is the smallest value of t such that z(t * , z 0 ) � z under prior conditions. en, (44) and, thus, equation (47) hold for ∀t ∈ [t 0 , t * ] which, due to comparison principle [4] and (50), implies that z(t * , z 0 ) ≤ Z (t * , z 0 , z) < z. is contradiction shows that z(t, z 0 ) < z, ∀t ≥ t 0 , ∀z 0 < z/Z s 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 z(t, z 0 ) ≤ Z(t, z 0 , z), ∀t ≥ t 0 , z 0 < z/Z s , z ∈ (0, z max ) which implies that the trivial solution to (43) is stable. e 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). e solution to (43) is stable for all values of z 0 which are contingent by the following inequality: region of stability of the trivial solution to (2) includes all values of x 0 meeting condition (52).
Let us recall that (52) depends upon t 0 -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.
Hence, setting in (52) that Z s (t 0 , z) � 1 prompts more conservative but abridged estimation of the stability region for our planar system which can be written as follows: e last formula was derived also in Section 1.2 by application of a simplified inference.
Theorem 3. Assume that assumptions of eorem 1 are met, l 2 (t, z) is continuous in t and bounded in z function, and equation (47) is asymptotically stable for some z ∈ (0, z max ).
Proof. Literally, under these more conservative conditions, eorem 2 and, hence, (51) hold, which affirms the current statement since in this case lim 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 eorem 1 hold, l 2 (t, z) is continuous in t and bounded in z function, and either (48) or (49) holds as well. en, inequality (51) holds, where Z(t, z 0 , z), z(t, z 0 ), and x(t, x 0 ) 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 eorem 2 and eorem 3, assures this statement □ Corollary 2. Assume that conditions of eorem 1 hold, l 2 (t, z) is continuous in t and bounded in z function, and where z v is the maximal value of z for which the last inequality holds. en, the trivial solution of (2) is asymptotically stable. Inequalities (51) and (52) where m(t, z) � ‖G − (t)‖ + l 2 (t, z) which steers us to the following. □ Corollary 3. Assume that conditions of eorem 1 hold, l 2 (t, z) is continuous in t and bounded in z function, and ϕ(t 0 , z) � α 1 + c(t 0 , z) < 0, ∀z ∈ (0, z ϕ ), t 0 ∈ Τ, where z ϕ is the maximal value of z for which last inequality holds. en, the trivial solution to (2) is asymptotically stable and (51) and (52) hold, where in (52), z max should be exchanged on z ϕ .

□
Remark. Let us represent the last stability criterion in an abridged form. Pretend that ‖G − (t)‖ � g av (t 0 ) + g(t 0 , t), en, it follows from definition of c that c(t 0 , z) � g av (t 0 ) + l 2av (t 0 , z 0 ) which implies that stability condition of Lemma 2 can be written as follows: Mathematical Problems in Engineering α 1 + g av t 0 + l 2av t 0 , z < 0, t 0 ∈ Τ.
(54) 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 f(t, x) is a piece-wise polynomial in x since in this case functions L(t, z) and l 2 (t, z) can be defined in analytical form. In contrast, the utilities of either (6) or (8) require simulations of w(t) and subsequent estimation of the corresponding parameters or function from the simulated data. e 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.

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). ese 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]).
Firstly, let us present an analog to Lemma 1 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
where, as prior, Proof. As is mentioned in Lemma 1, sup since Z h (t 0 , t, z) ≥ 0 and Z F (t 0 , t, z) ≥ 0 which yields (56) e next statement assesses the boundedness of solutions to the nonhomogeneous equation (41) and, in turn, (1). (55) and (56) hold, and l 2 (t, z) is continuous in t and bounded in z function. en, (i) ‖x(t, x 0 )‖ ≤ ∞, ∀t ≥ t 0 if these solutions are emanated from the interior of the ellipsoid that is defined as follows:
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 lim t⟶∞ Z h (t, x 0 ) � 0. Next, two corollaries present simplified sufficient conditions under which the presumptions of eorem 4 are met.
□ Corollary 5. Assume that the conditions of eorem 1 and inequalities (53), (55), and (56) hold, and l 2 (t, z) is continuous in t and bounded in z function. en, the solutions to (1)
Note that in this statement, it is assumed that in (58) and 59, z 1 � min(z ϕ , z ρ ).
Proof. In fact, in this case, (47) is asymptotically stable [1] and the conditions of eorem 4 hold. us, (59) yields that lim 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. us, to simplify our further inferences, we bring an autonomous counterpart to (25) which can be written as follows: where, as prior, we assume that G s � sup (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 Section 1.2, to simplify the assessment of the behavior of solutions to (62), we consider its linear counterpart where − s(z) � − α * + G s + (z) 2 σ. us, (63) can be viewed as a nonhomogeneous analog to (29). Furthermore, let us write a nonhomogeneous counterpart to (27) as follows: us, due to comparison principle, Clearly, the general solutions to the linear scalar (63) and (64) take the form (29) and (27) with z 0 � 1 and Z i,F (t, z 0 , z) are components of forced solution to the corresponding equations. Let us recall that sup To estimate the right side of (56), we recall that (z)(t − τ)) which, with use of (55), lets us bound ρ(t 0 , z) as follows: 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. e former approach is substantially abridged due to the following.

Theorem 5.
Assume that equation (41) possesses a unique solution z(t, z 0 ), ∀t ≥ t 0 , z 0 ≥ 0, and z(t, z 0 ′ ) and z(t, z ′ ′ 0 ) are solutions to (41) with 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, E(z) ⊂ R n , as follows: Also, we assume that zE(z) ⊂ R n− 1 defines the boundaries of these ellipsoids and E − (z) � E(z) − zE(z).
is prompts the following. □ Theorem 6. Assume that equations (2) and (43) possess unique solutions for ∀x 0 ∈ H ⊂ R n and ∀z 0 � ‖V − 1 x 0 ‖, and that the trivial solution to (43) is either stable or asymptotically stable. en, 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). en, set E − (z) is contained in stability region of (2).

E(z) is included into the trapping region of equation
Proof. e proof of this statement follows from inequality (42), where x(t, x 0 ) and z(t, z 0 ) are solutions to equations (1) and (41), respectively.
us, the problem of estimating the trapping/stability regions about point x ≡ 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. en, the boundary of the trapping/stability region for either (1) or (2) is given by the following formula: z � ‖V − 1 x 0 ‖. In turn, the task of simulating the threshold value is markedly simplified since z(t, z 0 ) monotonically increases in z 0 , ∀t ≥ t 0 due to eorem 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.
(69) en, (43) can be written as where and is comprises the following.  (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. en, 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: where In turn, we assume in the remainder of this section that (73) and (74) possess unique solutions, 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 z 3 (t, z 0 ) ≤ z(t, z 0 ) ≤ z 2 (t, z 0 ), ∀t ≥ t 0 , where z(t, z 0 ) 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 lim t⟶∞ z i (t, z 0 ) � ∞, ∀z 0 ≥ 0 i � 2, 3. e last conditions imply that in these cases, our bounds, i e., z i (t, z 0 ), become overconservative if values of t become relatively large.
is, in turn, implies that the trivial solution to (2) is asymptotically stable and lim t⟶∞ ‖x(t, x 0 )‖ � 0, ∀x 0 ∈ E − (z c ), i.e., E − (z c ) is enclosed into stability region of the trivial solution to (2).
e 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 ‖x(t, x 0 )‖.
Let us, for instance, show that z c1 and z c2 are unstable and stable fixed solutions to (73). In fact, we get from (73) that _ z 2|z 2 �0 � F 0 > 0 which implies that continuous right side of (73) is a positive function for ∀z < z c2 since z c2 is a simple and smallest root of equation, κ + z 2 + L + (z 2 ) + F 0 � 0. Hence, z c2 is a stable and, in turn, z c1 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, en, lim t⟶∞ z(t, z 0 ) � ∞, ∀z 0 > z ⌢ 0 as well.
us, 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 L(t, z) in this equation can be written as follows: L(t, z) � L 1 (t)z + L n (t)z n + P(t, z), where functions, L 1 (t), L n (t) > 0, ∀t ≥ t 0 , 1 < n ∈ N, are continous and function P(t, z) ≥ 0, z ≥ 0, ∀t ≥ t 0 is continuous and Lipschitz in z. is enfolds the following Bernoulli equation which assumes analytical solutions: Clearly, z 4 (t, z 0 ) ≤ z(t, z 0 ), ∀t ≥ t 0 , and z 4 (t, z 0 ) should be close to z(t, z 0 ) if P(t, z) is sufficiently small. us, 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).

Simulations
is 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]. ese 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, lim t⟶∞ c(t) � ∞.
Next, let V � [v ij ], i, j � 1, ..., 4, be the eigenvector matrix for matrix A that is defined above, and in (33) As in Section 2, Finally, we set in (41), G − � G − Im(diag(G)) and As in Section 1.2, abs(q) � Let us show briefly how to derive listed above formula for L(t, z). Assume that Y i � 4 k�1 v ik y k , i � 2, 4. en, (79) can be written as follows: T . Subsequently, application of standard inequalities, e 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: 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 ‖x(t, x 0 )‖ in two consecutive time steps were numerically detected. Similar approach was used to approximate the threshold values in simulations of our scalar auxiliary equations. e running time in the latter computations practically does not depend upon dimension of the original system and is reduced further since ‖z(t, z 0 )‖ monotonically increases in z 0 for ∀t ≥ t 0 . In contrast, the running time required to estimate the boundary of trapping/stability regions increases with dimensions as m n , where m is the number of points that were taken to discretize a phase-space variable and n is the number of these coordinates. e results of simulation of the Van der Pol-like system with the indicated parameters are partly shown in Figure 1.  .7), respectively, which correspond to a coupled Van der Pol-like system. Note that in this case, x 0 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 F 0 � 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 Pollike system. In fact, simulations of (43) estimate the threshold value of z 0 , i.e., z, which, in turn, defines the ellipsoid imbedded in the stability region through application of the formula, z � ‖V − 1 x 0 ‖. 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.

Coupled System of Duffing-Like
Oscillators. For this system, (77) remains intact, but (78) should be altered as follows: which, in turn, requires some alteration in (80), i.e., v 2k ⟶ v 1k and v 4k ⟶ v 3k , 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 (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 F 01 � 0.01, F 02 � 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).
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≥t 0 ‖G − (t)‖ is a large value, and sup ∀t≥t 0 L(t, z) becomes sufficiently large for relevant values of z. e accuracy of our prior technique is practically nonsensitive to variations of sup ∀t≥t 0 ‖G − (t)‖ but decreases with growth of sup ∀t≥t 0 L(t, z) as well.
However, this prior technique also requires numerical simulation of w(t), p(t), and c(t) which is avoided in the current approach.

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 ⟶ ∞. e 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. e 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. e 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]. is 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.