How do I find all the zeros of a function?
525 views (last 30 days)
MathWorks Support Team on 18 Oct 2013
The example below describes one way to find zeros between 0 and 2*pi.
lb = 0; % Set a lower bound for the function.
ub = 2*pi; % Set an upper bound for the function.
x = NaN*ones(100,1); % Initializes x.
% Look for the zeros in the function's current window.
%compare to theory values
Changing the step value of the loop (0.01) to a larger value will result in faster execution, but less accuracy, and vice versa when the specified step value is smaller.
The tolerance for duplicates, in this example 1e-12, should be specified depending on the problem at hand, an inappropriate setting for a certain problem might result in some of the zeros not being found.
Alternately, you can use the SOLVE function from the Symbolic Math Toolbox.
More Answers (2)
Fabien M. on 29 Apr 2016
If your function is continuous on a given range [xmin , xmax], an efficient solution consists in using the rolle's theorem.
Here is a standalone matlab code to find all zeros of a function f on a range [xmin , xmax] :
% Inputs :
% f : function of one variable
% [xmin - xmax] : range where f is continuous containing zeros
% N : control of the minimum distance (xmax-xmin)/N between two zeros
if (y1*y2<=0) % Rolle's theorem : one zeros (or more) present
z=[z,fsolve(f,(x2*y1-x1*y2)/(y1-y2))]; % Linear approximation to guess the initial value in the [x1,x2] range.
Using the Newton algorithm can fruitfully replace the fsolve function.
Here is an example to get all zeros between 0 and 100 of x-> besselj(40,x+cos(x)) :
>> AllZeros(@(x) besselj(40,x+cos(x)),0,100)
47.556 52.881 55.876 61.103 66.011 68.267 73.267 75.460 80.112 84.706 86.838 91.643 93.521 98.321
Dennis Yeh on 6 Dec 2014
Edited: Walter Roberson on 8 May 2015
EDIT: I suggested this solution in the past, but found that it does not work well.
For loops are slow!
Given a discrete set of data instead of a function, the following is one possible solution:
% Find all intersections of y1 and y2, where x is the independent variable:
x = linspace(0,10,10000); y1 = sin(x); y2 = cos(x);
coeff = polyfit(x,(y2-y1),100); % Perform a degree 100 polynomial curve fit on the difference between y1 and y2.
possible_zeros = sort(unique(abs(roots(coeff)))); % Roots of the polynomial curve fit - the absolute value is to convert complex roots, and the unique() function turns each complex pair to only 1 point
error = polyval(coeff,possible_zeros); % Evaluate the polynomial curve fit at every root - answers close to zero correspond to an intersection
tolerance = 1e-6; % Change this value to be larger if intersections are skipped, and smaller if fake intersections are returned
actual_zeros_indices = find(abs(error) < tolerance); % Find the indices of the vector "possible_zeros" that correspond to actual zeros
actual_zeros = possible_zeros(actual_zeros_indices); % Plug the indices into possible_zeros to find the real zeros
f1 = figure(1);
plot(x,y1,x,y2); % Plot the original functions
plot(x,(y2-y1),'b',x,polyval(coeff,x),'g',actual_zeros,polyval(coeff,actual_zeros),'r*') % Plot the difference between y1 and y2, the curve fit, and the zeros