Lirac<[email protected]> asked a question about mata function polyroot():
>I am running Stata 10.1 on Mac OS 10.4. The following issue is perplexing me:
>
>-------
>mata:
>mata clear
>polynomial1 = (1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-30,29)
>polyeval(polynomial1,1)
>polyroots(polynomial1)
>-------
>
>The polyeval line makes it clear that x=1 is a root for this polynomial.
>However, 1 is not a reported result for polyroots.
The first two results reported by polyroots(polynomial1) are
(1-1.1729e-08i, 1+1.1729e-08i). The results are correct.
Notice polynomial1 can be factored into
(x-1)*(29*x^29 - x^28 - x^27 - x^26 -... - 1),
and the second term is 0 at x=1. Hence polynomial has a factor (x-1)^2,
i.e., 1 is a double root of polynomial1. Mata correctly identifies this.
Due to numerical round-off, Mata computes the degree two factor as
(x-1)^2+e, where e is a small positive number, instead of
exact (x-1)^2, and consequently Mata reports two complex conjugate roots
close to 1 instead of two 1.
>What is odd is that
>if I add or subtract a zero term to the polynomial, polyroots() will
>then correctly report x=1 as a root. i.e., the following works:
>
>polynomial2 = (1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-30,29)
>polynomial3 = (1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-30,>29)
>
>Does anybody know what is going on here?
>
Both polynomial2 and polynomial3 have only (x-1) factor. Mata correctly
reports 1 as a simple root.
--Hua
([email protected])
*
* For searches and help try:
* http://www.stata.com/help.cgi?search
* http://www.stata.com/support/statalist/faq
* http://www.ats.ucla.edu/stat/stata/
