Нахождение минимума функции методом Фибоначчи. Полное решение

program pro;
const
  f:array [0..14] of integer=(1,1,2,3,5,8,13,21,34,55,89,144,233,377,610);
var L0,l,a,b,c,E,kn,fy,fz:real;
    ak,bk,y,z:array [0..20] of real;
    i,n,k:integer;
    Label Metka;
begin
     write('Vvedite a,b,c funkcii  ');
     readln(a,b,c);
     write('Vvedite a i b  ');
     readln(ak[0],bk[0]);
     write('Vvedite l  ');
     readln(l);
     L0 := abs(bk[0]-ak[0]);
     write('Vvedite E  ');
     readln(E);
     writeln('------------------');
     
     {1}
     writeln('1)  L.0 = [',ak[0],';',bk[0],'] , l=',l,' , E=',E);
     kn := abs(L0)/l;
     for i:=0 to 14 do
         if f[i] >= kn then
                     begin
                          n:=i;
                          break;
                     end;
     writeln('    F.',n,' = ',f[n],' > |L.0|/l = ',L0,'/',l,' = ',kn,' , N=',n);
     
     {2}
     write('2)  F.0=F.1=1');
     for i:=2 to n do
         write(';  F.',2,'=',f[i]);
     writeln();
     
     {3}
     writeln('3)  k=0');
     k:=0;
     
     {4}
     y[0] := ak[0] + f[n-2]/f[n]*(bk[0]-ak[0]);
     z[0] := ak[0] + f[n-1]/f[n]*(bk[0]-ak[0]);
     y[0] := round(y[0]*1000)/1000;
     z[0] := round(z[0]*1000)/1000;
     writeln('4)  y.0=a.0+F.',n-2,'/F.',n,'(b.0-a.0)=',ak[0],'-',f[n-2],'/',f[n],'*(',bk[0],'-',ak[0],')=',y[0]);
     writeln('    z.0=a.0+F.',n-1,'/F.',n,'(b.0-a.0)=',ak[0],'-',f[n-1],'/',f[n],'*(',bk[0],'-',ak[0],')=',z[0]);
     
     
     
     
      Metka:
     {5}
     fy := a*y[k]*y[k]+b*y[k]+c;
     fz := a*z[k]*z[k]+b*z[k]+c;
     fy := round(fy*1000)/1000;
     fz := round(fz*1000)/1000;
     writeln('5^',k,')  f(y.',k,')=f(',y[k],')=',fy);
     writeln('      f(z.',k,')=f(',z[k],')=',fz);
     
     {6}
     if fy<=fz then
        begin
            writeln('6^',k,')  f(y.',k,')<=f(z.',k,')');
            ak[k+1]:=ak[k];
            bk[k+1]:=z[k];
            z[k+1]:=y[k];
            y[k+1]:=ak[k+1]+f[n-k-3]/f[n-k-1]*(bk[k+1]-ak[k+1]);
            y[k+1] := round(y[k+1]*1000)/1000;
            writeln('      a.',k+1,'=a.',k,'=',ak[k]);
            writeln('      b.',k+1,'=z.',k,'=',z[k]);
            writeln('      z.',k+1,'=y.',k,'=',y[k]);
            writeln('      y.',k+1,'=',y[k+1]);
        end
        else
        begin
            writeln('6^',k,')  f(y.',k,')>f(z.',k,')');
            ak[k+1]:=y[k];
            y[k+1]:=z[k];
            bk[k+1]:=bk[k];
            z[k+1]:=ak[k+1]+f[n-k-2]/f[n-k-1]*(bk[k+1]-ak[k+1]);
            z[k+1] := round(z[k+1]*1000)/1000;
            writeln('      a.',k+1,'=y.',k,'=',y[k]);
            writeln('      y.',k+1,'=z.',k,'=',z[k]);
            writeln('      b.',k+1,'=b.',k,'=',bk[k]);
            writeln('      z.',k+1,'=',z[k+1]);
        end;
        
     {7}
     if k<>n-3 then
        begin
             writeln('7^',k,')  k=',k,'!=N-3=',n-3);
             writeln('      L.',k,'=[',ak[k+1],' ; ',bk[k+1],']');
             k:=k+1;
             writeln('      k=',k,'     => n.5');
             Goto Metka;
        end
        else
        begin
             writeln('7^',k,')  k=',k,'=N-3=',n-3);
             writeln('      Pust y.',n-1,' = y.',n-2,' = z.',n-2,' = ',z[n-2]);
             y[n-1] := y[n-2];
             z[n-1]:=y[n-1] + E;
             writeln('      z.',n-1,'=y.',n-1,'+E = ',z[n-1]);
             fy := a*y[n-1]*y[n-1]+b*y[n-1]+c;
             fz := a*z[n-1]*z[n-1]+b*z[n-1]+c;
             writeln('      f(y.',n-1,')=',fy);
             writeln('      f(z.',n-1,')=',fz);
             if fy<=fz then
                begin
                     writeln('      f(y.',n-1,') <= f(z.',n-1,')');
                     ak[n-1]:=ak[n-2];
                     bk[n-1]:=z[n-1];
                     writeln('      a.',n-1,' = a.',n-2,' = ',ak[n-1]);
                     writeln('      b.',n-1,' = z.',n-2,' = ',bk[n-1]);
                end
                else
                begin
                     writeln('      f(y.',n-1,') > f(z.',n-1,')');
                     ak[n-1]:=y[n-1];
                     bk[n-1]:=bk[n-2];
                     writeln('      a.',n-1,' = y.',n-1,' = ',ak[n-1]);
                     writeln('      b.',n-1,' = b.',n-2,' = ',bk[n-1]);
                end;
             writeln('      x*e[',ak[n-1],' ; ',bk[n-1],']');
             writeln('      x*=~(',ak[n-1],' + ',bk[n-1],')/2=~',(ak[n-1] + bk[n-1])/2);
        end;
     
end.