矩阵的QR分解
procedure QRBKSB(A:matrx2; N:integer; Q:matrx2; B:array of real;
var X:array of real);
var
I,J:integer; SUM:real;
begin
For I:=1 To N do
begin
Sum:=0;
For J:=1 To N do
Sum:=Sum + Q[I, J] * B[J];
X[I]:=Sum;
end;
For I:=N DownTo 1 do
begin
Sum:=X[I];
For J:=I + 1 To N do
Sum:=Sum - A[I, J] * X[J];
If A[I, I] = 0 Then ShowMessage('A is singular matrix.');
X[I]:=Sum / A[I, I];
end;
end;
procedure QRDCMP(var A:matrx2; M, N:integer;var Q:matrx2);
var
I,J,K:integer; S,T,Sgn,H,F:real;
begin
For I:=1 To M do
begin
For J:=1 To M do
Q[I, J]:=0;
Q[I, I]:=1;
end;
For K:=1 To M - 1 do
begin
S:=0;
For I:=K To M do
S:=S + Abs(A[I, K]);
If S <> 0 Then
begin
T:=0;
For I:=K To M do
begin
A[I, K]:=A[I, K] / S;
T:=T + A[I, K] * A[I, K];
end;
if A[K,K] >= 0 THEN
Sgn:=1
else
Sgn:=-1;
T:=-Sqrt(T) * Sgn;
A[K, K]:=A[K, K] - T;
H:=-T * A[K, K];
For J:=K + 1 To N do
begin
F:=0;
For I:=K To M do
F:=F + A[I, K] * A[I, J];
F:=F / H;
For I:=K To M do
A[I, J]:=A[I, J] - A[I, K] * F;
end;
For J:=1 To M do
begin
F:=0;
For I:=K To M do
F:=F + A[I, K] * Q[I, J];
F:=F / H;
For I:=K To M do
Q[I, J]:=Q[I, J] - A[I, K] * F;
end;
A[K, K]:=T * S;
For I:=K + 1 To M do
A[I, K]:=0;
end;
end;
end;