function sol = sec(fun,a,b)

fa = feval(fun,a);
fb = feval(fun,b);
flag = 0;
k = 0;

kmax = 100;
tol = 1e-4;

while (flag == 0)

   x = a + fa/(fa-fb)*(b-a);
   F = feval(fun,x);
   crit = fa*F;

   if (crit>0)
       a = x;
       fa= F;
   else
       b = x;
       fb= F;
   endif

   k = k+1;
   if abs(a-b)<tol, flag =1; endif
   if abs(F)<tol, flag =1; endif
   if k>kmax, flag =1; endif

endwhile

sol(1) = x;
sol(2) = F;
sol(3) = k;

endfunction