Code covered by the BSD License  

Highlights from
MATLAB en zijn Symbolic Math Toolbox

image thumbnail

MATLAB en zijn Symbolic Math Toolbox

by

 

Leergang voor geïnteresseerden in de Computeralgebra - Systemen MATLAB en zijn Symbolic Math Toolbox

ZeroBisect(f, a, b, aTol, rTol, maxIter, trace)
function [xZero, fZero, nIter] = ZeroBisect(f, a, b, aTol, rTol, maxIter, trace)

%
%                 ************
% Funnctie-file : ZeroBisect.m
%                 ************

% Bepaling van het nulpunt van de functie f met behulp van de bisectie methode
%
%  Input parameters : 
% -------------------
%
% --- f          : string die de functie definieert
% --- [a, b]     : het interval waar het nulpunt zich in bevindt
% --- aTol, rTol : respectievelijk absolute- en relatieve tolerantie
% --- maxIter    : maximaal aantal iteraties
% --- trace      : trace-functie : 1 'aan', 0 'af'
%
% Output parameters :
% -------------------
%
% --- xZero      : gevonden nulpunt
% --- fZero      : f(xZero)
% --- nIter      : benodigd aantal iteraties
%

fa  = feval(f, a);
sfa = sign(fa);
fb  = feval(f, b);
sfb = sign(fb);

if sfa == sfb
   error('Nulpuntsbepaling onmogelijk : geen nulpunts-interval')
end

if     sfa == 0
   xZero = a;
   fZero = 0.0;
   nIter = 0;
elseif sfb == 0
   xZero = b;
   fZero = 0.0;
   nIter = 0;
else
   if trace
      disp(' ')
      disp(sprintf ...
         ('**** traceer - optie ****                     a = %15.12e', a))
      disp(sprintf ...
         ('                                              b = %15.12e', b))
      disp(' ')
      disp(sprintf ...
         (' k             a                    b                   f(M)'))
   end
   err = abs(b - a);
   for k = 1 : maxIter
      err = err / 2;
      M   = (a + b) / 2;
      fM  = feval(f, M);
      sfM = sign(fM);
      if trace
         if rem(k, 10) == 1
            disp(' ')
         end % if rem  
         disp(sprintf('%3d  %+15.12e  %+15.12e  %+15.12e', k, a, b, fM))
      end % if trace
      if (sfM == 0) | (err < abs(M) * rTol + aTol)
         xZero = M;
         fZero = fM;
         nIter = k;
         break; % binnen tolerantie
      else
         if sfM == sfa
            a  = M;
            fa = fM;
         else
            b  = M;
            fb = fM;
         end
      end       % buiten tolerantie
   end % for k
   if k == maxIter
      warning('Nulpuntsbepaling niet binnen de vereiste nauwkeurigheid')
   else
      disp(' ')
      disp('Bisectie-methode convergeerde binnen de vereiste nauwkeurigheid')
   end % if k
end %else

%
%                ************
% Functie-file : ZeroBisect.m
%                ************
%

Contact us