Library
|
Your profile |
Software systems and computational methods
Reference:
Sklyar A.Y.
Numerical methods for finding the roots of polynomials with real and complex coefficients
// Software systems and computational methods.
2024. № 3.
P. 64-76.
DOI: 10.7256/2454-0714.2024.3.71103 EDN: KTJPCE URL: https://en.nbpublish.com/library_read_article.php?id=71103
Numerical methods for finding the roots of polynomials with real and complex coefficients
DOI: 10.7256/2454-0714.2024.3.71103EDN: KTJPCEReceived: 23-06-2024Published: 05-10-2024Abstract: The subject of the article is the consideration and analysis of a set of algorithms for numerically finding the roots of polynomials, primarily complex ones based on methods for searching for an approximate decomposition of the initial polynomials into multipliers. If the numerical finding of real roots usually does not cause difficulties, then a number of difficulties arise with finding complex roots. This article proposes a set of algorithms for sequentially finding multiple roots of polynomials with real roots, then real roots by highlighting intervals that potentially contain roots and obviously do not contain them, and then complex roots of polynomials. To find complex roots, an iterative approximation of the original polynomial by the product of a trinomial by a polynomial of a lesser degree is used, followed by the use of the tangent method in the complex domain in the vicinity of the roots of the resulting trinomial. To find the roots of a polynomial with complex coefficients, we propose a solution to an equivalent problem with real coefficients. The implementation of the tasks is carried out by step-by-step application of a set of algorithms. After each stage, a group of roots is allocated and the same problem is solved for a polynomial of lesser degree. The sequence of the proposed algorithms makes it possible to find all the real and complex roots of the polynomial. To find the roots of a polynomial with real coefficients, an algorithm is constructed that includes the following main steps: determining multiple roots with a corresponding decrease in the degree of the polynomial; allocating a range of roots; finding intervals that are guaranteed to contain roots and finding them, after their allocation, it remains to find only pairs of complex conjugate roots; iterative construction of trinomials that serve as an estimate of the values of such pairs with minimal the accuracy sufficient for their localization; the actual search for roots in the complex domain by the tangent method. The computational complexity of the proposed algorithms is polynomial and does not exceed the cube of the degree of the polynomial, which makes it possible to obtain a solution for almost any polynomials arising in real problems. The field of application, in addition to the polynomial equations themselves, is the problems of optimization, differential equations and optimal control that can be reduced to them. Keywords: polynomials, root finding, iterative methods, numerical methods, numerical algorithms, algebraic equation, conjugate complex roots, recursive algorithms, roots of polynomials, root localizationThis article is automatically translated. You can find original text of the article here. Introduction There are many tasks of a very different nature in which the definition of the roots of polynomials is required. Algebraic equations arise in the study of equilibrium states of complex thermodynamic and mechanical systems, they often appear in aerodynamics, in flight mechanics. For example, the speed of the fastest climb of an airplane is determined from an algebraic equation of the eighth degree. Algebraic equations also arise when performing various geometric calculations – determining the points of intersection and conjugation of curved contours, when designing smooth surfaces, well-streamlined bodies, and many other tasks. Methods for solving equations to the second degree have been known since the time of ancient Greece. In the XVI century, analytical expressions for polynomials of degree 3 (Cordano formula) and degree 4 were obtained by Scipio del Ferro, Tartaglia, Ferrari [3]. In the early 20s, Abel and then Galois in the 30s of the XIX century proved that such formulas for equations of the nth degree in the general case for any n ≥ 5 obviously cannot be found. Various methods are currently used for the numerical approximate solution of equations of higher degrees, such as the Lobachevsky method, the Hitchcock method, the Gorner scheme for dividing a polynomial into a binomial and a square trinomial [1, 4] and others. Another type of algorithms based on obtaining iterative formulas, by extracting simple and quadratic multipliers from polynomials, followed by comparing the entries of polynomials with the remainder when the roots are approximate, and without remainder when the values of the roots are exact, is considered in articles by Chye Yong Un and A.B. Shein [5,6,7]. A number of algorithms for solving polynomial problems and an extensive bibliography are given by G. P. Kutischev [8]. It is also worth noting the approaches for solving them in [9, 10, 11]. Note that the very existence of a large number of different methods generally indicates that there is not a single "completely satisfactory" one. 1. Removing multiple roots A polynomial of degree n has exactly n roots and is represented as In general, not all roots of x k are different. Finding the roots of a polynomial encounters certain difficulties in cases where it has multiple roots. Let there be a multiple root x* with multiplicity m, then P n(x) can be represented as Therefore, P n(x) and P'n(x) will have a common divisor (x-x*)m-1. Then, using Euclid's algorithm, we can find the largest divisor Q(x) of the polynomials P n(x) and P'n(x). And the original polynomial P n(x) can be represented as The polynomial H(x) will have the same roots as P n(x), but will not have multiple roots. Thus, the initial task is reduced to finding the roots of a polynomial whose roots are all different. 2. Determining the range of roots Consider the equation (a 0 will be assumed to be equal to 1) According to the well-known theorem [1,2], all roots x k (k = 1, 2, ..., n) of a polynomial in the complex plane lie in a ring r <| xk |< R, where At the same time, we should immediately note that this statement can be somewhat strengthened. Let's replace the variable x=y/t. In this case, the original equation with respect to y will take the form Given the nature of the substitution, you can write By changing the value of t, it is possible to reduce the estimate of the upper bound. Let the maximum value take the modulus of the coefficient at y m. Then the score for this coefficient will be no less than t* is located from Thus, it can be argued that the estimate lies between the values corresponding to the current value of t and t*. At the same time, when t changes, other coefficients also change. Thus, we get restrictions on the intersection of the lines of their change. If t r does not lie within the interval t and t*, then it can be ignored. From the rest, select the value closest to t. If there are none, then the optimal t is r=t*. If the signs of derivatives if they match, then further improvement is impossible. Optimal t r=t*. If the signs of the derivatives match, then we take m equal to k and repeat the procedure. Consider the derivative of the k coefficient at the point t Consider as an initial polynomial of the fourth degree t4-10t3+35t2-50t+24=(t-1)(t-2)(t-3)(t-4) The initial upper bound is 51, m=1, t=1, Points of intersection
And finally we get 3. Search for valid roots Let's consider separately the search for positive and negative roots. Let the initial polynomial have the form (the highest coefficient a 0, without reducing generality, can be assumed to be equal to 1). Let's introduce b k=(a k+|ak|)/2 and c k=(a k-|ak|)/2. For a k>0 b k=a k, c k=0 for a k≤0 b k=0. c k=-a k. In these designations The polynomials P+(x) and P-(x) for x>0 are continuous non-negative monotonically increasing functions. Using the introduced polynomials P+(x) and P-(x) and the upper bound of the values of the roots of R, we will look for roots in the range [A, B], where A=0, B=R. The lower bound can be refined to the value r* by obtaining the range [r*, R*] instead of the range [r, R], but this does not matter in principle. Let's take a closer look at the algorithm for finding roots. Calculate the values of P+(x) and P-(x) at the ends of the interval and call the Root function to calculate the roots with the parameters A, B, P+(A), P-(A) and P+(B), P-(B). The Root function returns either the found root (value greater than 0), or -1 if the polynomial P(x) has no roots in this interval. The algorithm of the Root function has the following form. 1. Calculate the values of the initial polynomial at the ends of the interval P(A)=P+(A)-P-(A) and P(B)= P+(B)-P-(B) 2. If P(A)=0, we return the found root A, if P(B)=0, we return the found root B, if P(A)P(B), we call the standard STROOT procedure for finding roots in a given interval [A, B] (for example, by dichotomy or chords) and return, the value found for it. 3. Calculate Q 1=P+(A)-P-(B) – the lower bound of the values of P(x) in the interval [A, B] and Q 2=P+(B)-P-(A) – the upper bound of the values of P(x) in the interval [A, B]. 4. If Q 1 Q 2≥0, then there are no roots in this interval and return the value -1. 5. Calculate C=(A+B)/2. Calling the Root function to calculate the roots with the parameters A, C, P+(A), P-(A) and P+(C), P-(C). 6. If the calculated value of x is >0, then we return the found root of x. 7. Calling the Root function to calculate the roots with the parameters C, B, P+(C), P-(C) and P+(B), P-(B). We return the found root x (if there are no roots, then -1). The above recursive algorithm allows you to find the root of a polynomial in the specified range or make sure that there are no valid roots inside it. After finding the root u, we proceed to search for the next root, replacing the original polynomial P n(x) with a polynomial of lesser degree P n-1(x), P n(x)=(x-u)P n-1(x). The polynomial P n-1(x) has all remaining roots P n(x). Repeating the above algorithm, we find all the positive roots. To find negative roots, replace the original polynomial P n(x) with the polynomial (-1)n P n(-x), the roots of which have the same values, but with the opposite sign. Having found all the positive roots of the second polynomial, we will thereby find all the real roots of the original polynomial. Since there are no other real roots, the remaining roots represent pairs of complex conjugate numbers and the remaining polynomial has an even degree. 4. Search for complex roots Let the polynomial P 2 n(x) with real coefficients have no real roots, then it is represented as Consider the representation of a polynomial of degree 2 n in the form To standardize calculations for all m, we introduce q-2=q-1=q n-1=q n≡0, q n-2=1. We will find the roots based on optimizing the approximation of the original polynomial Pn(x) by the polynomial Rn(x). Let's introduce the function For given values p 0, p 1, the coefficients q i are found from the minimization requirement F(p,q).
Or
For given values of q i, the coefficients p 0, p 1 are found from the minimization requirement F(p,q).
Or
To find the polynomial x 2+p 1 x+p 0, you can use the following iterative algorithm. 1. Set the initial value of the coefficients p 0, p 1. 2. Based on (4.5), we calculate the coefficients q i. The problem is reduced to solving a system of linear equations (SLA) for a 5-diagonal matrix. Given the type of matrix, the problem is reduced to a number of preparatory operations of difficulty O(n 2) and solving a SLOUGH of dimension 3×3. 3. The resulting set of coefficients q i is used to obtain adjusted values p 0, p 1 in accordance with (4.6). 4. We estimate the magnitude of the calculation error in accordance with (4.3). If the error is less than the specified threshold, then we finish the algorithm. Otherwise, we proceed to step 2 of the algorithm. As a result, we get the polynomial x 2 +p 1 x+ p 0 giving the first pair of roots. The coefficients of qi give the polynomial which can be used by the above algorithm to obtain the next pair of roots. Note that the convergence rate of this algorithm is low. To speed up, you can use the processing of previously obtained values. Let's assume that at the ith and i+1st steps we got the values p i,1, p i+1,1 and p i,0, p i+1,0 with residuals F(p,q) – Fa, Fb, respectively. Calculate the values p c 1=p i,1+2(p i+1,1-p i+1,1) and p c 0= p i+1,0+2(p i+1,0-p i,0) and the corresponding Fc value. If F c<F b, replace the points and the values in them a=b, b=c and repeat the calculation. If F c<F b, then by quadratic interpolation we find the point c*, taking it as a suboptimal value and return to the main algorithm. At the optimum point, both (4.5 and 4.6) must be performed simultaneously In matrix form, the aggregate (4.5, 4.6) is represented as At |A|≠0, the algorithm converges to the solution of the problem. In addition, from the point of view of practical implementation, the hybrid algorithm looks much more efficient. Several iterations of the basic algorithm lead us to the vicinity of the root. Considering that in the absence of multiple roots |P'n(x)| is guaranteed not to turn to 0, for P n(x)=0, a scheme using derivative values (Newton's method) or its more subtle modifications [12]) provides sufficiently fast convergence. Let the approximation s=x 2+p 1 x+p 0 be obtained. Its root will be The value of P n(z) in the neighborhood of z 0 is represented as This allows you to calculate the approximation of the z* root. At z 0, close to z*, the algorithm converges quite quickly. As a result, we get for s=x 2+p 1 x+p 0 The resulting expression is used to obtain the polynomial Q n-2(x) from P n(x)=(x 2+p 1 x+p 0)Q n-2(x), after which we look for the following roots, polynomials that coincide with the roots of Q n-2(x). 5. Search for roots of polynomials with complex coefficients Consider the equation (a 0 will be assumed to be equal to 1). The coefficients a k are complex numbers. Let its roots be x 1, x 2, ..., x n, then the roots of the equation will be Let's construct a polynomial R 2 n(x)=P n(x)Q n(x). This polynomial has roots
The quantities b k are representable as sums and products of complex conjugate numbers and, therefore, are real numbers. Thus, the polynomial R 2 n(x) is a polynomial with real coefficients. Accordingly, the search for the roots of the polynomial P n(x) with complex coefficients is reduced to the search for the roots of the polynomial R 2 n(x) with real coefficients. Note that all its real roots will be multiples, and some of the obtained complex conjugate roots will be extraneous. Conclusion The sequence of the proposed algorithms makes it possible to find all the real and complex roots of the polynomial. To find the roots of a polynomial of degree n with real coefficients, an algorithm is constructed that includes the following main steps: · definition of multiple roots, · highlighting a range of roots; · finding intervals guaranteed to contain roots (in time O(ln(R/H)), where R is an estimate of the total range of roots, and H is the length of the interval guaranteed to contain the root); · iterative construction of trinomials that serve as an estimate of the values of pairs of complex conjugate roots (the number of iterations does not necessarily have to meet the requirements for the accuracy of the solution, it is enough to fulfill the requirement of root isolation, after which the solution is achieved by traditional methods, for example, the Newton method). To find the roots of a polynomial of degree n with complex coefficients, an auxiliary polynomial of degree 2 n with real coefficients is constructed, containing all the roots of the original polynomial of degree n with complex coefficients, for which the above algorithm is used. The computational complexity of the algorithm is quite acceptable for calculating the roots of polynomials of almost any degree found in real problems. The software implementation of the proposed algorithms is performed in C++ in the Windows operating system environment. Given the computational nature of the algorithms, transferring them to a different operating environment does not cause any problems. References
1. Kurosh, A.G. (1968). Курс высшей алгебры [The course of higher algebra]. Moscow: Nauka.
2. Samarsky, A. A., & Gulin, A. V. (1989). Численные методы [Numerical methods]. Moscow: Nauka. 3. Stillwell, J. (1989). Mathematics and Its History. New York: Springer. 4. Tynkevich, M. A., & Pimonov, A. G. (2017). Введение в численный анализ [Introduction to Numerical Analysis]. Kemerovo: KuzGTU. 5. Chee, Yong Un, & Shein, A.B. (2012). Метод нахождения корней многочленов. I [Method of finding the roots of polynomials. I]. Информатика и системы управления, 4(34), 88-96. 6. Chee, Yong Un, & Shein, A.B. (2013). Метод нахождения корней многочленов. II [Method of finding the roots of polynomials. II]. Информатика и системы управления, 1(35), 108-118. 7. Chee, Yong Un, & Shein, A.B. (2013). Метод нахождения корней многочленов. III [Method of finding the roots of polynomials. III]. Информатика и системы управления, 3(37), 110-122. 8. Simon Telen. Polynomial Equations: Theory and Practice. Michal Kočvara; Bernard Mourrain; Cordian Riener. Polynomial Optimization, Moments, and Applications. Springer, pp. 215-240. 9. Simon Telen. Polynomial Equations: Theory and Practice. Michal Kočvara; Bernard Mourrain; Cordian Riener. Polynomial Optimization, Moments, and Applications. Springer, pp. 215-240. 10. B. Mourrain & J. P. Pavone. Subdivision methods for solving polynomial equations. Journal of Symbolic Computation, 44(3), 292-306, 2009. 11. Berthomieu, C. Eder, & M. Safey El Din. msolve: A library for solving polynomial systems. In Proceedings of the 2021 on International Symposium on Symbolic and Algebraic Computation, pages 51-58, 2021. 12. Statsenko, I. V. (2020). Исследование скорости сходимости одного обобщенного ньютоновского метода и классического метода ньютона в процедуре уточнения корней многочлена [Study of the Convergence Rate of a Generalized Newtonian Method and the Classical Newtonian Method in the Procedure for Refining the Roots of a Polynomial]. Точная наука, 78, 2-9.
First Peer Review
Peer reviewers' evaluations remain confidential and are not disclosed to the public. Only external reviews, authorized for publication by the article's author(s), are made public. Typically, these final reviews are conducted after the manuscript's revision. Adhering to our double-blind review policy, the reviewer's identity is kept confidential.
Second Peer Review
Peer reviewers' evaluations remain confidential and are not disclosed to the public. Only external reviews, authorized for publication by the article's author(s), are made public. Typically, these final reviews are conducted after the manuscript's revision. Adhering to our double-blind review policy, the reviewer's identity is kept confidential.
|