function fzero(a:double, b:double, tol:double):double
// a,b - границы интервала, в котором находится нуль функции
// tol - точность нахождения нуля
var x,c,d,fa,fb,fc,dxalfa,beta,gama;
  if a>b then begin c=a; a=b; b=c end;
  x=a; fa=fun(x); if fa=0 then goto fin;
  x=b; fb=fun(x); if fb=0 then goto fin;
  if (fa>0)=(fb>0) then goto fin;
ab:
  dx=(b-a)/2; c=a+dx;
// a,c,b - три равноотстоящие точки
  x=c; if (c-a)<=(abs(tol)) then goto fin;
  fc=fun(x); if fc=0 then goto fin;
  alfa=4*fc/(fb-fa); beta=(fa+fb)/(fb-fa);
  gama=-alfa/(1+sqrt(abs(1+alfa*(alfa-2*beta))));
  d=c+gama*dx;
// d - нуль квадратичной функции, полученной путем интерполяции
// по трем точкам
  if (fc>0)=(fb>0) then
  begin b=c; fb=fc end else begin a=c; fa=fc end;
  if (d-a)<(b-d) then c=d+(d-a) else c=d+(d-b);
  if (c<=a) or (c>=b) then goto ab;
  x=c; fc=fun(x); if fc=0 then goto fin;
  if (fc>0)=(fb>0) then
  begin b=c; fb=fc end else begin a=c; fa=fc end;
  goto ab;
fin:
  fzero=x;
  end;