A comparison of possible exponential polynomial approximations to get commensurate delays

The paper is aimed at a comparative simulation study on three prospective ideas how to approximate a general exponential polynomial by another one having all its exponents in the exp-function as integer multiples of some real number. This work is motivated by spectral properties of neutral time-delay systems (NTDS) and the contemporary state of the knowledge about the spectrum of NTDS with commensurate delays which are characterized by the latter family of exponential polynomials. The three ideas are, namely, those: Taylor series expansion, the interpolation in points given by dominant roots estimates and the special extrapolation technique presented by the authors recently. The goal is to match dominant parts of both the spectra as close as possible. However, some properties from the so called strong stability point of view can not be, in principle, preserved. The presented simulation example demonstrates the accuracy and efficiency of all the methods.


Introduction
The characteristic quasipolynomial of a linear-time invariant time delay system (TDS) gives the fundamental information about the system's dynamics; for instance, its zeros coincide (under some conditions) with system poles (i.e.eigenvalues) [1,2].Much effort has been made to analyse the infinite spectrum of TDS of retarded as well as neutral type (NTDS), mainly when the endeavour to decide about system stability and its dependence on particular delays' values [3][4][5][6][7][8].Compared to retarded ones, NTDS are much more advanced, tricky and intricate regarding spectral properties [1,4,9].Namely, positions of vertical strips of poles are sensitive to infinitesimal delay changes, which give rise to the notion of strong stability [10,11] that is affected i.a. by the rational dependence of delays [9].
These infinite vertical strips constitute the so called essential spectrum of a NTDS and they are unambiguously given by roots of the associated exponential polynomial [1].While there are many analytical results on the (essential and whole) spectrum of systems with commensurate delays in the literature, see e.g.[6,7], it is extremely difficult to find any exact laws for roots loci in the case of non-commensurate delays.
The goal of this work should be to establish a bridge between these two spaces by design of some possible ways how to approximate an exponential polynomial with non-commensurate delays by that with commensurate ones -which can be then analyzed easier.The accuracy of the approximation is measured by the matching of a subset of both the essential spectra with the most dominant poles.Roots loci are computed and displayed by means of the (advanced) Quasi-polynomial mapping Rootfinder (QPmR) [12] which is available as the aqpmr function in MATLAB.
We compare three methods in this contribution.As first, Taylor's series expansions of approximating and approximated exponential polynomial in the dominant root estimation point can be made, which is equivalent to the derivatives' equalities up to an appropriate order.As second, linear and quadratic extrapolation procedures based on the Taylor's series expansion again can be used.As third, the interpolation idea is another natural way how to cope with the problem.We present and compare two possibilities of the selection of points for the interpolation.
All these techniques require a sufficiently accurate estimation of dominant roots of the approximated exponential polynomial, which can be done e.g. by using our gridding algorithm [13].Another problem is the selection of the value of the base delay for the commensuracy; hence, we provide the reader with a possible solution and give the comparison with some other options.Last but not least, obtained approximated roots do not constitute complex pairs due to complexvalued coefficients.
A simulation example preformed in the MATLAB/Simulink environment provides the reader with the comparison case study and it clearly indicates that some of the above introduced ideas are good enough for the approximation in question.
Throughout the paper, & , 5 , 1 and 0 1 denote the sets of complex, real, integer and natural numbers, respectively.For ∈ s & , s Re denotes the real part of s.

NTDS spectrum
Let the transfer function ( ) ( ) of a NTDS have no common zeros and poles.Consider the corresponding monic characteristic quasipolynomial of a NTDS as ( ) ( ) Define the spectrum of the NTDS and its essential spectrum as respectively.Then it holds that [4,9]: • If there exists a nonzero pair ij d for some i, j, then • There exists a vertical chain of poles, Re sup : • The value of γ is not continuous with respect to Ĳ .
Moreover, let us define some other useful notions.
For commensurate delays, it is easy to see that if for any 1 ∈ k .

