IDENTIFICATION OF GENERIC STABLE DYNAMICAL SYSTEMS TAKING A NONLINEAR DIFFERENTIAL APPROACH

. Identifying new stable dynamical systems, such as generic stable mechanical or electrical control systems, requires questing for the desired sys- tems parameters that introduce such systems. In this paper, a systematic approach to construct generic stable dynamical systems is proposed. In fact, our approach is based on a simple identiﬁcation method in which we intervene directly with the dynamics of our system by considering a continuous 1-parameter family of system parameters, being parametrized by a positive real variable (cid:96) , and then identify the desired parameters that introduce a generic stable dy- namical system by analyzing the solutions of a special system of nonlinear functional-diﬀerential equations associated with the (cid:96) -varying parameters. We have also investigated the reliability and capability of our proposed approach. To illustrate the utility of our result and as some applications of the nonlin- ear diﬀerential approach proposed in this paper, we conclude with considering a class of coupled spring-mass-dashpot systems, as well as the compartmental systems – the latter of which provide a mathematical model for many complex biological and physical processes having several distinct but interacting phases.


(Communicated by Peter Hinow)
Abstract. Identifying new stable dynamical systems, such as generic stable mechanical or electrical control systems, requires questing for the desired systems parameters that introduce such systems. In this paper, a systematic approach to construct generic stable dynamical systems is proposed. In fact, our approach is based on a simple identification method in which we intervene directly with the dynamics of our system by considering a continuous 1-parameter family of system parameters, being parametrized by a positive real variable , and then identify the desired parameters that introduce a generic stable dynamical system by analyzing the solutions of a special system of nonlinear functional-differential equations associated with the -varying parameters. We have also investigated the reliability and capability of our proposed approach.
To illustrate the utility of our result and as some applications of the nonlinear differential approach proposed in this paper, we conclude with considering a class of coupled spring-mass-dashpot systems, as well as the compartmental systems -the latter of which provide a mathematical model for many complex biological and physical processes having several distinct but interacting phases.
1. Introduction. The question as to whether a system of differential equations has stable solutions is of vital importance in engineering where it occurs in the investigation of mechanical, control, and electrical systems (for more details, see e.g. [23], [24], [27] and the references given therein). It may be rephrased in the form of an equivalent question as to whether a polynomial with real coefficients has all its roots in the open left half-plane {z ∈ C : Re(z) < 0}.
Most differential equations and systems of differential equations one encounters in practice are in fact nonlinear. But in applications one is often interested in the behavior of a dynamical system near an equilibrium state. To study the behavior of a nonlinear dynamical system near an equilibrium point, we can linearize the system and deal with the derived linear system. Hence, it is of our interest here to know for which constant N -vector of real coefficients α = (α 1 , . . . , α N ), a system modeled by the following linear differential equation is stable z (N ) + α 1 z (N −1) + α 2 z (N −2) + · · · + α N −1ż + α N z = F(t); (1) or equivalently, for which constant vector α of coefficients, the characteristic polynomial s N + α 1 s N −1 + α 2 s N −2 + · · · + α N −1 s + α N (2) associated with the homogeneous differential equation corresponding to (1), in which we have set the input F(t) = 0, is a Hurwitz polynomial. The important point to note here is that the constant coefficients α i 's appearing in the differential equation (1) are in fact (given by) the parameters of the system whose dynamics is modeled by (1).
Having designed a dynamical system that works for any long period of time is perfect; however, at some point, one might decide to modify the system to work in a different but stable way. For instance, the control parameters for most of the flight controller need to be tuned, as inappropriate parameters will cause the flight be suboptimal or even unstable. Since the simulation results are different from real time experiments, usually, at first, a range for control parameters are obtained through simulation, and then the parameters are further tuned by hanging the vehicle on a test bench. But, due to the existence of some sensor noises which are hard to be filtered, the design and implementation of a test bench that can verify a six degree-of-freedom (DOF) flight control and simulate the real six DOF flight -in order to be a convenient platform for parameters identification of a (generic) stable system-is actually a hard task to do. In fact, several test benches have been recently designed and implemented but none of them are able to meet the requirements of the stable real flight simulation; for more details on this issue, see e.g. [1], [11], [15], [34], [36], [37].
In this paper, we wish to investigate as to how one can systematically modify a system being modeled by (1) by varying (some or all of) the system parameters in order to construct a stable system. This is precisely the situation that occurs in the designing of new stable mechanical or electrical control systems. To this end, following the original ideas of a recent work by Calogero [4] -instead of a test-bench approach as mentioned above, we take a nonlinear differential approach which is based on a simple identification method. In fact, we intervene directly with the dynamics of our system by considering a continuous 1-parameter family of system parameters, being parametrized by a positive real variable , and then identify the desired parameters that introduce a generic stable dynamical system by analyzing the solutions of a special system of nonlinear functional-differential equations associated with the -varying system parameters. One of the advantages of our proposed approach over the above-mentioned test-bench approach is in its reliability: taking our approach iteratively, one may achieve any desired accuracy in the parameters identification of our (generic) stable dynamical system, just by using the derived solution of each iteration as the initial value of a subsequent iteration.
The paper is organized as follows. Sect. 2 provides a setup for the required notation and terminology and presents some standard facts on stability. In Sect. 3, our main result Theorem 3.1 is formulated and proved. We conclude the third section with some facts and remarks on the reliability of the nonlinear differential approach proposed in the present paper, and show the utility of Theorem 3.1. To illustrate the capability of our proposed approach, we consider some examples in Sect. 4. Finally, some other interesting topics for future work in this area are given in the conclusion (Sect. 5) of the paper.
2. Background: The setup. In this section we set up some of the required notation and terminology, and compile some standard facts on stability.
We begin by introducing the notion of a Hurwitz polynomial. A polynomial with real coefficients is said to be a Hurwitz polynomial if all of whose roots lie in the open left half-plane {z ∈ C : Re(z) < 0}. Accordingly system that is modeled by the differential equation (1), regardless of its input function F(t), is called a stable system if its corresponding characteristic polynomial (2) is Hurwitz.
Let us denote by the monic characteristic polynomial of degree N ≥ 2, as in (2), where s is the complex variable, components of the N -vector α = (α 1 , . . . , α N ) are the coefficients of the polynomial (3), and x = {x 1 , . . . , x N } is the unordered array of the roots of polynomial (3). We say that a system being modeled by (1), or its corresponding characteristic polynomial (3), is generic if the coefficients and roots are generic complex numbers and the roots are simple (i.e., the roots are all different). In fact, (3) gives two representations of the same characteristic polynomial p N (s; α, x), as this monic polynomial can be identified by assigning either its N coefficients or its N roots; indeed, any coefficient α k can be expressed in terms of the N roots x n via the following standard formula There is a simple necessary condition for Hurwitz polynomials, being referred to as the coefficients criterion, which states: the coefficients of a monic Hurwitz polynomial are all positive. It is actually a simple matter to check this. Let p(s) be a monic polynomial of a given degree N ≥ 2, and suppose that x = {z 1 , . . . , z a ,z 1 , . . . ,z a , y 2a+1 , . . . , y N } is the set of its roots, where z i 's are the complex roots with non-zero imaginary parts,z i 's are the complex conjugate of z i 's, and y j 's are the real roots of the monic polynomial p(s). Hence whose roots must read Re(z i ) < 0 and y j < 0 since p(s) is Hurwitz. From the right-hand side, one can see that the coefficients are obtained by multiplications and additions of positive numbers; thus they are all positive.
2.1. Routh-Hurwitz criteria. One of the most important criteria that gives necessary and sufficient condition for the characteristic polynomial (with real coefficients) to be a Hurwitz polynomial is known as the Routh-Hurwitz criterion.
According to a standard fact, if the characteristic polynomial (2) is Hurwitz, then any solution to the homogeneous differential equation associated with (1) -in which we have set the input function F(t) = 0-converges to zero; this just amounts to saying that the corresponding homogeneous system models an exponential decay. In fact, Routh-Hurwitz criteria for differential equations are analogous to the Jury conditions for difference equations. The Routh-Hurwitz criterion can be formulated as follows.
Proposition 1 (The Routh-Hurwitz criterion). Given the characteristic polynomial where the coefficients α i are all real constants, define the N Hurwitz matrices using the coefficients α i of the polynomial as where α j = 0 for j > N . All of the roots of the polynomial p(s) are negative or have negative real part (i.e., p(s) is Hurwitz) if and only if the determinants of all the Hurwitz matrices are positive: For a proof of the Routh-Hurwitz criterion, we refer the reader to the classical work of Gantmacher [10].
For a polynomial to be a Hurwitz polynomial, Liénard and Chipart [20] established necessary and sufficient conditions which require about half the amount of computations required in the Routh-Hurwitz criterion (cf. Proposition 1 above).
Proposition 2 (Liénard-Chipart). Let the monic characteristic polynomial p(s) of degree N ≥ 2, and the Hurwitz determinants ∆ i be as in Proposition 1. Then the following statements are equivalent.

