InraPorc

{ ******************************************************************

  Numerical gradient

  ****************************************************************** }

procedure Gradient(X, G : TVector);

const

  Eta = 1.0E-4;     { Relative increment }

var

  I     : Integer;  { Loop variable }

  Temp  : Float;    { Temporary variable }

  Delta : Float;    { Increment }

  Xm    : Float;    { X - Delta }

  Xp    : Float;    { X + Delta }

  Fm    : Float;    { F(X - Delta) }

  Fp    : Float;    { F(X + Delta) }

begin

  for I := 1 to Nvar do

    begin

      if X[I] <> 0.0 then

        Delta := Eta * Abs(X[I])

      else

        Delta := Eta;

      Xp := X[I] + Delta;

       Xm := X[I] - Delta;

      Temp := X[I];

      X[I] := Xm;

      Fm := Func(X);

      X[I] := Xp;

      Fp := Func(X);

      G[I] := (Fp - Fm) / (2.0 * Delta);

      X[I] := Temp

    end;

end;