uses graph, crt;

const eps = 0.00001;         { точность решения задачи }
      Dim = 2;             { размерность }
      omega = 1;

type  vect = array [1..Dim] of real;

var dF, x, xk, New_X, dx, dg, gk: vect;
    H, dH: array [1..Dim, 1..Dim] of real;
    xtrace: array[0..200] of vect;
    isolinies, step, sminx, smaxx, sminy, smaxy, xbegin, ybegin,
    offsety, iter, o, Error, is: integer;
    lambda, kc, kx, ky, da,
    minx, maxx, miny, maxy, minz, maxz: real;
    Com: ( Quit, ZoomIn, ZoomOut, Left, Right, Up, Down, DBetter, DWorst,
           Nothing, WBigger, WSmaller, Repaint, AddIsol, DelIsol );
    num: string;
    Digits: array [1..12,1..5,1..4] of integer;
    Ch: text;

Procedure gr_init; forward;

Function F( var x: vect ): real;        { минимизируемая функция }
begin
  F:=100*(x[2]-x[1]*x[1])*(x[2]-x[1]*x[1])+(1-x[1])*(1-x[1]);
 { F:=x[1]+x[2]*x[2]; }
end;

procedure Make_dF;            { частные производные }
begin
  dF[1]:=-400*(x[2]-x[1]*x[1])*x[1]-2*(1-x[1]);
  dF[2]:=200*(x[2]-x[1]*x[1]);
{  dF[1]:=1;
  dF[2]:=2*x[2];    }
end;

procedure Init;               { выбор начальной точки и }
var i, j, k: integer;         { начального приближения матрицы Гессе }
begin
  gr_init;
  for i:=1 to Dim do
  begin
    for  j:=1 to Dim do
      H[i,j]:=0;
    x[i]:=0;
    H[i,i]:=1;
    dg[i]:=1;
  end;
  x[1]:=10;
  x[2]:=10;
  xtrace[0][1]:=x[1];
  xtrace[0][2]:=x[2];
  Com:=Nothing;
  iter:=0; Error:=0;
end;

procedure put_digit( x, y: integer; Num: integer );
var st: string;
begin
  str( Num, st );
  setcolor( BLACK );
  if Num>10 then outtextxy( x, y, '██' )
  else outtextxy( x, y, '█' );
  setcolor( YELLOW );
  outtextxy( x, y, st );
end;

procedure gr_init;
var i, j: integer;
    z: real;
    grDriver : Integer;
    grMode   : Integer;
    ErrCode  : Integer;
begin
  step:=3; isolinies:=100;
  sminx:=00; smaxx:=450;
  sminy:=00; smaxy:=450;
  minx:=-2 ; maxx:=12 ;
  miny:=-2 ; maxy:=12 ;
  kx:=(smaxx-sminx)/(maxx-minx);
  ky:=(smaxy-sminy)/(maxy-miny);
  xbegin:=1;
  ybegin:=1;
  maxz:=-1e10; minz:=1e10;
  i:=sminx;
  while( i<smaxx ) do
  begin
    i:=i+3;
    x[1]:=minx+(i-sminx)/kx;
    j:=sminy;
    while( j<smaxy-step ) do
    begin
      j:=j+3;
      x[2]:=miny+(j-sminy)/ky;
      z:=F(x);
      if z>maxz then maxz:=z;
      if z<minz then minz:=z;
    end;
  end;
  kc:=isolinies/(maxz-minz);
  grDriver := Detect;
  InitGraph(grDriver,grMode,'d:\tp\bgi');
  ErrCode := GraphResult;
  if ErrCode <> grOk then
  begin
    WriteLn('Graphics error:',
            GraphErrorMsg(ErrCode));
    halt(1);
  end;
end;

procedure Make_dg;
var i: integer;
begin
  for i:=1 to Dim do
  begin
    dg[i]:=-dF[i];
    x[i]:=New_x[i];
  end;
  Make_dF;
  for i:=1 to Dim do
    dg[i]:=dg[i]+dF[i];
