1. warning: divede by zero

Hello. I have a problem with a function to find the zero of a function by the method of secants. this is the text of my function:

function[pf,y1,err]=secant(f,p0,p1,delta,epsilon,max1)
%--------------------------------------------------------------------------
% Questa function permette di trovare lo zero di una funzione col metodo
% delle secanti
%
%parametri in ingresso:
% f= funzione
% p0= primo valore iniziale
% p1= secondo valore iniziale
% delta= tolleranza per p1
% epsilon= tolleranza per y1
% max1= numero max di iterazioni
%
%parametri in uscita:
% p1= soluzione
% y1= valore di f in p1
% err= errore stimato
%--------------------------------------------------------------------------

y0=feval(f,p0);
y1=feval(f,p1);
for k=1:max1
df=(y1-y0)/(p1-p0);
if df==0
dp=0;
else
dp=y1/df;
end
p2=p1-dp;
y2=feval(f,p2);
err=abs(dp);
relerr=err/(abs(p2)+eps);
p0=p1;
y0=y1;
pf=p2;
y1=y2;
if (err<delta)|(relerr<delta)|(abs(y2)<epsilon)
return
end
end

I tried to run it with a simple function (x ^ 2-16) giving initial values p0=3 and p1=5 but it says that the line 23 ( df=(y1-y0)/(p1-p0); )divided by zero! and then I returned as a result pf = 5 ... line 23 seems to me ok ... where is the error?

ps. I apologize for my bad English but I'm Italian ...Ascolta
Trascrizione fonetica

2. Originally Posted by sty86
Hello. I have a problem with a function to find the zero of a function by the method of secants. this is the text of my function:
Code:
function[pf,y1,err]=secant(f,p0,p1,delta,epsilon,max1)
%--------------------------------------------------------------------------
% Questa function permette di trovare lo zero di una funzione col metodo
% delle secanti
%
%parametri in ingresso:
% f= funzione
% p0= primo valore iniziale
% p1= secondo valore iniziale
% delta= tolleranza per p1
% epsilon= tolleranza per y1
% max1= numero max di iterazioni
%
%parametri in uscita:
% p1= soluzione
% y1= valore di f in p1
% err= errore stimato
%--------------------------------------------------------------------------

y0=feval(f,p0);
y1=feval(f,p1);
for k=1:max1
df=(y1-y0)/(p1-p0);
if df==0
dp=0;
else
dp=y1/df;
end
p2=p1-dp;
y2=feval(f,p2);
err=abs(dp);
relerr=err/(abs(p2)+eps);
p0=p1;
y0=y1;
pf=p2;
y1=y2;
if (err<delta)|(relerr<delta)|(abs(y2)<epsilon)
return
end
end
I tried to run it with a simple function (x ^ 2-16) giving initial values p0=3 and p1=5 but it says that the line 23 ( df=(y1-y0)/(p1-p0); )divided by zero! and then I returned as a result pf = 5 ... line 23 seems to me ok ... where is the error?
I don't know the programming language that you are using, but one thing seems clear to me: you are iterating a fixed number of times. So if, after a number of iterations, your p1 and p0 get so close together that they cannot be distinguished given the floating point accuracy that your run-time system is using, then you get that p1-p0=0, and that's it. Your division is effectively of the form 0/0.

So, I think you need to insert a test for y1=y0 (or p1=p0) before you do that division on line #23.
P.S: I should add that I don't quite understand your code and don't feel like I really want to understand it in detail.