NTDS stability
Regarding exponential stability, this notion coincides with the finite-dimensional case, i.e. a NTDS is exponentially stable if 1 < α (including the limit of infinite vertical chains); however, other habitual types of stability (asymptotic, H , etc.) are much more complicated, see e.g.[6,14].Moreover, there exists a specific stability notion for NTDS -strong exponential stability -expressing that it remains 0 < γ under small delay perturbations; note that isolated roots with γ are not considered in this definition.The system is strongly (exponentially) stable if and only if 1 : Moreover it holds for NTDS with rationally independent delays that 0 > ξ implies strong instability [15].Hence, from this point of view, the rational independency is a rather stronger notion compared to the rational dependency.Notice that condition (4) does not include delay values.
Whereas the value of γ is not continuous with respect to delay values, this relation holds for the safe upper bound estimation on γ introduced e.g. in [12] as ( ) In fact, the value of c calculated from (5) gives the relevant information about the exponential stability of ( )

Problem formulation and possible solutions
As introduced above, there exist many results on the (essential) spectrum of NTDS with commensurate delays but only a few on non-commensurate ones.Therefore, our intention is to approximate with a general vector Ĳ by ( ) ( ) where max k means the degree of commensuracy.Let us denote the spectrum of (6) as , and, analogously, the corresponding values γˆ, ξ ˆ, c ˆ. Three possible ways how to cope with this task follow.

Extrapolation method
This idea stems from the discrete formulation of ( ) Assume the substitution and the shifting operator q corresponds to the variable 1 − z from the z-transform.Hence, one can write for some where ϑ means the non-integer fractional part.
Then the Taylor's series expansion of the factor ( ) can be made.In [16] we have presented analytic result for the linear ); both of them are used in this study yet not displayed herein the paper.
The dominant roots have the decisive impact to the dynamics.From one point of view, the rightmost poles coincide with the dominant ones.However, this is not the best conception for our case since there are infinitely many rightmost roots, see (3).Another idea is to measure the dominancy by the distance from the zero (i.e. the absolute value) which, however, is applicable to stable systems only.Thus, if & ∈ 0 s expresses the dominant root, then ( )

=
can also be expressed as where derivatives are calculated in a suitable point 0 s , for instance, again in the dominant root.
Whenever there exists a number max k of unknown coefficients (if linear algebraic equations (9) are independent) to get the unambiguous solution.

Interpolation method
Unknown coefficients of a prescribed ( ) can be determined by the interpolation with Let us choose two options.As first, consider the infinite chain of the rightmost roots (3).As second, try to take roots estimates with the minimum modulus (except for those in the chain with a high imaginary part), see the last paragraph of 3.1.
Due to (3), the exact imaginary parts can not be taken; otherwise, the set (10) would have no solution.Therefore, small deviations r i in the imaginary parts are made.

Discussion on methods
The crucial task is the option of the base delay 0 τ .From the z-transform point of view, it has the meaning of the sampling period, the recommended value interval of which for a second order finite-dimensional system reads for some suitable selected real 0 λ .
Another issue to be touched is that all the above methods result in complex-valued coefficients of the eventual exponential polynomial.This i.a.implies that roots are not symmetrical to the real axis.Thus, in the example, we benchmark the option to take 0 Re s rather than 0 s in ( 7)- (10).Assume an exponential term in ( 6), ( ) , and its linear approximation as where for some integer k .Since the linear extrapolation method is based on the Taylor's series expansion, it can easily be deduced that this method coincides with (9) whenever the approximating exponential polynomial is constructed as in (14), see details in [16].
Apparently, strong stability condition (4) is not affected only if ; however, it can be verified that it is satisfied only for 0 0 = s , i.e. 1 0 = z [16].Another problem is that there is the transition from non-commensuracy to commensuracy that yields rational dependency.Therefore, condition (4) becomes less strict than in the case of rational independency.
Last but not least, the approximation might be improved by the iterative use of the particular method via the re-calculation of the leading root estimation 0 , 0 A s s ≈ .
Simply, the discretization described in the subsection 3.1 is used to get a polynomial and it is followed by the use of the function to (8) to get the eventual A s , 0 .