Monic characteristic polynomials.
In the spirit of a paper by Hardy [12] on the roots of a class of integral functions, Hutchinson [17] formulated an elegant necessary and sufficient condition for a polynomial to belong to a remarkable proper subset H of the set of Hurwitz polynomials, all of whose members have just simple negative real roots. In fact, Hutchinson's derived result is as follows.
Proposition 3 (Hardy-Hutchinson). Let p(s) be any monic characteristic polynomial of degree N ≥ 2 as in Proposition 1. The polynomial p(s) belongs to the class H of those Hurwitz polynomials having just simple negative real roots if and only if the coefficients α i of p(s), setting α 0 = 1, read the following inequalities

IDENTIFICATION OF GENERIC STABLE DYNAMICAL SYSTEMS 4545
The constraints (5) will be referred to as the Hardy-Hutchinson criterion. Assuming all the coefficients α i of the monic polynomial p(s) to be positive, one can readily deduce from (5) that Orlando's formula. In his work on the Routh-Hurwitz problem, Orlando [28] obtained a nice and simple formula which expresses the Hurwitz determinant ∆ N −1 in terms of the roots x 1 , x 2 , . . . , x N of the monic characteristic polynomial p(s) (of degree N ≥ 2 as in Proposition 1): When N = 2, this formula reduces to the well-known formula for the coefficient α 1 in the quadratic polynomial . . x N , we obtain from (7) what will be referred to as the Orlando's formula: Using the Orlando's formula (8) corresponding to a monic characteristic polynomial of degree N ≥ 2, we introduce a hypersurface sitting in R N (to be used in Sect. 3) as follows. Let u = (u 1 , . . . , u N ) be any vector in R N , and ρ = {ρ 1 , . . . , ρ N } be an arbitrary set of numbers that can be thought of as the roots of the monic polynomial under consideration. We can then define the hypersurface Π u ρ to be the locus of the set of points q = (q 1 , · · · , q N ) ∈ R N , with q k < u k for any k ∈ {1, · · · , N }, that satisfy the following relation of Orlando type We call Π u ρ the Orlando's hypersurface associated with u and ρ.
3. The result: Identification of generic stable systems. To formulate our main result, we introduce an additional independent continuous variable , which is assumed to be a positive real parameter. Now, similar to (3), we define two continuous 1-parameter families of auxiliary monic polynomials with -varying coefficients as follows which, for any value • of our parameter, can be seen as the monic characteristic polynomial of a system being modeled by a differential equation analogous to (1) with the coefficients vector α( • ) = (α 1 ( • ), . . . , α N ( • )), and where the components ∆ k ( ) of the N -vector ∆( ) = (∆ 1 ( ), . . . , ∆ N ( )) are the Hurwitz determinants associated with the monic characteristic polynomial p N associated with the parameter , as in (10), and r( ) = {r 1 ( ), . . . , r N ( )} is the set of roots of the monic polynomial ϕ N corresponding to , as in (11). As mentioned in Sect. 1, the coefficients appearing in the differential equation (1) are in fact (given by) the parameters of the system whose dynamics is modeled by (1); hence, in practice, there should always be some (say, upper) bounds for the coefficients α k in (1). Let us suppose that the components of the N -vector u = (u 1 , . . . , u N ) ∈ R N give the upper bounds, namely that cf. the relations (6) which give a necessary condition for a monic polynomial with positive coefficients to belong to the class H, as in Proposition 3, of those Hurwitz polynomials having just simple negative real roots.
We can now formulate our main result.
Proof. This theorem provides a solution to the problem of questing for constant coefficients α k that make a system being modeled by (1) a stable system, or equivalently, make the corresponding monic characteristic polynomial with -varying coefficients, as in (10) or (11), at some point, a Hurwitz polynomial. Taking the -derivative of the auxiliary monic polynomial ϕ N (s; ∆( ), r( )), as in (11), yields setting in which s = r m ( ) for m = 1, 2, . . . , N gives our system (13) In fact, the above system of equations is a system of nonlinear functional-differential equations in terms of the roots r m ( ) of the monic polynomial ϕ N (s; ∆( ), r( )), as some member of a continuous 1-parameter family of such polynomials, whose coefficients ∆ k ( ) are the Hurwitz determinants associated with the monic polynomial p N (s; α( )). Since determinant is a multi-linear function that sends an ordered array of columns of a square matrix to a real number, letting H k,i ( ) denote the i th column of the Hurwitz matrix H k ( ), one can evaluate the -derivative∆ k ( ) of the Hurwitz determinant as followṡ If for some λ ∈ (0, T ) the solution r( ) = {r m ( ); m = 1, . . . , N } of (13) reads Re (r m (λ)) < 0 for all m, then the characteristic polynomial ϕ N (s; ∆( ), r( )) at = λ is by definition a Hurwitz polynomial. From this we deduce that a system being modeled by (1) with the components of the N -vector ∆(λ) = (∆ 1 (λ), . . . , ∆ N (λ)) as its coefficients introduces a generic stable system. In addition, since ϕ N (s; ∆(λ), r(λ)) is a Hurwitz polynomial, the coefficients criterion implies ∆ k (λ) > 0 for any k ∈ {1, . . . , N }. By the Routh-Hurwitz criterion (see Proposition 1), this just amounts to saying that at = λ the polynomial p N (s; α( )) is Hurwitz, which is equivalent to the fact that the system being modeled by (1) and having its coefficients given by the components of the constant vector α(λ) = (α 1 (λ), . . . , α N (λ)) is a stable system. This finishes the proof.
We put off, to some concluding remarks below, discussing both the question as to whether this procedure always works and the utilities of Theorem 3.1.
It is worth pointing out here that we can attain an even better control over the parameters identification procedure to identify our (generic) stable systems, if we appropriately modify the proposed algorithm in Theorem 3.1 as follows (cf. the concluding remark of [5] as a commentary written by Calogero on [4]).