end;

procedure Make_H;
var i, j: integer;
    delta, k: vect;
    s: real;
begin
  Make_dg;
  for i:=1 to Dim do
  begin
    delta[i]:=0;
    for j:=1 to Dim do
      delta[i]:=delta[i]+H[i,j]*dg[j]
  end;
  s:=0;
  for i:=1 to Dim do
  begin
    k[i]:=dx[i]-delta[i];
    s:=s+k[i]*dg[i];
  end;
  if s<>0 then
  begin
    for i:=1 to Dim do
      for j:=1 to Dim do
        H[i,j]:=H[i,j]+k[i]*k[j]/s;
  end
  else Error:=1;
end;

function gold: real;
var xx1, xx2: vect;
    x1, x2, a, b, delta: real;
    i: integer;
begin
  a:=0; b:=1.5*da; delta:=0.000001;
  while  abs(b-a) > 1e-5 do
  begin
    for i:=1 to Dim do
    begin
      xx1[i]:=x[i]+dx[i]*(a+b-delta)/2;
      xx2[i]:=x[i]+dx[i]*(a+b+delta)/2;
    end;
    if ( F( xx1 ) < F( xx2 ) ) then b:=(a+b+delta)/2
    else a:=(a+b-delta)/2;
  end;
  gold:=b; da:=b;
end;