Example
Consider the following exponential polynomial to be approximated ( ) ( ) We are going to present several numerical test below.Let the initial dominant root estimation simply be 0 0 = s .Since the imaginary part is "unknown" in this case, set the default value 9 .0 0 = τ .In the interpolation method, let us set 0 . Corresponding results without and with iterations are displayed in Figure 1 and Figure 2, respectively, in the form of the dominant subset of A Σ .In all figures of the paper, the legend is as follows: Roots of ( 15) are denoted by circles (ż); those of the ( ) obtained by using the linear and the quadratic extrapolation method are denoted by crosses (x) and plusses (+), respectively, A Σ from the Taylor's series expansion is denoted by squares (Ƒ), and the interpolation method with the vertical chain of dominant roots i s , 0 and the roots .The iterative procedure has resulted in almost exact estimation of these roots (but also some more not included in the original spectrum); however, the dominant rightmost roots of a Σ have not been approached at all (except for the last submethod).Note that results for the interpolation method with roots i s , 0 in a vertical chain are not displayed since a part of the spectrum is quite far in the right half-plane.
As the second experiment, we attempt to get better initial estimation by the simple linear extrapolation ( 14) where where almost all the (sub)methods have given a very close spectrum estimation; however, it has not been proved that the iterative procedure should have given better results in terms of our example.Notice that the linear extrapolation method and Taylor series expansion for ( 14) yield the same spectra, as mentioned above.given by ( 16), iterative procedure.in the complex plane is not as multifarious as in the previous experiment, positions of the roots give a very good estimation of e Σ .
Eventual forms of ( ) and some spectral and stability measures (see section 2) from the last two experiments for 3 .0 0 = τ are given to the reader in Table 1.In Table 1, method "I" and "I-R" means the linear extrapolation with complex-and real-valued coefficients, respectively; quadratic extrapolation method is denoted as "II"; the Taylor's series expansion is denoted as method "III"; and the interpolation methods with i s , 0 in a vertical chain and i s , 0 with the minimum modulus are denoted as "IV" and "V", respectively.Finally, the last one submethod giving real-valued coefficients is denoted as "V-R".The original .All the methods give almost identical spectral abscissa that is, surprisingly, better than the original one.This is, however, due to the identified dominant root in the selected subset of the complex plane located left from the exact spectral abscissa.Another interesting result is given by method "II": A relatively high value of ξ and a low value of c .This is becase of a higher value of max k .
In general, strong stability becomes worse under the approximation.

Conclusions
Several possible methods of the characteristic exponential polynomial approximation for systems with neutral delays have been proposed and compared by an example.The goal has been to get commensurate delays for a simpler spectral properties determination.We have i.a.found that the simple linear extrapolation procedure gives satisfactory results and almost keeps basic spectral and (strong exponential) stability measures compared to the original; moreover, it coincides with Taylor series expansion with the linear delay expansion.Both the approximated and approximating spectra are matched best via the quadratic extrapolation method and the interpolation method that uses the dominant roots with the minimum modulus.On the contrary, the interpolation in points given by the vertical strip of roots has given poor results.An idea how to obtain real-valued coefficients has also been given to the reader.
In the future research, the natural extension of this work may lie in a complex commensurate or a finitedimensional approximation of the whole characteristic quasipolynomial.
respectively, of the order l .The identity

Figure 2
Figure 2. e Σ , A Σ for 9 .0 0 = τ minimum modulus give A Σ denoted by triangles ( ∇ , Δ ), respectively.Figures1 and 2clearly display very poor results where the non-iterative procedure has yielded the acceptable estimation of roots about 77 that last two values are the root estimate in the vertical chain and that with the second minimum modulus, respectively.The corresponding z-values via (non-iterative and iterative procedures are displayed in Figures3 to 8. Clearly, the best one has been obtained for 3

Table 1 .
Resulting exponential polynomials and their measures.