

%  Speichern sie die Funktionen vor Gebrauch jeweils in 
%  bisect.m bzw. newton.m ab
%

function m =bisect(fun,a,b)
% calculate zero of function by bisection
%           x=bisect(fun,a,b)
% input: a<b mit f(a)f(b) < 0
%        fun ... character string with name of file which evaluates function  
% output: m approximation to zero with |m-xi|< tol/2
% termination: b-a < eps
% tol set below

tol = 1e-10;

fa=feval(fun,a);fb=feval(fun,b);
if fa*fb>=0
   error('f(a)f(b) nicht negativ!')
end
n=0;
while b-a > tol
   n=n+1;
   m=(a+b)/2;
   fm=feval(fun,m);
   if fm*fb <=0
      a=m;
      fa=fm;
   else
      b=m;
      fb=fm;
   end
end
  

function [x,X] = newton(fun,x,n)
% Newton's method to calculate the root of a  differentiable function
% usage: X = newton(fun,x,n)
% input: fun ... character string with filename evaluating f and df
%        x ... initial guesss for root
%        n ... maximumnumber of iterations
%output: x ... approximation to root
%        X ... iteration history
% Note: fun has to return  f and its derivative
%       [f,df] = fun(x);

maxiter=n;
tol=1e-10;
term = inf;
X=x;
iter = 1;

while term > tol & iter <= maxiter
   %  [y1,y2] = fun(x); 
   [y1,y2] = feval(fun,x);
   if abs(y2) < 1e-14
      disp('derivative close to zero')
      return
   end
   xold = x;
   x= x- y1/y2;
   X=[X;x];
   term = abs(x-xold);
   iter=iter+1;
end
if iter >= maxiter
    disp('keine Konvergenz')
end
return
