Текст подпрограммы и версий
zf31r_p.zip  zf31e_p.zip  		Тексты тестовых примеров
tzf31r_p.zip  tzf31e_p.zip 

Подпрограмма:  ZF31R (модуль ZF31R_p)

Назначение

Решение системы  N нелинейных уравнений с  N неизвестными методом Ньютона.

Математическое описание

Пусть дана система  N нелинейных уравнений с  N неизвестными. 

             fi ( x1, x2, ..., xN )  =  0 ,   i = 1, 2, ..., N

Для решения заданной системы подпрограмма ZF31R на  k - м шаге итерационного процесса Ньютона вычисляет решение следующей линейной системы методом Гаусса с выбором ведущего элемента по столбцу: 

                   N
                   ∑  αi j δx j( k )  =  β i ,
                  j =1

где  αi j = ∂fi /∂xj - якобиан системы, и  βi = - fi, а затем вычисляет (k + 1) - е приближение к вектору решения на основе ранее найденного  k - го приближения: 

            x ik+1  =  x i(k) + δx i(k) ,   i = 1, 2, ..., N

Итерационный процесс по  k повторяется до тех пор, пока не будет выполнен один из двух критериев сходимости: 

                    N
(1)               ∑   | δx i(k) |  ≤  EPSX
                  i =1
                    N
(2)               ∑   | f i ( x1(k+1), x2(k+1), ..., xN(k+1))  | ≤  EPSF
                  i =1

Если же ни одно из этих двух условий не выполнено за ITMAX итераций, то ZF31R прекращает свою работу и выдает диагностическое сообщение.

Н.С.Бахвалов. Численные методы. Изд - во "Наука", 1979.

Использование 

procedure ZF31R(FUN :Proc_F3A; var ROOT :Array of Real; N :Integer;
                EPSX :Real; EPSF :Real; ITMAX :Integer;
                var ALPHA :Array of Real; var BETA :Array of Real;
                var NLEAD :Array of Integer; var IERR :Integer);


   Параметры
FUN - имя подпрограммы вычисления  fi (x1, x2, ..., xN) и
∂fi /∂xj ǀ (x1, x2, ..., xN) ,  i=1,2,...,N,  j=1,2,...,N; оформляется в виде подпрограммы с тремя формальными параметрами:
X - вещественный вектор длины N, содержащий текущее приближение к решению исходной системы;
ALPHA - двумерный вещественный вектор размеров N на N , в который помещаются значения элементов Якобиана системы в точке X, т.е.
∂fi /∂xj ǀ (x1, x2, ..., xN) ;
BETA - вещественный вектор длины N, в который помещаются значения  - f i ( x1, x2, ..., xN );
ROOT - вещественный вектор длины N, содержащий вычисленное решение системы; перед началом работы подпрограммы вектор ROOT должен содержать начальное приближение к решению;
N - заданное число уравнений системы (тип: целый);
EPSX - первый критерий сходимости (1) (тип: вещественный);
EPSF - второй критерий сходимости (2) (тип: вещественный);
ITMAX - целая переменная, значение которой перед началом работы подпрограммы должно быть положено равным максимальному числу итераций, ориентировочно требуемых для обеспечения сходимости в соответствии с заданными критериями;
ALPHA - двумерный вещественный массив размеров N на N, используемый в качестве рабочего;
BETA - вещественный вектор длины N, используемый в качестве рабочего;
NLEAD - целый вектор длины N, используемый в качестве рабочего;
IERR - целая переменная, служащая для сообщения об ошибках, обнаруженный в ходе работы подпрограммы; при этом:
IERR=65 - когда решение системы не может быть найдено в пределах заданного числа итераций;
IERR=66 - когда якобиан заданной системы имеет особенность. 
   Версии
ZF31E - решение системы N нелинейных уравнений с N неизвестными методом Ньютона в режиме расширенной (Extended) точности; при этом параметры ROOT, EPSX, EPSF, ALPHA, BETA и формальные параметры подпрограммы FUN должны иметь тип Extended. 
   Вызываемые подпрограммы
       ASG9R -
       ASG9E   		решение системы Ax = b или ATx = b методом Гаусса с выбором ведущего элемента по столбцу; используются в подпрограммах ZF31R и ZF31E соответственно;
      UTZF10 -
      UTZF11   		подпрограммы выдачи диагностических сообщений при работе подпрограмм ZF31R и ZF31E соответственно. 
   Замечания по использованию
  В подпрограммах ZF31R и ZF31E имеется глобальная запись (структура данных) _ZF31RR с элементом elm1, который полагается равным количеству итераций, выполненных при решении системы.

Пример использования 

Unit TZF31R_p;
interface
uses
SysUtils, Math, LStruct, Lfunc, UtRes_p, FZF31R_p, ZF31R_p;

function TZF31R: String; 

implementation

function TZF31R: String;
var
N,ITMAX,_i,IERR :Integer;
EPSX,EPSF :Real;
ALРНА :Array [0..8] of Real;
ВЕТА :Array [0..2] of Real;
NLEAD :Array [0..2] of Integer;
const
RООТ :Array [0..2] of Real = ( 1.0,1.0,1.0 );
begin
Result := '';
N := 3;
EPSX := 1.E-5;
EPSF := 1.E-5;
IТМАХ := 30;
ZF31R(FZF31R,ROOT,N,EPSX,EPSF,ITMAX,ALPHA,BETA,NLEAD,IERR);
Result := Result + #$0D#$0A;
for _i:=0 to 2 do
 begin
  Result := Result + Format('%15.6f ',[ROOT[_i]]);
  if ( ((_i+1) mod 3)=0 )
   then Result := Result + #$0D#$0A;
 end;
Result := Result + #$0D#$0A;
Result := Result + Format('%5d %5d ',[_ZF31RR.elm1,IERR]) + #$0D#$0A;
UtRes('TZF31R',Result);  //вывод результатов в файл TZF31R.res
exit;
end;

end.

Unit FZF31R_p;
interface
uses
SysUtils, Math, { Delphi }
LStruct, Lfunc;

procedure FZF31R(X :Array of Real; var ALPHA :Array of Real;
                 var BETA :Array of Real);

implementation

procedure FZF31R(X :Array of Real; var ALPHA :Array of Real;
                var BETA :Array of Real);
begin
BETA[0] := -(X[0]+Exp(X[0]-1)+IntPower(X[1]+X[2],2)-27);
BETA[1] := -(X[0]*Exp(X[1]-2)+IntPower(X[2],2)-10);
BETA[2] := -(X[2]+Sin(X[1]-2)+IntPower(X[1],2)-7);
ALPHA[0] := 1+Exp(X[0]-1);
ALPHA[3] := 2*(X[1]+X[2]);
ALPHA[6] := 2*(X[1]+X[2]);
ALPHA[1] := Exp(X[1]-2);
ALPHA[4] := X[0]*Exp(X[1]-2);
ALPHA[7] := 2*X[2];
ALPHA[2] := 0.0;
ALPHA[5] := Cos(X[1]-2)+2*X[1];
ALPHA[8] := 1.0;
end;

end.

 Результаты:  ROOT  =  ( 1.0, 2.0, 3.0 ) 
                       ITER  =  7 ,     IERR  =  0