// Формирование столбца L-матрицы for j:=i+l to N do begin Rtmp := 0;
if(i >
0) then for k:=0 to i-1 do begin Rtmp2 := Rtmp;
Rtmp := Rtmp + Matr[j, k] * Matr[k,i];
if (Abs(Rtmp) <
Abs(Rtmp2) * 1E-10) then Rtmp := 0;
end;

Matr[j, i] := (Matr[j, i] - Rtmp) * Matr[i, i]; if (Abs(Matr[j, i]) < Abs(Rtmp) * 1E-10)

then Matr[j, i] := 0;
end;
// Формирование (i+1)-ой строки U-матрицы for j:=i+l to N do for k:=0 to i do begin Rtmp := Abs(Matr[i+1, j]);
Matr[i+1, j] := Matr[i+1, j] - Matr[i+l,k] * Matr[k,j];
if(Abs(Matr[i+1, j]) <
Rtmp*le-10) then Matr[i+1, j] := 0;
end;
end;
end;
function LUVector(const B: TVector;

const Matr: TMatrix; const V: Tint) : TVector; // Решение системы линейных уравнений с уже проведенным LU-разложением // (результаты в Matr и V) и с заданным вектором правых частей В var i, j, N: word; Rtmp: real; Vz: TVector; begin N: = High(Matr); if (N <> High(Matr[0])) then begin

Messages(1 Функция LUVector :1#131 матрица не квадратная, 1#13 + 'ее размерность ' + IntToStr(N+l) + ' х ' + IntToStr(High(Matr[0])+1) ) ;
exit ;
end;

if(N <> High(V)) then begin Messages('Функция LUVector :' +

#13'размер вектора не соответствует размеру матрицы');
exit;
end;
SetLength(Vz, N+l);
SetLength(Result, N+l);

// Прямой ход

Vz[0] := B[0];

for i:=l to N do begin

Rtmp := 0;
for j:=0 to i-1 do Rtmp := Rtmp + Matr[i, j] * Vz[j];
Vz[i] := B[i] - Rtmp;
end;

// Обратный ход

Vz[N] := Vz[N] * Matr[N, N];

for i:=N-l downto 0 do begin

Rtmp := 0;
for j:=i+l to N do Rtmp := Rtmp + Matr[i, j] * Vz[j];
Vz[i] := (Vz[i] - Rtmp) * Matr[i, i];
end;
// Восстановление начальной последовательности переменных for i:=0 to N do Result[V[i]] := Vz[i];
end;
function LUSystem(const A: TMatrix;
const VRight: TVector): TVector;
// Решение системы линейных уравнений с помощью LU-разложения без // выбора главного элемента с сохранением исходной матрицы var Matr: TMatrix;

V: Tint; // массив индексов с учетом перестановок

Vz : TVector; i, j, N: word; Rtmp: real; begin N:= High(A); if(N <> High(A[0])) then begin

Messages(1 Функция LUSystem :' +

#13'матрица не квадратная,1#13+

'ее размерность ' + IntToStr(N+l) + ' х 1 +

IntToStr(High(A[0])+1)) ;
exit ;
end;

SetLength(Matr, N+l, N+l); for i:=0 to N do for j:=0 to N do

Matr[i, j] := A[i, j] ;
SetLength(Vz, N+l);
SetLength(Result, N+l);
LU(Matr, V);
if (ErrorMatrix) then exit;

// Прямой ход

Vz[0] := VRight[0];

for i:=l to N do begin

Rtmp := 0;
for j:=0 to i-1 do Rtmp := Rtmp + Matr[i, j] * Vz[j];
Vz[i] := VRight[i] - Rtmp;
end;

// Обратный ход

Vz[N] := Vz[N] * Matr[N, N];

for i:=N-l downto 0 do begin

Rtmp := 0;
for j:=i+l to N do Rtmp := Rtmp + Matr[i, j] * Vz[j];
Vz[i] := (Vz[i] - Rtmp) * Matr[i, i];
end;
// Восстановление начальной последовательности переменных for i:=0 to N do Result[V[i]] := Vz[i];
end;
function LUMaxSystem(const A: TMatrix;
const VRight: TVector): TVector;
// Решение системы линейных уравнений с помощью LU-разложения // с выбором главного элемента с сохранением исходной матрицы var Matr: TMatrix;

⇐ Предыдущая страница| |Следующая страница ⇒

Приемы программирования в Delphi на основе VCL



Новости за месяц

  • Сентябрь
    2020
  • Пн
  • Вт
  • Ср
  • Чт
  • Пт
  • Сб
  • Вс