roots(x) vs root= interp1(y,x,0)

1 view (last 30 days)
Matt Thomas
Matt Thomas on 24 May 2022
Commented: Matt Thomas on 24 May 2022
I tried to use the root = interp1 command to find the root of my function and it does not find the correct root, however, the roots command does. what is the difference?
x= linspace(-100,100,100);
y=85.*x- 91.*(x.^2)+ 44.*(x.^3)- 8.*(x.^4)+(x.^5)-26;
plot(x,y)
grid on
p=[1 -8 44 -91 85 -26]
r=roots(p)
thanks

Answers (3)

the cyclist
the cyclist on 24 May 2022
It's not clear, at least to me, what you mean by "root = interp1 command". But, my best guess is that if you were somehow trying to use interp1 to find roots, you would run into the problem that it does linear interpolation by default, and root-finding for your equation will be non-linear.
Can you post the code that you thought should work, but did not give the result you expected?
  1 Comment
Matt Thomas
Matt Thomas on 24 May 2022
root=interp1(y,x,0) rather than roots(p)
the interp1 command yielded 0.9xxx rather than 0.5xxx

Sign in to comment.


Paul
Paul on 24 May 2022
The query points in x are spaced too far apart for the linear interpolation used by interp1 to be accurate relative to the actual curve. Plotting the curve around the real root shows why the original data yields an apparent root at r = 0.9xxxx when the actual root is around r = 0.5xxx
p = [1 -8 44 -91 85 -26];
x1 = linspace(-100,100,100);
y1 = polyval(p,x1);
x2 = linspace(-100,100,1000);
y2 = polyval(p,x2);
plot(x1,y1,'-x',x2,y2,'-o'),grid
xlim([-1.1 2.1])
roots(p)
ans =
2.5907 + 4.4333i 2.5907 - 4.4333i 1.1308 + 0.7012i 1.1308 - 0.7012i 0.5570 + 0.0000i
[interp1(y1,x1,0) interp1(y2,x2,0)]
ans = 1×2
0.9713 0.5702
  2 Comments
Torsten
Torsten on 24 May 2022
You need more points with a smaller distance between the x-values.
The maximum distance between the x-values is an estimate of the maximum error in determining the root of the function.
x= linspace(-3,3,1000);
y=85.*x- 91.*(x.^2)+ 44.*(x.^3)- 8.*(x.^4)+(x.^5)-26;
x0 = interp1(y,x,0)

Sign in to comment.


John D'Errico
John D'Errico on 24 May 2022
Edited: John D'Errico on 24 May 2022
You need to understand the mathematics. Lost you will be without that (to misquote someone named Yoda.)
For this polynomial:
syms x
y=85.*x- 91.*(x.^2)+ 44.*(x.^3)- 8.*(x.^4)+(x.^5)-26;
fplot(y,[-2,5])
grid on
There seems to be a root near zero. Thus a root is where the function has a value of y=0. But when it is flat like that, this is often an indication there may be multiple roots, or at least several roots that are very close to each other. Or you may have a complex root relatively close to that point too. And the scale on the y axis makes it difficult to be sure.
grid on
p=[1 -8 44 -91 85 -26];
roots(p)
ans =
2.5907 + 4.4333i 2.5907 - 4.4333i 1.1308 + 0.7012i 1.1308 - 0.7012i 0.5570 + 0.0000i
However, roots tells us there is only one real root of the polynomial. It lies at x = 0.5570. while the other roots are complex.
What happens when you tried to use interp1? What does interp1 do here? Actually, interp1 is a terribly bad idea. When you flip x and y like that, for a function that has a nearly zero slope, then flipping y and x produces a function with nearly infinite slope. And using an INTERPOLATION function to approximate something with a nearly infinite slope is a BAD IDEA.
You used
x = linspace(-100,100,100);
so 100 points, that varied from -100 to 200. That means a spacing of roughly 2 between each point. Do you see how much curvature you have in that function? Then using interp1, which by default uses LINEAR interpolation is just completely misunderstanding what is happening. The only points that are near zero in x at all were:
x(49:53)
ans = 1×5
-3.0303 -1.0101 1.0101 3.0303 5.0505
yfun = matlabFunction(y)
yfun = function_handle with value:
@(x)x.*8.5e+1-x.^2.*9.1e+1+x.^3.*4.4e+1-x.^4.*8.0+x.^5-2.6e+1
plot(yfun(x(49:53)),x(49:53),'o-')
grid on
the problem is, if you use interp1 here, it uses only linear interpolation. And that means it tries to interpolate the connect-the-dots curve drawn there. Using interp1 here to solve for a root of this function is again, just a flat out bad idea. I know that people have told you to do that online. But they were wrong, at least to tell you to just swap x and y, and then to use interp1, if done without any understanding of what is happening.) That is certainly true if you use such a coarsely sampled grid in x.
CAN we do a better job, using interp1 more intelligently? Well, yes. By using a more finely sampled grid, we might do this:
xbetter = linspace(-10,10,100);
ybetter = yfun(xbetter);
xinterproot = interp1(ybetter,xbetter,0)
xinterproot = 0.5698
xbetter = linspace(-10,10,1000);
ybetter = yfun(xbetter);
xinterproot = interp1(ybetter,xbetter,0)
xinterproot = 0.5572
xbetter = linspace(-10,10,10000);
ybetter = yfun(xbetter);
xinterproot = interp1(ybetter,xbetter,0)
xinterproot = 0.5570
So you see that by use of a finely sampled grid, the LINEAR interpolation can do a better job. However, we could also have used a spline interpolant with interp1. This is still a dangerous thing to do with interp1 using a spline, because splines themselves have serious problems.
xbetter = linspace(-10,10,100);
ybetter = yfun(xbetter);
xinterproot = interp1(ybetter,xbetter,0,'spline')
xinterproot = 0.5549
xbetter = linspace(-10,10,1000);
ybetter = yfun(xbetter);
xinterproot = interp1(ybetter,xbetter,0,'spline')
xinterproot = 0.5570
which you can see that even on the moderately coarsely sample grid, did a decent job. I did not need to go to as finely sampled a grid to get 4 digits correct there, when using a spline.
Just use roots, which will produce the correct solution, here accurate to full double precision accuracy. Or, you could have used solve or vpasolve. Since this is a 5th degree polynomial, we cannot compute an algebraic form of the roots of a 5th degree polynomial anyway, so vpasolve is correct here.
vpasolve(y)
ans = 
  1 Comment
Matt Thomas
Matt Thomas on 24 May 2022
thanks for that nice reply with all the examples. Yeah, i am very familar with the math. Im in my 13 semester of math at the moment. However, i have never used matlab or excell or any type of programing. Im old school and do iterations by hand, usually 4 of them.
Usually what i do is take the deriv of the function and use that to find my slopes of 0, concave and convex points, then just graph it. But had to use matlab for this one haha.
i have sent this thread to my class so they can see what is going on. most of them are really good at matlab, except me haha

Sign in to comment.

Tags

Community Treasure Hunt

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

Start Hunting!