{function gold: real;
label L5, L2, L10, L11, L32, L40, _end;
var a, r, t, s: real;
    temp_x, xx: vect;
    i, k: integer;
begin
  a:=F(x);
  is:=0;
  for i:=1 to Dim do
    if dx[i]<>0 then goto L32;
  write('\n Вектоp напpавления поиска нулевой!\n');
  is:=1;
  goto _end;

 L32:
  k:=0;
  for i:=1 to Dim do
    xx[i]:=x[i]+da*dx[i];
  r:=F(xx);
  if a<r then goto L5;

 L2:
   for i:=1 to Dim do
     xx[i]:=xx[i]+da*dx[i];
  t:=F(xx);
  if r<t then goto L10;
  da:=da*2;
  r:=t;
  goto L2;

 L5:
  k:=k+1;
  da:=da/2;
  t:=r;
  for i:=1 to Dim do
    xx[i]:=x[i]+da*dx[i];
  r:=F(xx);
  if r <= a then goto L10;
  if k <= 20 then goto L5;
  is:=1;
  goto _end;

 L10:
  s:=a+t-2*r;
  if s=0. then goto _end;

 L11:
  s:=da*(a*3+t-4*r)/(s*2);
  da:=s;
  gold:=da;
  goto _end;

  is:=0;

  for i:=1 to Dim do
    if dx[i]<>0. then goto L32;
       write('\n Вектоp напpавления поиска нулевой!\n');
  is:=1;
  goto _end;

           {
 is:=0;
 da:=0.1;
 a:=F(x);
 k:=0;
 for i:=1 to Dim do
   x[i]:=temp_x[i]+da*dx[i];
 r:=F(x);
 if(a<r) then goto L5;

L2:
  for i:=1 to Dim do
    x[i]:=x[i]+da*dx[i];
 t:=F(x);
 if(r<t) then goto L10;
 da:=2*da;
 r:=t;
 goto L2;

L5:
 k:=k+1;
 da:=da/2;
 t:=r;
 for i:=1 to Dim do
   x[i]:=temp_x[i]+da*dx[i];
 r:=F(x);
 if( r <= a ) then goto L10;
 if( k <= 20 ) then goto L5;
 is:=1;
 goto L40;

L10:
 s:=a+t-2*r;
 if( s = 0 ) then goto L40;

L11:
 s:=da*(a*3+t-4*r)/(s*2);
 gold:=s;
L40:

 _end:
end; }

procedure Make_New_X;
var i, j: integer;
begin
  for i:=1 to Dim do
    xk[i]:=x[i];
  Make_dF;
  for i:=1 to Dim do
  begin
    dx[i]:=0;
    for j:=1 to Dim do
      dx[i]:=dx[i]-H[i,j]*dF[j];
  end;
  lambda:=gold;
  for i:=1 to Dim do
    dx[i]:=dx[i]*lambda;
  for i:=1 to Dim do
    New_x[i]:=x[i]+dx[i];
end;

function Norma_dF: real;
var i: integer;
    s: real;
begin
  s:=0;
  for i:=1 to Dim do
    s:=s+dF[i]*dF[i];
  Norma_dF:=s;
end;

procedure mputpixel( i, j: real; c: integer );
begin
  if c<16 then setfillstyle( 1, c+1 )
  else setfillstyle( 1, 0 );
  bar( sminx+round((i-minx)*kx),
       round((smaxy-(j-miny)*ky)),
       sminx+round((i-minx)*kx)+(step-1),
       round((smaxy-(j-miny)*ky)-(step-1)) );
end;

procedure paintx;
var i: integer;
begin
  setviewport( sminx, sminy, smaxx, smaxy, ClipOn );
  setcolor( LIGHTRED );
  for i:=1 to iter do
  begin
    line( round((xtrace[i-1][1]-minx)*kx),
          round((smaxy-sminy)-(xtrace[i-1][2]-miny)*ky),
          round((xtrace[i][1]-minx)*kx),
          round((smaxy-sminy)-(xtrace[i][2]-miny)*ky) );
    mputpixel( xtrace[i-1][1], xtrace[i-1][2], YELLOW );
  end;
  for i:=1 to iter do
  begin
    put_digit( round((xtrace[i][1]-minx)*kx),
               round((smaxy-sminy)-(xtrace[i][2]-miny)*ky), i );
  end;
  setviewport( 0, 0, getmaxx, getmaxy, ClipOn );
end;

procedure mline( x1, y1, x2, y2: real; c: integer );
begin
  setlinestyle( 0, 1, 3 );
  setcolor( c );
  line( sminx+round((x1-minx)*kx), round((smaxy-(y1-miny)*ky)),
        sminx+round((x2-minx)*kx), round((smaxy-(y2-miny)*ky)) );
end;

Procedure ParamsOut;
var num: string[9];
begin
  setcolor( WHITE );
  setfillstyle( 9, BLUE );
  offsety:=sminx+80;
  bar( o, offsety, o+130, offsety+70 );
  rectangle( o, offsety, o+130, offsety+70 );
  outtextxy( o+4, offsety+10, 'Параметры рис.' );
  str( isolinies, num );
  outtextxy( o+4, offsety+25, concat('Изол.:', num ) );
  str( step, num );
  outtextxy( o+4, offsety+35, concat('Точка:', num ) );
  str( minz:8:3, num );
  outtextxy( o+4, offsety+45, concat('Fmin: ', num ) );
  str( maxz:8:3, num );
  outtextxy( o+4, offsety+55, concat('Fmax: ', num ) );
end;

procedure Paint;
var i, j, c, c1: integer;
    x_temp: vect;
    z: real;
begin
  kx:=(smaxx-sminx)/(maxx-minx);
  ky:=(smaxy-sminy)/(maxy-miny);
  maxz:=-1e10; minz:=1e10;
  for i:=1 to Dim do
    x_temp[i]:=x[i];
  i:=sminx;
  while( i<smaxx ) do
  begin
    i:=i+3;
    x[1]:=minx+(i-sminx)/kx;
    j:=sminy;
    while( j<smaxy ) do
    begin
      j:=j+3;
      x[2]:=miny+(j-sminy)/ky;
      z:=F(x);
      if z>maxz then maxz:=z;
      if z<minz then minz:=z;
    end;
  end;
  kc:=isolinies/(maxz-minz);

  setfillstyle( 7, DARKGRAY );
  bar( sminx, sminy, smaxx, smaxy );
  c:=0;
  i:=sminx;
  while( i<smaxx-step ) do
  begin
    x[1]:=minx+(i-sminx)/kx;
    j:=sminy;
    while( j<smaxy-step ) do
    begin
      x[2]:=miny+(j-sminy)/ky;
      z:=F(x);
      c1:=round(( z-minz )*kc);
      if c1<>c then
        begin
          mputpixel( x[1], x[2], YELLOW );
          c:=c1;
        end;
        j:=j+step;
    end;
    i:=i+step;
  end;
  for i:=1 to Dim do
    x[i]:=x_temp[i];
  mline( minx,miny,maxx,miny, WHITE );
  mline( maxx,miny,maxx,maxy, WHITE );
  mline( maxx,maxy,minx,maxy, WHITE );
  mline( minx,maxy,minx,miny, WHITE );

  if ( 0>minx ) and ( 0<maxx ) then
    mline( 0, miny, 0, maxy, LIGHTGREEN );
  if ( 0>miny ) and ( 0<maxy ) then
    mline( minx, 0, maxx, 0, LIGHTGREEN );
  o:=smaxx+xbegin+20;
  setcolor( WHITE );
  setfillstyle( 9, BLUE );
  offsety:=sminx;
  setcolor( WHITE );
  bar( o, offsety, o+130, offsety+70 );
  rectangle( o, offsety, o+130, offsety+70 );
  outtextxy( o+4, offsety+10, 'Параметры окна' );
  str( minx:8:3, num );
  outtextxy( o+4, offsety+25,  concat('Xmin:', num ) );
  str( maxx:8:3, num );
  outtextxy( o+4, offsety+35, concat('Xmax:', num ) );
  str( miny:8:3, num );
  outtextxy( o+4, offsety+45, concat('Ymin:', num ) );
  str( maxy:8:3, num );
  outtextxy( o+4, offsety+55, concat('Ymax:', num ) );

  ParamsOut;
  offsety:=offsety+80;
  bar( o, offsety, o+130, offsety+070 );
  rectangle( o, offsety, o+130, offsety+070 );
  outtextxy( o+4, offsety+10, ' Алгоритм' );
  str( eps:8:5, num );
  outtextxy( o+4, offsety+25, concat('eps1:', num ) );
  str( eps:8:5, num );
  outtextxy( o+4, offsety+35, concat('eps2:', num ) );
  str( xtrace[0][1]:8:5, num );
  outtextxy( o+4, offsety+45, concat('x0  :', num ) );
  str( xtrace[0][2]:8:5, num );
  outtextxy( o+4, offsety+55, concat('x0  :', num ) );
  offsety:=offsety+70;
  bar( o, offsety, o+130, offsety+70 );
  rectangle( o, offsety, o+130, offsety+70 );
  outtextxy( o+4, offsety+10, ' Результаты' );
  str( F(x):9:4, num );
  outtextxy( o+4, offsety+25, concat('Fmin:', num ) );
  str( x[1]:9:4, num );
  outtextxy( o+4, offsety+35, concat('xmin:', num ) );
  str( x[2]:9:4, num );
  outtextxy( o+4, offsety+45, concat('ymin:', num ) );
  str( iter, num );
  outtextxy( o+4, offsety+55, concat('Итер:', num ) );

  offsety:=offsety+80;
  bar( o, offsety, o+130, offsety+130 );
  rectangle( o, offsety, o+130, offsety+130 );
  outtextxy( o+4, offsety+10, ' Помощь ' );
  outtextxy( o+4, offsety+25, 'z: Увеличить');
  outtextxy( o+4, offsety+35, 'x: Уменьшить');
  outtextxy( o+4, offsety+45, 'a: Точка больше ');
  outtextxy( o+4, offsety+55, 's: Точка меньше ');
  outtextxy( o+4, offsety+65, 'c: Окно больше');
  outtextxy( o+4, offsety+75, 'v: Окно меньше');
  outtextxy( o+4, offsety+85, 'd: Изол. больше');
  outtextxy( o+4, offsety+95, 'f: Изол. меньше');
  outtextxy( o+4, offsety+105, ' : Перерисуй');
  outtextxy( o+4, offsety+115, 'q: Выход');
  paintx;
  if Error=1 then
  begin
    setfillstyle( 1, RED );
    bar( sminx, smaxy+5, smaxx, smaxy+15 );
    setcolor( YELLOW );
    outtextxy( sminx+2, smaxy+6, 'Error!!!' );
  end;
end;

Procedure DoBrojden;
var i: integer;
   as: real;
begin
  da:=2;
  repeat
    Make_New_X;
    iter:=iter+1;
    for i:=1 to Dim do
      xtrace[iter][i]:=x[i];
    Make_H;
  until (Norma_dF < eps) or (Error=1);
  as:=Norma_dF;
end;

Procedure Menu;
var ch: char;
begin
  case readkey of
    'q': Com:=Quit;
    'z': Com:=ZoomIn;
    'x': Com:=ZoomOut;
    'a': Com:=DBetter;
    's': Com:=DWorst;
    'c': Com:=WBigger;
    'v': Com:=WSmaller;
    ' ': Com:=Repaint;
    'd': Com:=AddIsol;
    'f': Com:=DelIsol;
     #0: begin
           ch:=readkey;
           case ch of
             'K': Com:=Left;
             'M': Com:=Right;
             'H': Com:=Up;
             'P': Com:=Down;
           end;
         end
    else
      Com:=Nothing;
  end;
end;

Procedure ProcessCommand;
var hx, hy: real;
begin
  case Com of
    Left:    begin hx:=maxx-minx; maxx:=minx+hx/2; minx:=maxx-hx; end;
    Right:   begin hx:=maxx-minx; minx:=maxx-hx/2; maxx:=minx+hx; end;
    Up:      begin hy:=maxy-miny; miny:=maxy-hy/2; maxy:=miny+hy; end;
    Down:    begin hy:=maxy-miny; maxy:=miny+hy/2; miny:=maxy-hy; end;
    ZoomIn:  if (maxx-minx)>0.08 then
             begin
               hx:=(maxx-minx)/2; hy:=(maxy-miny)/2;
               maxx:=minx+hx; maxy:=miny+hy;
             end;
    ZoomOut: begin
               hx:=(maxx-minx)*2; hy:=(maxy-miny)*2;
               maxx:=minx+hx; maxy:=miny+hy;
             end;
    Dbetter: if step>1 then step:=step-1;
    Dworst:  if step<10 then step:=step+1;
    WBigger: begin
               hy:=smaxy-sminy;
               if sminy+hy*1.5<getmaxy
                 then begin
                   smaxy:=smaxy+round(hy*0.5);
                   smaxx:=smaxx+round(hy*0.5);
                 end;
             end;
    WSmaller:begin
               hy:=smaxy-sminy;
               if hy/1.5>100
                 then begin
                        smaxy:=sminy+round(hy/1.5);
                        smaxx:=sminx+round(hy/1.5);
                      end;
             end;
    AddIsol: Isolinies:=round(Isolinies*1.5);
    DelIsol: if Isolinies>1 then Isolinies:=round(Isolinies/1.5);
    Repaint: cleardevice;
    Nothing:
  end;
end;

begin
  Init;
  DoBrojden;
  Paint;
  repeat
    if ( Com=Left ) or ( Com=Right ) or
       ( Com=Up ) or( Com=Down ) or
       ( Com=ZoomIn ) or ( Com=ZoomOut ) or
       ( Com=Repaint )
    then Paint;
    if ( Com=AddIsol ) or ( Com=DelIsol ) or
       ( Com=Dbetter ) or ( Com=Dworst )
    then ParamsOut;
    Menu;
    ProcessCommand;
  until Com=Quit;
  closegraph;
end.
