I want to find the real roots of a polynomial that are 0<=roots<=1, for example:
>> rts=roots([16 -32 116 0 -100])
rts =
0.9050 + 2.6262i
0.9050 - 2.6262i
1.0000 + 0.0000i
-0.8100 + 0.0000i
I should be able to find that rts(3) satisfies these requirements. However,
>> myrt=rts(logical(~imag(rts) & real(rts)>=0 & real(rts)<=1))
myrt =
Empty matrix: 0-by-1
This is because real(rts(3))<=1 yields false, which is a bug so far as I can tell since complex(1,0)<=1 yields true. If I open the variable rts and double click rts(3) (which shows that it really is stored just like complex(1,0)), I can then execute real(rts(3))<=1 and it produces the correct result. What gives?
How do I find that rts(3) is the correct answer when I write a script? I am using MATLAB R2015b.

 Réponse acceptée

James Tursa
James Tursa le 21 Nov 2016
Modifié(e) : James Tursa le 21 Nov 2016

0 votes

R2016b:
>> rts=roots([16 -32 116 0 -100])
rts =
0.904997916303652 + 2.626226923351172i
0.904997916303652 - 2.626226923351172i
1.000000000000000 + 0.000000000000000i
-0.809995832607303 + 0.000000000000000i
>> rts(logical(~imag(rts) & real(rts)>=0 & real(rts)<=1))
ans =
1.000000000000000
>> num2hex(rts(3))
ans =
3feffffffffffffc
>> polyval([16 -32 116 0 -100],rts(3))
ans =
-8.526512829121202e-14
>> polyval([16 -32 116 0 -100],1)
ans =
0
R2011a:
>> rts=roots([16 -32 116 0 -100])
rts =
0.904997916303650 + 2.626226923351169i
0.904997916303650 - 2.626226923351169i
1.000000000000000
-0.809995832607303
>> rts(logical(~imag(rts) & real(rts)>=0 & real(rts)<=1))
ans =
Empty matrix: 0-by-1
>> num2strexact(rts(3))
ans =
1.000000000000000444089209850062616169452667236328125
>> num2hex(rts(3))
ans =
3ff0000000000002
>> polyval([16 -32 116 0 -100],rts(3))
ans =
8.526512829121202e-014
>> polyval([16 -32 116 0 -100],1)
ans =
0
So it looks like the roots function (which uses eig in the background) has different behavior depending on the version. The versions missed being exactly 1.0 by 2-3 least significant bits (one version less than and the other version greater than). You will have to adjust your code to account for this.

2 commentaires

Bryce Inman
Bryce Inman le 22 Nov 2016
Thank you for the reply! It's frustrating how misleading matlab's 'long' format can be. And I still wonder, why would double clicking on the variable cause the equality to evaluate correctly?
James Tursa
James Tursa le 22 Nov 2016
Not sure what you are saying here. How do you get real(rts(3))<=1 to evaluate as true? You double clicked on rts in the Workspace list and then that brings up the Variable Editor. Then what did you do?

Connectez-vous pour commenter.

Plus de réponses (1)

Walter Roberson
Walter Roberson le 21 Nov 2016

1 vote

You will find that real(rts(3))-1 is not 0.
When I test on my Mac, in R2015b, rts(3) is 1+2*eps but in R2016b, rts(3) is 1-2*eps

1 commentaire

Bryce Inman
Bryce Inman le 22 Nov 2016
Good point, I should have subtracted 1 to see the error.

Connectez-vous pour commenter.

Catégories

En savoir plus sur MATLAB dans Centre d'aide et File Exchange

Produits

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by