1 { Version 990830. Copyright � Alexey A.Chernobaev, 1996-1999 }
2
3 unit SLS_Iter;
4 {
5 ������� ������� �������� ��������� ������������ ������� ����������� �������.
6
7 �������: Ax = f;
8 ������������ �������: (x(n+1) - x(n))/t(n+1) + Ax(n)=f;
9 �������: r(n)=Ax(n) - f;
10 t(n+1)=(r(n), Ar(n)) / (Ar(n), Ar(n)).
11
12 ����������� ������� ����������: A > 0.
13 }
14
15 interface
16
17 {$I VCheck.inc}
18
19 uses
20 ExtType, Aliasv, Aliasm, F32v, F64v, F80v, F32m, F64m, F80m, VectErr;
21
SolveLinearSystemIternull22 function SolveLinearSystemIter(A: TFloatMatrix; f: TFloatVector;
23 MaxIter: Integer; Precision: Float; x: TFloatVector): Bool;
24 { �� �����:
25 A - ������� ������������� (����������); f - ������ ������ �����;
26 MaxIter - ������������ ���������� ��������, ������� ����������� �����������;
27 Precision - �������� (����� ���������� ��� |r(n)| <= Precision).
28
29 A � f ������ ���� ���������� �����������.
30
31 �� ������: ��������� ����� True, ���� ���������� ��������� ��������, � False,
32 ���� ��������� ����� �� ���������� ���������� �������� MaxIter;
33 x - ������������ �������. }
34
35 implementation
36
SolveLinearSystemIternull37 function SolveLinearSystemIter(A: TFloatMatrix; f: TFloatVector;
38 MaxIter: Integer; Precision: Float; x: TFloatVector): Bool;
39 var
40 xn, Arn, rn: TFloatMatrix;
41 I, m: Integer;
42 begin
43 {$IFDEF CHECK_MATH}
44 if (A.RowCount <> A.ColCount) or (A.RowCount <> x.Count) or
45 (A.RowCount <> f.Count) then TFloatMatrix.Error(SErrorInParameters, [0]);
46 {$ENDIF}
47 Result:=False;
48 m:=A.RowCount;
49 xn:=TFloatMatrix.Create(m, 1, 0); { x(0):=0 }
50 rn:=TFloatMatrix.Create(m, 1, 0);
51 Arn:=TFloatMatrix.Create(m, 1, 0);
52 Precision:=Sqr(Precision);
53 try
54 rn.Vector.SubVector(f); { r(0):=-f }
55 for I:=0 to MaxIter - 1 do begin
56 if rn.Vector.SqrSum <= Precision then begin
57 Result:=True;
58 Break;
59 end;
60 Arn.MatrixProduct(A, rn); { Ar(n):=A*r(n) }
61 xn.Vector.AddScaled( - Arn.Vector.DotProduct(rn.Vector) /
62 Arn.Vector.DotProduct(Arn.Vector), rn.Vector); { x(n+1):=-t(n+1)*r(n)+x(n) }
63 rn.MatrixProduct(A, xn); { r(n):=A*x(n) }
64 rn.Vector.SubVector(f); { r(n):=A*x(n)-f }
65 end;
66 x.Assign(xn.Vector);
67 finally
68 xn.Free;
69 rn.Free;
70 Arn.Free;
71 end;
72 end;
73
74 end.
75