Evaluation of Computing Symmetrical Zolotarev Polynomials of the First Kind

This report summarize and compares various methods for computing the symmetrical Zolotarev Polynomial of the first kind and its spectrum with each other. Suitable criteria are suggested for the comparison. The best numerical stability shows the method employing Chebyshev polynomial recurrence. In case of the polynomial spectrum computation the best method is the one using the difference backward recursion introduced by M. Vlček. Both methods are able to generate the polynomial of high degree up to, at least, 2000, using 32-bit IEEE floating point arithmetics.


Introduction
Zolotarev polynomials (ZP) are important generalization of Chebyshev polynomials.The practical utilization of ZP was not possible till recently.The reason was the absence of method able to generate ZP of high enough degree, as was emphasized in [1].The ZP play, among others, the key role in analytical design of the optimal FIR filters [2], [3].Particularly, the symmetrical ZP of the first kind (ZP1S) are utilized in the design of almost equiripple half-band FIR filters [4].The utilization of ZP1S for nonstational signal analysis also seems to be promisable [5].Recently ZP were employed for design of circular antenna array [6].In spectral analysis the degree of the ZP1S is related to the number of samples of analyzed signal segment, and it is also directly related to the order of a filter.Therefore, utilizing ZP1S in both applications requires to generate polynomials with relatively high degree, typically, degree of hundreds.There are two typical tasks in computing ZP: One is the computation of actual ZP, while the ZP spectral coefficients computation.No method had been able to generate the spectrum of polynomial with satisfying high degree [2], until the recursive algorithm was developed [7].Nevertheless, based on findings in [7], and others, there are several methods of the computation of Chebyshev polynomials or ZP or their coefficients (spectra) e.g.[8][9][10][11][12].However, the usability boundaries of these methods and a detailed numerical comparison with each other have not yet been systematically explored.It is also unclear which degree of ZP1S or its spectrum can be achieved.
The goal of this article is to explore and compare selected methods for computation of ZP1S and its spectrum.In task of actual ZP1S computation we decided to use two often used methods: Chebyshev recursion [13] and power expansion with binomial coefficients [13].For computation of ZP1S spectral coefficients we employed backward recursion [7].We developed an empirical technique determining if the results is not a valid polynomial.Using this technique we distinguish between selected ZP1S computational methods in range of ZP1S parameters.We also explore dynamic range during the computation and its stability.The paper is organized as follows: an introduction is in Sec.1; the Chebyshev polynomials are briefly described in Sec.2; the symmetrical Zolotarev polynomial of the first kind is outlined in the Sec.3; Section 4 summarizes various methods of the polynomial computation, the methods are evaluated in Sec.5; a brief discussion and conclusion is in Sec.6 and 7, respectively.

Chebyshev Polynomials
Before describing the Zolotarev polynomials we briefly summarize the Chebyshev polynomials.The Chebyshev polynomial T n (x) of the first kind is defined by the relation [13] T n (x) = cos(nθ) when x = cos(θ). (1) Variable x on interval ⊂ [−1, 1] corresponds to variable θ on interval ⊂ [0, π].Nonlinear transformation of variable x onto variable θ converts T n (x) into function cos nθ as illustrated in Fig. 1a and Fig. 1b, respectively.
The Chebyshev polynomial T n (x) of degree n can be generated iteratively by following forward recursion [13] T  with initial condition Substituting for x = cos θ in (2) results in well-known trigonometric identity Nevertheless, this formula can not be used for computing ZP1S.
The solution of the recursion (2) results in a formula which describes the Chebyshev polynomial of the first kind of degree n in terms of power of x as [13] T where coefficients Thus, the Chebyshev polynomial T n (x) of order n can be also written as power expansion with coefficients given by binomial expansion.

