![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/182002/image.jpeg)
finding multiple roots using newton raphson
77 views (last 30 days)
Show older comments
Isaac Al-rai
on 10 Feb 2018
Answered: Radu Trimbitas
on 19 May 2022
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/170194/image.jpeg)
Hello everyone, I am being asked in a homework question to find the instants a function y(t)=4*exp(-0.3t)sin(3t+0.25) crosses zero (the roots) in the interval 0<t<4. I have developed a code that uses Newton Raphson to find roots for functions. Here is that function:
function Xs=NewtonRoot(Fun,FunDer,Xest,Err,imax)
% NewtonRoot: finds the root of Fun=0 near the point Xest using Newton's
% method.
%Fun: Name of a user-defined funtion that calculates Fun for a given x.
% FunDer: Name of a user-defined function that calculates the derivative of
% Fun for a given x.
% Xest: Initial estimate of the solution.
% Err: Maximum error
% imax: Maximum number of iterations
% Output variable
%Xs Solution
for i=1:imax
Xi=Xest-Fun(Xest)/FunDer(Xest);
if abs((Xi-Xest)/Xest)<Err
Xs=Xi;
break
end
Xest=Xi;
end
if i==imax
fprintf('Solution was not obtained in %i iterations.\n',imax)
Xs=('No answer');
end
end
And then I create a function for both y(t) and y'(t). Here is y(t):
function y=problem3fun(x)
% This is the function for damping oscillator
y=4*exp(-0.3*x)*(sin((3*x)+0.25));
end
and y'(t):
if true
function y=problem3funDer(x)
y=-exp((3*x/10)*(6*sin((3*x)+0.25))-60*cos((3*x)+0.25)/(5));
endfunction y=problem3funDer(x)
y=-exp((3*x/10)*(6*sin((3*x)+0.25))-60*cos((3*x)+0.25)/(5));
end
So I think what my problem is that I do not know how to display multiple roots. I have a few guesses (Xest) but I can only at this point guess one root and thats it. On the command line:
Tsolution=NewtonRoot(@problem3fun,@problem3funDer,1,0.001,20)
Tsolution =
1.0000
So my question is how can I make my function work to find multiple roots of this y(t).
0 Comments
Accepted Answer
John D'Errico
on 10 Feb 2018
Edited: John D'Errico
on 10 Feb 2018
Newton's method cannot be used to find multiple roots. You basically find one root, and you are done.
Then you can start at a new point. Hopefully Newton's method will converge to another root. Or it might just find the first one.
So a simple approach for 1-d problem is to evaluate the function at a bunch of points, looking to see where there might be zero crossings. A root MUST exist in the area of any zero crossing. But you might miss two roots, if they are close to each other. So sample sufficiently tightly so this will not happen. But wherever your tests see a crossing, start Newton's method there, and you SHOULD probably find another root.
There is no magic in the above - just careful programming.
Now, could you use homotopy (continuation) methods to find multiple roots? Well yes, in some cases, yes. But that is wildly beyond just using Newton's method. And it involves an ODE solver, which is NOT at all Newton's method. I don't think this is the goal of your assignment, unless you wanted to invest a lot of effort.
Ok, there are other tricks. You can use a root deflation scheme, so as you find a root, you modify the function, so the root you just found is no longer a root. This can get tricky too, if you are not careful.
The point is, you cannot simply just modify Newton's method to find multiple roots. The simple approach is as I suggested, just look for sign changes in a sequence of function evaluations.
yfun = @(x) 4*exp(-0.3*x).*(sin((3*x)+0.25));
xtest = linspace(0,4,1000);
ytest = yfun(xtest);
indp2n = strfind(sign(ytest),[1,-1]);
indn2p = strfind(sign(ytest),[-1,1]);
xstart = (xtest([indp2n,indn2p]) + xtest([indp2n,indn2p]+1))/2;
Note my use of .* in yfun. That will be necessary.
ezplot(yfun,[0 4])
grid on
hold on
plot(xstart,yfun(xstart),'ro')
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/182002/image.jpeg)
Now, just start your Newton scheme at each of those points. It will converge in only a few iterations from each start point, since Newton is quadratically convergent, and those are darn good starting values.
2 Comments
John D'Errico
on 10 Feb 2018
Edited: John D'Errico
on 11 Feb 2018
Part of the problem may be that you got the derivative wrong? ;-)
syms x
y = 4*exp(-0.3*x)*(sin((3*x)+0.25))
y =
4*exp(-(3*x)/10)*sin(3*x + 1/4)
diff(y,x)
ans =
12*exp(-(3*x)/10)*cos(3*x + 1/4) - (6*exp(-(3*x)/10)*sin(3*x + 1/4))/5
subplot(2,1,1)
ezplot(diff(y,x),[0 4])
subplot(2,1,2)
ezplot(-exp((3*x/10)*(6*sin((3*x)+0.25))-60*cos((3*x)+0.25)/(5)),[0 4])
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/187768/image.jpeg)
The lower plot is the derivative you used. I know they look very much alike. But they are just slightly different, which may cause failure of Newton's method. You might have other problems.
If you STILL had problems, I'd probably learn to use the debugger, checking every step. Watch to see what happens. But I'll bet fixing your derivative might just help. I wonder if the derivative that you wrote may have had some mis-placed parens when you derived it.
More Answers (1)
Radu Trimbitas
on 19 May 2022
The standard solution, if you know the multiplicities in advance is to use the reccurence relation
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1004660/image.png)
where m is the multiplicity. This modification has the order of convergence equal to 2. The ordinary Newton for multiple roots has order 1. If the multiplicity m is unknown, you can estimate it using
![](https://www.mathworks.com/matlabcentral/answers/uploaded_files/1004665/image.png)
You may use this estimation in the previous relation.
0 Comments
See Also
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!