|
Текст подпрограммы и версий aln1r_p.zip |
Тексты тестовых примеров taln1r_p.zip |
Вычисление приближения к псевдорешению системы линейных алгебраических уравнений на множестве неотрицательных векторов.
Для системы линейных алгебраических уравнений:
A * X = B
ищется вектор X с неотрицательными компонентами, минимизирующий невязку
|| A * X - B ||
Подробное описание алгоритма смотри в [1] (алгоритм NNLS).
| 1. | Lawson C.L., Hanson R.J. "Solving Least Squares Problem", Prentice - Hall Inc., Englewood Cliffs, New Jersey, 1974. |
procedure ALN1R(var A :Array of Real; MDA :Integer;
M :Integer; N :Integer; var B :Array of Real;
var X :Array of Real; var RNORM :Real;
var W :Array of Real; var ZZ :Array of Real;
var INDEX :Array of Integer; var IERR :Integer);
Параметры
| A - | вещественный двумерный массив размера MDA на N, содержащий при входе в подпрограмму в первых M строках матрицу A (I, J) размерности M на N. По окончании работы подпрограммы содержит произведение матриц Q * A, где Q - ортогональная матрица, размера M на M, вычисленная данной подпрограммой; |
| MDA - | число строк массива A (тип: целый); |
| M, N - | число строк и столбцов матрицы A, соответственно (тип: целый); |
| B - | при входе в подпрограмму содержит вещественный вектор длины M - правую часть уравнения. В результате работы подпрограммы содержит вектор Q * B; |
| X - | вещественный вектор длины N - в результате работы подпрограммы содержит решение; |
| RNORM - | вещественная переменная, содержащая по окончании работы подпрограммы евклидову норму вектора невязки; |
| W - | вещественный рабочий вектор длины N; по окончании работы подпрограммы содержит вектор двойственного решения, удовлетворяющего следующим соотношениям (см.[1]): |
| W(I) = 0. - | для всех I из множества P; |
| W(I) ≤ 0. - | для всех I из множества Z; |
| ZZ - | вещественный рабочий вектор длины M; |
| INDEX - |
целый рабочий вектор длины N; по окончании работы
подпрограммы определяет множества P и Z следующим образом (см.[1]): INDEX(1) до INDEX(NSETZ) = множество P INDEX(IZ1) до INDEX(N) = множество Z IZ1 = NSETZ + 1 = NPP1; |
| IERR - | диагностический параметр - код завершения программы со значениями: |
| IERR= 0 - | решение получено успешно; |
| IERR=65 - | неправильная размерность параметров; |
| IERR=66 - | требуется более чем 3*N итераций. |
Версии: нет
Вызываемые подпрограммы: нет
Замечания по использованию
| 1. | Подпрограмма ALN1R обращается к вспомогательным подпрограммам с именами: ALN1R1, ALN1R2, ALN1R3, ALN1R4. |
Unit TALN1R_p;
interface
uses
SysUtils, Math, { Delphi }
Lstruct, Lfunc, UtRes_p, ALN1R_p;
function TALN1R: String;
implementation
function TALN1R: String;
var
_i,IERR :Integer;
RNORM :Real;
X :Array [0..2] of Real;
W :Array [0..2] of Real;
ZZ :Array [0..1] of Real;
INDЕХ :Array [0..2] of Integer;
const
A :Array [0..5] of Real = ( 1.0,2.0,1.0,2.0,1.0,1.0 );
B :Array [0..1] of Real = ( 2.0,4.0 );
A1 :Array [0..5] of Real = ( 1.0,1.0,1.0,-1.0,0.0,0.0 );
B1 :Array [0..1] of Real = ( 1.0,3.0 );
begin
Result := ''; { результат функции }
ALN1R(A,2,2,3,B,X,RNORM,W,ZZ,INDEX,IERR);
Result := Result + Format('%s',[' ТЕСТ ПОДПРОГРАММЫ ALN1R' + #$0D#$0A]) + #$0D#$0A;
Result := Result + Format('%s',[' СИСТЕМА УРАВНЕНИЙ: ' + #$0D#$0A]) + #$0D#$0A;
Result := Result + Format('%s',[' X + Y + Z = 2' + #$0D#$0A +
' 2*X + 2*Y + Z = 4' + #$0D#$0A]) + #$0D#$0A;
Result := Result + Format('%s',[' IERR: ']);
Result := Result + Format('%15d ',[IERR]) + #$0D#$0A;
Result := Result + Format('%s',[' РЕШЕНИЕ: ']);
Result := Result + #$0D#$0A;
for _i:=0 to 2 do
begin
Result := Result + Format('%15.5f ',[X[_i]]);
if ( ((_i+1) mod 3)=0 )
then Result := Result + #$0D#$0A;
end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' НЕВЯЗКА: ']);
Result := Result + Format('%15.5f ',[RNORM]) + #$0D#$0A;
Result := Result + Format('%s',[' W: ']);
Result := Result + #$0D#$0A;
for _i:=0 to 2 do
begin
Result := Result + Format('%15.5f ',[W[_i]]);
if ( ((_i+1) mod 3)=0 )
then Result := Result + #$0D#$0A;
end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' INDEX: ']);
Result := Result + #$0D#$0A;
for _i:=0 to 2 do
begin
Result := Result + Format('%15d ',[INDEX[_i]]);
if ( ((_i+1) mod 3)=0 )
then Result := Result + #$0D#$0A;
end;
Result := Result + #$0D#$0A;
ALN1R(A1,2,2,3,B1,X,RNORM,W,ZZ,INDEX,IERR);
Result := Result + Format('%s',[' СИСТЕМА УРАВНЕНИЙ:' + #$0D#$0A]) + #$0D#$0A;
Result := Result + Format('%s',[' IERR: ']);
Result := Result + Format('%15d ',[IERR]) + #$0D#$0A;
Result := Result + Format('%s',[' РЕШЕНИЕ: ']);
Result := Result + #$0D#$0A;
for _i:=0 to 2 do
begin
Result := Result + Format('%15.5f ',[X[_i]]);
if ( ((_i+1) mod 3)=0 )
then Result := Result + #$0D#$0A;
end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' НЕВЯЗКА: ']);
Result := Result + Format('%15.5f ',[RNORM]) + #$0D#$0A;
Result := Result + Format('%s',[' W: ']);
Result := Result + #$0D#$0A;
for _i:=0 to 2 do
begin
Result := Result + Format('%15.5f ',[W[_i]]);
if ( ((_i+1) mod 3)=0 )
then Result := Result + #$0D#$0A;
end;
Result := Result + #$0D#$0A;
Result := Result + Format('%s',[' INDEX: ']);
Result := Result + #$0D#$0A;
for _i:=0 to 2 do
begin
Result := Result + Format('%15d ',[INDEX[_i]]);
if ( ((_i+1) mod 3)=0 )
then Result := Result + #$0D#$0A;
end;
Result := Result + #$0D#$0A;
UtRes('TALN1R',Result); { вывод результатов в файл TALN1R.res }
exit;
end;
end.
Рассмотрим решение системы уравнений:
x + y + z = 2
2*x + 2*y + z = 4
Результаты:
IERR: 0
решение: 2.00000e+00 .00000e+00 .00000e+00
невязка: 8.11612e-08
w: 0.00000e + 00 .00000e + 00 -3.62964e-08
INDEX: 1 2 3
Рассмотрим решение системы уравнений:
x + y + 0 * z = 1
x - y + 0 * z = 3
Результаты:
IERR: 0
решение: 2.00000e+00 .00000e+00 .00000e+00
невязка: 1.41421e+00
W: .00000e + 00 -2.00000e + 00 .00000e + 00
INDEX: 1 2 3