Interpretation of Chebyshev Polynomial and its Spectrum
For the purpose of this study it is valuable to introduce the following interpretation of Chebyshev polynomial.If we regard the continuous-time variable x (or theta) as "time" then Chebyshev polynomial T n (x) (or cos(nθ)) can be regarded as a signal waveform.Chebyshev polynomial cos(nθ) can be converted to a spectral domain using Fourier series concept.Because we deal with an even symmetrical polynomial the Fourier transform degrades into cosine transform.The resulting spectrum of Chebyshev polynomial cos(nθ) of degree n is given where T 0 = π for this case.Using the principle of the orthogonality we obtain only one nonzero spectral coefficient on frequencies 2π T 0 .The Chebyshev polynomial and its corresponding spectrum is shown in Fig. 1b and Fig. 1c, respectively.This interpretation will be generalized later for Zolotarev polynomials.
It should be noted, it is possible to introduce the opposite interpretation, i.e. spectral coefficient a , = 0, 1, ..., n represent the impulse response of a filter in time domain and Chebyshev polynomial can be considered as the frequency response of the respective filter.This approach is used in filter design using ZP [3], [4], [7], [14].

Symmetrical Zolotarev Polynomial of the First Kind
The ZP1S can be expressed using the Chebyshev polynomial of the first kind as where p is the degree of Chebyshev polynomial as well as of the ZP1S, see [7].The ZP1S is defined on standard interval for polynomial approximation w ⊂ [−  The ZP are members of elliptical functions family of which mathematical discipline is rather difficult.To understand the ZP1S nature with ease we can approach to the Zolotarev polynomial as to trigonometric function.Employing this approach the ZP1S is expressed by the Chebyshev polynomial of the first kind: by using the substitution w = cos θ for ZP1S, denoted as Z p (w, k )|w = cos θ, we obtain function, which can be interpreted as a cosine with the central lobe elevated, see Fig. 3a.Similar approach is used also in [7] for the approximation of FIR filters by a ZP.For w = cos θ and k = 0 the Z p (w, k ) degenerates to the trigonometric cosine of the degree p.This can be easily justified: by substituting k = 0 in (6) we obtain Z p (w, k ) = (−1) p 2 T p (w), and by using the substitution w = cos θ we obtain cos pθ.
The global extreme of the ZP1S is influenced by both k and p: its absolute value gets larger when k or p get larger.This phenomena become obvious by looking at Fig. 3. Compare Fig. 3a with Fig. 3b, where the degree is the same and k is 0.10 or 0.11, respectively.Furthermore, compare Fig. 3a with Fig. 3c, where the k is the same and the degree is 16 and 20, respectively.Recent analysis shows the ZP1S are orthogonal in both disjoint intervals w ⊂ for fixed k and degrees p = 2, 4, 6, . . . .These findings are yet to be published.
The equation ( 6) can be rewritten using Use of this equation for future generation of ZP1S is beneficial compared with use of (6) for three reasons.Firstly, it lowers the upper index of the sum in equation (4a) to the half.Secondly, the number of iterations of recurrence (2) is also reduced to the half.Thirdly, the argument of T p 2 (x) is free of the square root witch is rather problematic computational operation.For all three reasons, the use of equation (7) enables to enumerate the ZP1S with less errors compared with equation (6).
It is also worthy to note that the ZP1S can be alternatively expressed as the linear combination of the Chebyshev polynomials T n (x) as [7], [15].
By using the substitution w = cos θ in (8) we obtain The a(2 ) coefficients2 are, in fact, the coefficients of Fourier series3 [18].Therefore, the a(2 ) coefficients can be regarded as the spectral coefficients of Z p (w, k ).Equation ( 9) expresses the spectrum of ZP which involves generally p 2 spectral lines.The expansion ( 8) can be generally used for the approximation of an arbitrary continuous function over a finite interval [−1, 1].One possible solution how to generate the coefficients a( ) can be found in e.g.[8].Another very effective and simple approach, in which we are interested, can derived from the linear differential equation for ZP [7] resulting in backward recursion with time-varying coefficients.

Discretization
This section introduces a discretization of continuoustime variables x and θ which is needed for computing ZP1S on a computer.All above defined formulas shall be used in a discrete-time manner requiring proper sampling of continues-time variables x or θ.The discretization can be done by letting x = iT, where i = 0, 1, 2, ... and T is the sample step for variable x.Similarly, the variable θ is discretized by letting θ = Ti, where i = 0, 1, 2, .... From this point forward the discussion will be supposing that all schemes of ZP1S are discrete-time functions.Let us remind that when a proper sample step is used then the behaviour of a discrete-time method is close to a behavior of original continuous-time one.We use sample step approximately

An Approach to FIR Design Using Zolotarev Polynomial
Almost equiripple half-band filter originates from ZP1S differential equation [4]; However, in this article we show simplified possible approach to FIR filter design employing ZP1S.ZP1S can form frequency response of a pass band FIR filter with passband around central frequency.A FIR filter transfer function can be written in general form [16] as 2Since the ZP1S is even symmetrical function only even coefficients a are nonzero.3In this case the Fourier cosine series.where h(v) is the impulse response, N = 2M + 1 is the filter length and M is the filter order.The filter transfer function employing ZP1S polynomial is given [7] as where A is the scaling coefficient.By comparing (10) with (11) the impulse response can assembled from a coefficients as Such an impulse response is even symmetrical and forms a linear phase type I. FIR filter.The coefficients a(v) are direct output of the backward recursion algorithm [7].The variable w in ( 11) is related to complex variable z by transformation w = 1 2 (z + z −1 ).The transformation degrades to w = cos(θ) in case of frequency response H (z = e jθ ): The k is a function of θ S which is symmetrical stopband edge.Having in mind the transformation of w in case of frequency response w = cos(θ) the k is simply given as where θ S1 is left stopband edge.The stopband ripple is gives as The scaling coefficient A is virtually maximum value of ZP1S.The maximal value is localized in the center of the polynomial and can be derived from ( 8) by substituting w = 0 as The derivation of degree equation p = f (δ, θ S ) is out of scope of this paper.The example of symmetrical pass band FIR filter using ZP1S is in Fig. 4. Also a high pass filter can be formed by discarding zero samples of filter impulse response.A low pass filter can be easily converted from a high pass one.Various filter types based on Zolotarev polynomials are described in [3], [4], [7], [14] and others.

Possible Schemes of ZP1S Computation
This section defines possible schemes of ZP1S computation.Any ZP1S of order p can be computed using two basic approaches.Firstly, direct computation in "time" domain: by using power expansion with binomial coefficients (4a) or by implementing the forward recursion (2).Both methods in "time" generate a Chebyshev polynomial of degree p along with the transformation of parameter x onto θ.The parameter x is transformed as in terms of (6) and as in terms of (7).
There are other possibilities of computing ZP1S, as mentioned before.For example, using the factorized form given by ZP1S zeros [7], or rather complicated methods using Remez algorithm or other approaches based on least deviation and polynomials with weights in "time" domain [10] and [9].But we are interested in very effective recursive method [7] promising the computation of ZP spectra of high degree and in its counterpart in the "time" domain (2).The aim is to compare the numerical behaviour of these selected methods.

ZP1S Computation in "Time" Domain
This subsection briefly discusses and compares the numerical behaviour of linear time variant (LTV) recursion ( 2) and the power expansion with binomial coefficients (4a), both used for computing a Chebyshev polynomial of degree p.

Computation Using Expansion of T n (x)
This method uses binomial coefficients given by equation (4a) to generate Chebyshev polynomials T n (x).Due to use of binomial coefficients the computation demands rather big dynamic range [17] and then poor computational precision can be expected.In consequence only the low degrees of ZP1S can be generated using this technique.

Computation by Using LTV Recursion
This method of computing Chebyshev polynomials based on using recursion (2) with the substitution (18).The discretized version of equation ( 2) can be written as with the initial condition where N is the number of points at which ZP1S is sampled.
Parameter N is given by used sample step T and the definition scope T 0 of variable x or θ: N = T 0 /T.In our case we use N = 10 4 .
Recursive schemes similar to (2) no matter if they are time-invariant or time-varying are known to be very effective and stable when used for limited number of samples even when fixed-point arithmetics is used [18][19][20].But to our knowledge, no systematic study of numerical behaviour of this type of LTV recursion has not yet been performed for the special case of computing of ZP1S of very high degree.If recursions (2) or (3) are used then the zeros of resulting Chebyshev polynomials are almost equidistantly spread in the intervals [−1, 1] or [0, π].In the latter case of T n (x) = cos(nθ) interval between zeros is precisely equidistant.However, if the Chebyshev recursion (2) is applied to compute ZP1S, especially for higher values of k or higher ZP1S degree, its zeros are pushed towards the edges of both the definition scopes [−1, 1] and [0, π].Hence, higher errors in computing ZP1S can be expected.The recursive schemes (2) or (3) are known to be very effective and stable when used for limited 4For discrete-time functions by discrete time cosine transform.number of samples.The scheme is stable even when fixedpoint arithmetics is used.Detailed numerical study of the recursive equation of this type can be found in e.g.[19].

Conclusion.
We briefly summarize comparison of methods of T n (x) computing using either the Chebyshev recursion (2) or using power expansion with binomial coefficients (4a).The latter method is expected to has worse numerical behavior.The evaluation of the binomial coefficients by equation (4b) causes severe problems due to high dynamic range.We can see in (4a) the values of coefficients c (n)  k can be up to twice greater than the values of binomial coefficients n−k k , which results in additional increasing of dynamical range.Nevertheless, there are some modified schemes of computing binomial coefficients (4a) which offers better numerical performance, e.g.[21]; however, we did not implemented these schemes.Due to described reasons we favor the use of Chebyshev recursion (2).

Generation of ZP1S Using Spectral Coefficients
The algorithm to compute ZP1S spectral coefficients is derived by Vlcek [7].The algorithm for computing ZP1S is just one of a family, of which algorithms are able to generate various kinds of ZP.The algorithm for ZP1S computation is a backward recursion with time-varying coefficients, and it can be described by where The result is set of spectral coefficients a , = 1, . . ., p.The spectral coefficients are transformed into ZP1S using the inverse cosine series5.
Tab. 1. Recursive evaluation of the coefficients a(2 ) for symmetrical Zolotarev polynomial of the

Evaluation of ZP1S Generation Methods Numerical Behavior
In this section we explore numerical behavior of the ZP1S computational methods described in previous section.Firstly we show that both LTV recursion (2) and ( 21) are stable.Secondly, we explore dynamic range of the methods.And at thirdly, we apply empirical criterion onto the methods results.
In case of recursion stability determination we employ following criteria.The first stability criterion states: if cumulated error throughout recurrence iterations is bounded by a linear function than the recursive system is considered stable.A recurrence system has to have a limited number of iterations.This approach is based on the one described in [19].As a second stability criterion we chose more analytical approach.We consider each iteration of LTV system as a LTI system.We show that for every iteration of LTV system every particular LTI system is stable.
The empirical criterion is based on fundamental theorem of algebra: the number of zeros of every polynomial is uniquely given by its degree.We facilitate this polynomial property by counting number of zeros in waveform generated by tested method.If a number of zeros differs from given polynomial degree than the generated waveform is definitely not the desired ZP1S polynomial.This criterion is not able to validate a method, but it is capable to point at one which is definitely unusable in a subset of ZP1S parameter space.

Stability of Chebyshev LTV Forward Recursion
In this section we show illustrative results confirming the stability of Chebyshev LTV forward recursion (2).We use previously drafted stability criterion: we adopt and slightly modify the approach suggested in [19].We evaluate accumulative sums of quantization errors of quantized recurrence.The computation using 32-bit floating-point IEEE number format is considered as "precise".We compare the computation of quantized recurrence in 32-bit fix point arithmetics with "precise" one.For the recursion (2) error is defined as where Tn [i] represents the quantized computation version, while T n [i] represents the "precise" one.The variable n is an index of recurrence iteration.The p stands for the given polynomial order and N for the number of samples taken from respective definition scope.The cost function given by the sum of errors is estimated according to The stability criterion is as follows: if the function J[n] is bounded by a liner function then the recursion in question is stable for limited number of iterations within the finite interval 1, p , where p < ∞.
The accumulative error J[n] for Chebyshev forward recurrence for k = 0.1 and degree p = 100 is shown in Fig. 5.The function J[n] is bounded, that is, any linear function with properly chosen slope gives greater values than J[n]; therefore, the recursion can be considered stable.We verified this criterion to be valid up to degree 2000 with various values of k .
Another interesting property of computed ZP1S waveform by the Chebyshev recursion method reveals the shape of cumulated error throughout θ, "time" dimension.This error which is than cumulated is given by The cumulative sum of error between "precise" computing T n [i] and quantized computing Tn [i] is shown in Fig. 6 for different polynomial degree within the interval [3,100].One can see that the error rapidly grows at the centre of definition scope, where the main lobe of ZP1S is placed.The error at the centre of definition scope gets bigger as polynomial degrees or values of k increases.For higher polynomial degree are errors the bigger as can be seen in Fig. 6 when we examine curves from bottom to up.This error increasing is the consequence of the central lobe maximum value and of ZP1S zeros pushing towards the edges of the definition scope.For higher polynomial degree and constant elliptical modulus k the main lobe of ZP1S is higher and wider pushing zeros more to the ends of definition scope as can be seen in Fig. 3a and 3c.

Stability of LTV Backward Recursion
The backward LTV recursion is given by ( 21) and Tab. 1.To analyze its stability we employ the first of described stability criteria.For recursion (21) we use the error of spectral coefficients given by where the â[ ] represents the quantized version with 32-bit fix-point arithmetics, and the a[ ] represents version using 32-bit IEEE floating point arithmetics, which is again considered as "precise".The p stands for the given polynomial degree and the n for iteration index of the recursion.The cost function is now given by the cumulative sum of errors (25) as The stability of backward LTV recursion ( 21) can be verified similarly to the stability of forward Chebyshev recursion (2).Again, if J[n] is a bounded (or linear) function of the index iteration n then the algorithm stability is ensured for limited number of iterations within the finite interval 3, p + 2 , where p < ∞.The shape of accumulated error J[n] (26) is shown in Fig. 7.This function is clearly bounded and has a similar shape as function J[l] in previous Fig. 5.The stability was experimentally proved for different polynomial degrees within the interval [10,2000].Furthermore, one can see that the recursion (21) gives less maximum of accumulated error than the Chebyshev forward recursion ( 2): compare values 6 × 10 −7 with 4 × 10 −6 .Also the slope of function J[l] at the beginning and at the end is greater than the slope of function J[l].Thus, it can be concluded, the backward recursion ( 21) is more robust with respect to stability and gives more precise results than the forward Chebyshev recursion (2).It is interesting to note that the stability of recursion improves with increasing number of iteration index n.This is illustrated by decreasing the slope (or "saturation") of J[n] with increasing iteration index n.It can be concluded that recursion (1) is stable.Note that the cumulative error of the backward recursion increases as the polynomial degree increases in linear, not exponential, rate.The second stability criterion examination follows.The analysis of this type of recursion can be performed using state-space equation and transition matrix (e.q.[20]).The dependency of coefficients d[ ] on both m = p 2 and parameters is nonlinear.The dynamical ranges corresponding to time-varying coefficients d[ ] are significantly different.For both these reasons we approach to LTV system is to set of LTI systems fixed for each iteration rather than to use the analysis of transition matrix behaviour.
Analyzed backward LTV recursion can be interpreted as IIR filter of order 6 with time-varying coefficients.When we trace parameters of the recursion in each iteration then the stability of this recursion during the whole iteration process can be verified using pole positions given by the denominator of the transfer function.The idea behind this approach is very simple.We substitute the LTV system by linear timeinvariant one.That means for each fixed iteration instant the system is supposed to be LTI with transfer function where d x [ ] are time varying coefficients defined in (21).Resulting denominator of the transfer function represents a polynomial for each iteration index .Roots of this polynomial, poles of transfer function H (z), are inspected and pole positions and the maxima of pole modulus are checked.This evaluation process is repeated throughout whole iteration process.One chosen result is illustrated for ZP1S parameters of k = 0.1 and p = 100 in Fig. 8, which illustrates the trace of pole positions during the whole iteration process.All pole modulus during the iteration process are less then one: hence, ensuring the stability.The values of poles modulus decreases with iteration index.
The best test results (the lowest values of poles modulus: the greatest distances of poles from the unit circle) are obtained for low values of k : k ≤ 0.2.Nevertheless, for all permissible values of k ⊂ [0, 1] the solution is stable.The stability was experimentally proved for different polynomial degrees within the interval [10,2000].

Dynamic Range
Dynamic range needed for all selected method for the computing ZP1S was checked for all polynomial degrees within the interval [10,2000].The result is as expected.Dynamic range needed for binomial coefficients computing (4a) is much greater than for the LTV forward Chebyshev recursion (2) and also greater than backward LTV recursion (21).
For example, for degree N = 100 and k = 0.1, the range of binomial coefficients c (n)  k in eq.(4b) values is approximately 10 37 , while for LTV forward recursion is about 10 6 and for LTV backward recursion and its coefficients a(2 ) approximately 10 5 .When we compare the dynamic range of the floating-point 32-bit IEEE number format with the range needed for the computing binomial coefficients we can see that the 32-bit floating-point number format is not sufficient for safe computing these coefficients.On the other hand dynamic ranges for both LTV recursions are comparable and much less than the dynamic range of the 32-bit floating-point number format.This is the reason, why both recursions work properly while binomial coefficients does not.

Empirical Test Criterion
The empirical criterion tests if a generated waveform representing ZP1S satisfies the number of zeros.The number of zeros is uniquely given by the polynomial order p.Therefore, we can find the number of zeros of computed polynomial ZP1S and compare it with the polynomial order p.If the number of found zeros differs from the given polynomial degree then the tested waveform ZP1S can not be a proper Zolotarev polynomial.The simple and robust enough way of enumerating the number of zeros is the zero crossings count (ZCC) where y[i] is tested discrete generated polynomial.Using (28) we can define the zero count error (ZCE) as Thus, if the ZCE computed by ( 29) is bigger then zero than the generated waveform do not meet ZP1S properties.If the number of zeros of computed polynomial ZP1S is equal to the given polynomial order then this ZP1S might be the proper polynomial.Tests for all selected methods of ZP1S computation are using 32-bit IEEE floating number format.All waveforms are represented by N = 10 4 samples.
The results of ZCE evaluated for polynomial generated using power expansion with binomial coefficients (4a) is in Fig. 9.We can see that the method (4a) definitely fails if the polynomial degree is greater than 80. Tested method Fig. 9.The ZCE evaluated for a waveform generated using power expansion with binomial coefficients (4a); non black areas denote the ZP1S parameter space where the generated waveform definitely does not meet ZP1S properties.
Fig. 10.The ZCE evaluated for a waveform generated using the LTV backward recursion followed by the iDCT; non black areas denote the ZP1S parameter space where the generated waveform definitely does not meet ZP1S properties.
computes waveform based on (7).If the waveform is computed based on (6) the method fails if degree reaches ≈ 35.In case a waveform generated by LTV Chebyshev forward recursion (2) the results show correct number of zeros throughout all values of tested ZP1S parameter space: k ⊂ (0, 1) and The results of the ZCE evaluated for a waveform computed using the LTV backward recursion (21) followed by the iDCT are in Fig. 10.We can see that it should be possible to generate ZP1S of high degrees; however the interval of k is being severally limited with increasing degree.

Discussion
The study of selected ZP1S computation methods reveals that the one employing the Chebyshev LTV forward recursion ( 2) is the most suitable method for computing ZP1S: the empirical criterion shows the number of zeros matching the polynomial degree throughout whole parameter space k ⊂ (0, 1) and p = 3, • • • , 2000; furthermore, the study confirms the stability of the LTV forward recursion in the same ZP1S parameter space.
In case of method using the LTV backward recursion followed by the iDCT the empirical criterion shows much worse performance, compared with forward recursion.However, the analysis of stability shows the LTV backward recursion more stable than the forward one.This seemingly contradicts, but it can be explained as follows.Both the LTV backward recursion and the iDCT are not free of the error propagation caused by quantization used in the numerical computation of ZP1S on a digital system.The problem is, we are unable to separate errors of the LTV backward recursion from errors caused by iDCT.Thus, one can expect that this method has worse results than direct computing ZP1S by forward LTV recursion (2).Since the LTV backward recursion is stable, even more stable than forward, it is obvious that the iDCT introduces greater errors during the computation compared with the LTV backward recursion.The iDCT errors are the result of the accumulation of cosines weighed by the spectral coefficients a(2 ) in (9).More detailed analysis of iDCT errors is out of the scope of this paper; however, it could prove interesting.
In case of employing of the ZP for design of FIR filters we have shown that the methods based on backward recursion are stable and capable of generate a filter coefficients of high orders up to thousands.However, corresponding filter frequency response, which is bounded with a filter coefficients by inverse Fourier series, might be problematic in some cases.Therefore, it seems to be advisable to use the LTV forward Chebyshev type of recursion to generate a filter frequency response.
We intend to focus future research on errors of iDCT, and generally on reducing of errors of latter methods for computing ZP1S.

Conclusion
Numerical study of computing symmetrical Zolotarev polynomials of the first kind shows that both recursions, the forward recursion for direct computing ZP1S and the backward recursion for computing spectrum of ZP1S, are stable.Method employing forward LTV Chebyshev recursion is the one most suitable one for computing the polynomial.For computing spectrum of the polynomial the Vlcek backward recursion has the best performance with respect to errors.This study also reveals that the method of computing ZP1S by backward recursion followed by iDCT has much worse numerical behavioral than the one using forward recursion.
Miroslav VLČEK was born in Prague, The Czech Republic, in 1951.He received the graduate degree in theoretical physics from Charles University, Prague, in 1974, and the Ph.D. degree in communication engineering and the D.Sc.degree from the Czech Technical University (CTU), Prague, in 1979 and 1994, respectively.From 1974 to 1993, he was with the Department of Circuit Theory, Faculty of Electric Engineering, CTU, where he is currently a Professor of the-ory of electrical engineering.Since 1995, he has been the Head of the Department of Applied Mathematics, Faculty of Transportation Sciences, CTU.He was the Alexander-von-Humboldt Fellow at the University of Erlangen, Nürberg, Germany, in 1988 and 1997.He currently teaches courses in system theory and digital filter design.His scientific interests include filter design and digital signal processing, and theory of approximation and higher transcendental functions.

Fig. 1 .
Fig. 1.An illustration of mapping of T n (x) onto cos nθ by transformation x = cos θ together with its spectrum.

Fig. 2 .
Fig. 2. The symmetrical Zolotarev polynomial of the first kind Z p (w, k ) for p = 8 and k = 0.1.

Fig. 8 .
Fig. 8. Pole positions of LTI transfer function 27 for all iterations of LTV backward recursion for given ZP1S parameters: k = 0.1 and degree p = 200.