Remark 1.
Let L > 0 be as above and f k : [0, L] → R be arbitrary differentiable functions with f k (L) − f k (0) = 0 for any k ∈ {1, . . . , N }. Then we may modify the assignment (14) for the components of the N -vector α( ) as Applying this modification, we get a new version of our algorithm which works exactly as before but it has some advantages over the proposed procedure in Theorem 3.1; e.g., employing the modified version, we may deal with each component α k of the N -vector α independently and thus obtain an even better control over the whole mechanism of our approach.
By the following remark we give an affirmative answer to the question as to whether this procedure always works.

Remark 2.
The only possible mechanism whereby the solution r( ) of the dynamical system (13) might run into a singularity during its evolution is exactly when two different r m ( ) coincide, i.e., r m1 ( ) = r m2 ( ) at some point , causing the righthand side of (13) to blow up (cf. Sect. II of [4]). This collision might indeed happen, however it is not a generic phenomenon, but it can be avoided by carefully making the assignment of the initial data α(0); in fact, for this purpose, it is advisable to include in the set of initial values r(0) the complex conjugate of any chosen complex member, with non-zero imaginary part, of the set r(0), and to have at least one pair of such complex conjugates in r(0); then it suffices to let α(0) be any point on the Orlando's hypersurface Π u r(0) as above. By Proposition 3, the Hardy-Hutchinson criterion gives a necessary and sufficient condition for a polynomial to belong to a remarkable class H of those Hurwitz polynomials all of whose members have just simple negative real roots; hence, the class H of Hurwitz polynomials given by the Hardy-Hutchinson criterion is, in fact, a subclass of the generic Hurwitz polynomials which are the characteristic polynomials of the generic stable systems given by Theorem 3.1.
One (not recommended) way of finding the desired coefficients α k (and hence the desired parameters) of our system being modeled by (1) to get a generic stable system, instead of applying Theorem 3.1, is to derive them by inspection: employing either of the Routh-Hurwitz criterion (or equivalently, its refinement, the Liénard-Chipart criterion; see Propositions 1, 2), or the Hardy-Hutchinson criterion (see Proposition 3). But who knows where in the continuous forest of real numbers the desired real parameters live, especially when we are dealing with systems of higher orders. In fact, as the order N of the monic characteristic polynomial p N corresponding to our desired systems increases, the randomly constructed stable systems become more and more scarce, because the probability of randomly constructing a stable system is bounded above by ( 1 2 ) N , due to the coefficients criterion which is a necessary condition for p N to be Hurwitz.
We conclude with some important points on the differential approach that we have proposed in Theorem 3.1, comparing it with a "Newtonian differential approach".
Remark 3. To achieve any desired accuracy, one may use the method of our proposed approach of Theorem 3.1 iteratively, using the derived solution r(λ) = Figure 1. Schematic of a coupled spring-mass-dashpot system {r 1 (λ), . . . , r N (λ)} of each iteration as the initial value r(0) of a subsequent iteration (cf. Sect. II of [5]). In fact, following our differential approach, we can identify with any desired accuracy the system parameters of a generic stable system; i.e., a system whose characteristic polynomial has N simple roots r k (λ) with negative real parts -from the practical point of view, it is important to note here that to the dynamics of such a system, ignoring the input, there corresponds a homogeneous differential equation of the same order N which has N independent solutions of the following forms: exp(r k (λ)t), if r k (λ) is real, or exp(Re(r k (λ))t) cos(Im(r k (λ))t) and exp(Re(r k (λ))t) sin(Im(r k (λ))t) otherwise; i.e., N independent solutions in the form of exponential decays with bounded functional multiples.
In fact, one of the current well-known classical root-finding methods with a very high rate of convergence is the Newton-Raphson (also known as Newton's) method. It is of course possible to use Newton's method as part of a continuation algorithm to propose some other differential approach to derive the desired parameters in the construction of a stable system. But one of the advantages of our proposed approach (cf. Theorem 3.1) over the above-mentioned "Newtonian differential approach" is in its convergence: the differential approach based on the Newton's method may not even converge to the roots if, at least, one of the roots is very close to some stationary points. However, comparison of the actual effectiveness of our proposed approach with any other differential method (being essentially of the same class as mentioned above) is beyond the scope of this relatively short paper; but, this issue could be of interest to and considered by specialists in numerical analysis.
The above-mentioned facts and remarks show the reliability and capability of the main result of this paper, Theorem 3.1.
We shall illustrate, in Sect. 4, the utility of our proposed nonlinear differential approach in the identification of the system parameters for constructing some generic stable dynamical systems.
Applying Newton's law of motion, it can be easily seen that the equations of motion for this system satisfy the following system of differential equations Deriving x 2 from equation (17), and then substituting x 2 and its derivatives into equation (18), we obtain the following fourth-order ODE for x 1 which models the dynamics of our coupled spring-mass-dashpot system, where Now, using our nonlinear differential approach proposed in Theorem 3.1, we intend to identify the system parameters, {m 1 , m 2 , k 1 , k 2 , b}, so that the resulting system modeled by (19) becomes a generic stable system. To do so, for simplicity, we let L = 1, f ( ) = sin( ), and set the upper-bounds for the coefficients u = (u 1 , u 2 , u 3 , u 4 ) := (10, 10, 10, 10), and assume the initial values r(0) = {r 1 (0), r 2 (0), r 3 (0), r 4 (0)} of the system of equations (13) to be as follows. Then, assigning α(0) = (α 1 (0), . . . , α 4 (0)) := (1.1666764676, 4, 3.5, 3) which is a point on the Orlando's hypersurface Π u r(0) as above, and running the integration of the system of equations (13) using a computer algebra system (like Maple), we arrive at the following polynomial s 4 + 0.8989502523s 3 + 1.265740825s 2 + 0.7236743319s + 0.09921684031 (20) that can be considered as the characteristic polynomial of some fourth-order generic stable systems, whose roots are as follows In order to identify the system parameters for constructing a generic stable coupled spring-mass-dashpot system, we put the characteristic polynomial associated with (19) to be identically equal to (20); hence, the desired system parameters are just given by the solutions of some system of algebraic equations. Accordingly, the desired system parameters read m 1 = 0.1232474328, m 2 = 0.337471835, Needless to say that changing the initial values r(0) of the system of equations (13) and accordingly taking α(0) from the new Orlando's hypersurface in our nonlinear differential approach, we may arrive at a new set of desired system-parameters, constructing another generic stable coupled spring-mass-dashpot system. Compartmental systems. Many complex biological and physical processes can be subdivided into several distinct but interacting phases, each of which being called a compartment in the overall process, where the modeling system is referred to as a compartmental system (also known as, a block diagram, usually among engineers); cf. [35]. To illustrate the capability of our nonlinear differential approach, here, we shall consider an open compartmental system, with 12 compartments, in which there exist an input to as well as an output from the system. Let Y 1 , Y 2 , . . . , Y 12 denote the 12 separate interacting compartments of our system, and y 1 (t), y 2 (t), . . . , y 12 (t) be respectively the concentration of a homogeneous substance in the compartments at time t. In addition, suppose that b in , b 1 , b 1 , b 2 , b 2 , . . ., b 11 , b 11 , and b out are the rates per volume at which the fluids containing the substance flow into and out of the compartments of the open compartmental system shown in Fig. 2. Also let c denote the constant concentration of the substance in the fluid flowing into compartment Y 1 .
Usually it is assumed that the volume of each compartment remains constant, and hence the sum of the input rates per volume into a compartment must be equal to the sum of the output rates per volume from the compartment. This amounts to say that the rates per volume must read So the system parameters are just b in , b 1 , b 2 , . . . , b 11 . Under the assumptions that the volume of each compartment remains constant, the time for material to flow from one compartment to another is negligible, and the fluid in each compartment has uniform concentration, the concentrations y i for i = 1, . . . , 12, satisfy the following system of differential equationṡ Taking the similar steps as done for the previous example in order to derive a single dynamical equation from the system of differential equations of motion, coupling (21) with the system of equations (22), (23), (24) we can get the following twelfth-order ODE for y 1 , which models the dynamics of our compartmental system, y and the rest of coefficients in (25) being just some lengthy algebraic functions in the constants b in , b 1 , b 2 , . . . , b 11 . Using our nonlinear differential approach proposed in Theorem 3.1, we intend now to identify the system parameters, {b in , b 1 , b 2 , . . . , b 11 }, so that the resulting system modeled by (25) becomes a generic stable system. To do so, for simplicity, we let L = 1, f ( ) = sin( ), and set the upper-bounds for the coefficients u = (u 1 , u 2 , . . . , u 12 ) with u i := 10 for all i, and assume the initial values r(0) = {r 1 (0), r 2 (0), . . . , r 12 (0)} of the system of equations (13) to be as follows. and running the integration of the system of equations (13) using our computer algebra system, we arrive at the following polynomial s 12 + α 1 s 11 + α 2 s 10 + · · · + α 11 s + α 12 , that can be considered as the characteristic polynomial of some twelfth-order generic stable systems, whose roots are as follows It is worth noting that the long-term behavior of our generic stable compartmental system being modeled by (25), with the above system parameters, can be given (independent of the initial conditions for the system) by the constant function F = cb in b 1 b 2 b 3 · · · b 11 for any positive real c.

5.
Conclusion. In this paper, we have proposed a nonlinear differential approach for identification of (generic) stable dynamical systems. This approach is based on a simple identification method in which we intervene directly with the dynamics of our system by considering a continuous 1-parameter family of system parameters, being parametrized by a positive real variable , and then identify the desired parameters that introduce a generic stable dynamical system by analyzing the solutions of a special system of nonlinear functional-differential equations associated with the -varying parameters. We have also indicated some of the advantages of our approach over both the test-bench approach (being widely used in the parameters identification of aerial vehicles/robots) and the classical Newtonian differential approach (as discussed above).
Due to the simple mechanism and exceptional agility of quadrotors 1 , many theoretical researches have been conducted on developing the control system for these aerial vehicles, which provided an ideal field for performing advanced tasks in the area of aerial robotics; for more details, see [3], [13], [21], [22], [25], and the references therein. Since quadrotors use fix-pitch propellers, the rotational inertia of the motor is a critical factor for the stability of the vehicle. The actuator dynamics of the motor-propeller can be simplified as a first order linear system. It has been shown in [26,30] that if the inertia of motor-propeller for a quadrotor increases and exceeds a certain value in its range, then the vehicle will become unstable. However, to the best knowledge of the author, there exists actually no reliable method to resolve this stability issue. Application of our proposed approach to resolve the stability issue, concerning the actuator dynamics of quadrotors, and possible modifications of this approach to carry out the parameters identification for some stable complex aerial vehicles will be the subject of future publications.