Home » Source Code » » Matrix.pas

## Matrix.pas ( File view )

From：
• By 2010-08-31
• View(s)：3
• Point(s)： 1
```			unit Matrix;

interface

uses
Windows, Messages, SysUtils, Classes, Graphics, Controls, Forms, Dialogs,
ExtCtrls,stdctrls, buttons,  grids, Menus, comctrls;

type
TMatrix = array of array of Extended;
TDataToFloatFun=function (Str:String; var Value:Extended):Boolean of Object;

//------------------------------------
Type
TMatrixS=Array of TMatrix;
var
MatrixErrStr:String;

procedure FreeMatrix(var Matrix:TMatrix);

//转置矩阵
function RollMatrix(Matrix:TMatrix; var Value:TMatrix):Boolean;
//求方阵行列式
function AbsMatrix(Matrix:TMatrix; var Value:Extended):Boolean;
//
function IsMatrixOdd(Matrix:TMatrix):Boolean;
//求方阵的逆
function AthWart(Matrix:TMatrix; var Value:TMatrix):Boolean;
//将Matrix1与Matrix2相加(Anti=False)或相减(Anti=True)
function AddMatrix(Matrix1,Matrix2:TMatrix; var Value:TMatrix; Anti:Boolean=False):Boolean;//加法
//将Matrix1(M*K)和Matrix2(K*N)相乘得Value(M*N)
function MulMatrix(Matrix1,Matrix2:TMatrix; var Value:TMatrix):Boolean;  //乘法

//用选主元的高斯消元法求A*X=B 的X方程解
function SysLin(A,B:TMatrix; var X:TMatrix):Boolean;
//用松弛迭代法解方程
function SorGS(A,B:TMatrix; Var X:TMatrix; W,V:Extended; N:Integer):Boolean;

implementation
uses Math,ZXDDialogs;

type
PMatCom=^TMatCom;
TMatCom=Array of Record
Matrix:TMatrix;
FunIndex:Byte;
NexFun:Byte;
Fun:PMatCom;
end;
procedure FreeMatrix(var Matrix:TMatrix);
//var i:Integer;
begin// for i:=0 to High(Matrix) do SetLength(Matrix[i],0);
SetLength(Matrix,0,0);
end;
//转置     行列式值  逆          加          减          乘
const _Roll=1; _Abs=2;  _athwart=3; _Add=128+1; _Sub=128+2; _Mul=128+3;
ReserveStrCount=12;
ReserveStr:Array [0..ReserveStrCount-1] of String=('ROLLMAT','ABSMAT','ATHWART',
'+','-','*','(',')','[',']','{
','
}');

function RollMatrix(Matrix:TMatrix; var Value:TMatrix):Boolean;
var
W,H,i,j:Integer;
begin
Result:=False;
H:=Length(Matrix); if H=0 then exit;
W:=Length(Matrix[0]); if W=0 then exit;
SetLength(Value,H,W); //for i:=0 to H-1 do SetLength(Value[i],W);
for i:=0 to H-1 do for j:=0 to W-1 do Value[j,i]:=Matrix[i,j];
Result:=True;
end;

procedure GetMatrix(Matrix:TMatrix; N,R:Integer; var Mat:TMatrix);
var
i,j:Integer;
begin
SetLength(Mat,R,R);// for i:=0 to R-1 do SetLength(Mat[i],R);
for i:=0 to N-1 do for j:=0 to R-1 do Mat[j,i]:=Matrix[j+1,i];
for i:=N+1 to R do for j:=0 to R-1 do Mat[j,i-1]:=Matrix[j+1,i];
end;
function GetMatValue(Matrix:TMatrix; R:Integer):Extended;
var i:Integer; M:TMatrix;
begin
if R=0 then begin Result:=Matrix[0,0]; exit; end;
Result:=0;
for i:=0 to R do  begin
GetMatrix(Matrix,i,R,M);
if i mod 2=0
then Result:=Result+Matrix[0,i]*GetMatValue(M,R-1)
else Result:=Result-Matrix[0,i]*GetMatValue(M,R-1);
end;
SetLength(M,0,0);// for i:=0 to R-1 do SetLength(M[i],0); SetLength(M,0);
end;
function AbsMatrix(Matrix:TMatrix; var Value:Extended):Boolean;
var
W,H:Integer;
begin
Result:=False;
H:=Length(Matrix); if H=0 then exit;
W:=Length(Matrix[0]); if W<>H then exit;
Value:=GetMatValue(Matrix,H-1);
Result:=True;
end;

procedure GetMatrixS(Matrix:TMatrix; R,M,N:Integer; var Mat:TMatrix);
var
i,j:Integer;
begin
for i:=0 to N-1 do  for j:=0 to M-1 do Mat[j,i]:=Matrix[j,i];
for i:=0 to N-1 do  for j:=M+1 to R do Mat[j-1,i]:=Matrix[j,i];
for i:=N+1 to R do  for j:=0 to M-1 do Mat[j,i-1]:=Matrix[j,i];
for i:=N+1 to R do  for j:=M+1 to R do Mat[j-1,i-1]:=Matrix[j,i];
end;
function ConCom(Matrix:TMatrix; var Value:TMatrix):Boolean; //求伴随矩阵
var
W,H,L,i,j:Integer;
X:Extended;
M:TMatrix;
begin
Result:=False;
H:=Length(Matrix);    if H=0 then exit;
W:=Length(Matrix[0]); if W<>H then exit;

L:=H-1; SetLength(M,L,L);// for i:=0 to L-1 do SetLength(M[i],L);
SetLength(Value,H,H);    // for i:=0 to L do SetLength(Value[i],H);

for i:=0 to L do  for j:=0 to L do
begin
GetMatrixS(Matrix,L,i,j,M);
AbsMatrix(M,X);
if (i+j)mod 2=0 then Value[i,j]:=X
else Value[i,j]:=-X;
end;
SetLength(M,0,0);// FreeMatrix(M);
Result:=True;
end;
function AthWart(Matrix:TMatrix; var Value:TMatrix):Boolean;//求逆矩阵
var
X:Extended;
i,j,R:Integer;
Temp:TMatrix;
begin
Result:=False;
if AbsMatrix(Matrix,X) then begin
if X=0 then begin MatrixErrStr:='该方阵是奇异的'; exit; end;
Result:=True;

R:=High(Matrix);  if R=0 then begin Value[0,0]:=1/X; exit; end;

SetLength(Value,R+1,R+1); //for i:=0 to R do SetLength(Value[i],R+1);
ConCom(Matrix,Temp);
for i:=0 to R do for j:=0 to R do Value[j,i]:=Temp[i,j]/X;
FreeMatrix(Temp);
end else MatrixErrStr:='该矩阵不是方阵';
end;

function AddMatrix(Matrix1,Matrix2:TMatrix; var Value:TMatrix; Anti:Boolean=False):Boolean;//加法
var
H,W,i,j:Integer;
begin
Result:=False; MatrixErrStr:='加法运算的矩阵不匹配';
H:=High(Matrix1);    if H<>High(Matrix2) then exit; if H<0 then exit;
W:=High(Matrix1[0]); if W<>High(Matrix2[0]) then exit; if W<0 then exit;
SetLength(Value,H+1,W+1);// for i:=0 to H do SetLength(Value[i],W+1);
if Anti
then for i:=0 to H do for j:=0 to W do Value[i,j]:=Matrix1[i,j]-Matrix2[i,j]
else for i:=0 to H do for j:=0 to W do Value[i,j]:=Matrix1[i,j]+Matrix2[i,j];
Result:=True;
end;
function MulMatrix(Matrix1,Matrix2:TMatrix; var Value:TMatrix):Boolean;  //乘法
var
H1,WH,W2,
i,j,k:Integer;
X:Extended;
begin
Result:=False; MatrixErrStr:='乘法运算的矩阵不匹配';
H1:=High(Matrix1); if H1<0 then exit;
WH:=High(Matrix1[0]); if WH<0 then exit;

if (H1=0) and (WH=0) then begin
WH:=High(Matrix2); if WH<0 then exit;
W2:=High(Matrix2[0]); if W2<0 then exit;
Result:=True; SetLength(Value,WH+1,W2+1); // for i:=0 to WH do SetLength(Value[i],W2+1);
X:=Matrix1[0,0];
for i:=0 to WH do for j:=0 to W2 do Value[i,j]:=Matrix2[i,j]*X;
end;

if WH<>High(Matrix2) then exit;
W2:=High(Matrix2[0]); if W2<0 then exit;
Result:=True; SetLength(Value,H1+1,W2+1); //for i:=0 to H1 do SetLength(Value[i],W2+1);
for i:=0 to H1 do for j:=0 to W2 do begin
X:=0; for K:=0 to WH do X:=X+Matrix1[i,K]*Matrix2[k,j];
Value[i,j]:=X;
end;
end;

type RCoVecE=Array of Integer; RCOVec=Array of Extended;

function TrianguleO(A:TMatrix;var AT:TMatrix;var Jx,Iy:RCOVecE;DimMat:Integer):Boolean;
var i, k     : integer;
im , jm  : integer;
ATkk     : extended;
C        : extended;
Max,Det  : extended;
SignDet  : integer;
ErrorMat : Integer;

PROCEDURE InitJxIy;  //原序号
var i:integer;
begin
for i:=0 to pred(DimMat) do
begin
Jx[i]:=i;
Iy[i]:=i;
end;
end;

PROCEDURE CopiAAT;
var i,j:integer;
begin
for i:=0 to pred(DimMat) do
for j:=0 to pred(DimMat) do
AT[i,j]:=A[i,j];
end;

PROCEDURE Reduit;     //消元
var i,j:integer;
begin
ATkk:=AT[k,k];
for i:=k+1 to pred(DimMat) do
begin
C   :=AT[i,k]/ATkk;
for j:=k+1 to pred(DimMat) do
AT[i,j]:=AT[i,j]-C*AT[k,j];
for j:=0 to k do
if j<>k
then AT[i,j]:=AT[i,j]-C*AT[k,j]
else AT[i,j]:=-C;
end;
end;

PROCEDURE ChercheMax(var Max:extended;var im,jm:integer); //寻找最大值
var i,j : integer;
M   : extended;
begin
Max:=abs(AT[k,k]);
im:=k;jm:=k;
for i:=k to pred(DimMat) do
for j:=k to pred(DimMat) do
begin
M:=abs(AT[i,j]);
if M>Max then
begin
Max:=M;
im:=i;jm:=j;
end;
end;
end;

PROCEDURE PermuteCol;  //交换列
var Ind : integer;
i   : integer;
x   : extended;
begin
for i:=0 to pred(DimMat) do
begin
x:=AT[i,k];
AT[i,k]:=AT[i,jm];
AT[i,jm]:=x;
end;
Ind   :=Jx[k];
Jx[k] :=Jx[jm];
Jx[jm]:=Ind;
SignDet:=-SignDet;
end;

PROCEDURE PermuteRow; //交换行
var Ind : integer;
j   : integer;
x   : extended;
begin
for j:=0 to pred(DimMat) do
begin
x:=AT[k,j];
AT[k,j]:=AT[im,j];
AT[im,j]:=x;
end;
Ind:=Iy[k];
Iy[k]:=Iy[im];
Iy[im]:=ind;
SignDet:=-SignDet;
end;

begin
InitJxIy;
CopiAAT;
k:=0;
SignDet:=1;    ErrorMat:=0;
while (k<=pred(DimMat)-1) and (ErrorMat=0) do
begin
ChercheMax(Max,im,jm);
if Max<>0
then
begin
if k<>jm then PermuteCol;
if k<>im then PermuteRow;
Reduit;
end
else begin ErrorMat:=1; MatrixErrStr:='最大值为0'; end;
k:=k+1;
end;
Det:=SignDet;
for i:=0 to pred(DimMat) do
Det:=Det*AT[i,i];
if Det=0 then Result:=False
else Result:=True;
end;
PROCEDURE ResouTrO(A:TMatrix;var X,Y:RCOVec;Jx,Iy:RCOVecE; DimMat:Integer);
var i, j     : integer;
PVX      : RCOVec;
Xi       : extended;

begin
SetLength(PVX,DimMat);
for i:=pred(DimMat) downto 0 do
begin
Xi:=Y[Iy[i]];
for j:=0 to i-1 do Xi:=Xi+A[i,j]*Y[Iy[j]];
for j:=i+1 to pred(DimMat) do Xi:=Xi-A[i,j]*PVX[j];
PVX[i]:=Xi/A[i,i];
end;
for i:=0 to pred(DimMat) do X[Jx[i]]:=PVX[i];
Setlength(PVX,0);
end;
function SysLin(A,B:TMatrix; var X:TMatrix):Boolean;
var AT       : TMatrix;
Jx,Iy    : RCOVecE;
X2,B2    : RCOVec;
i:Integer;
DimMat:Integer;
begin
Result:=False; MatrixErrStr:='不匹配的矩阵';
DimMat:=Length(A);
if DimMat=0 then exit; if DimMat<>Length(A[0]) then exit;
if Length(B)<>DimMat then exit; if High(B[0])<>0 then exit;

SetLength(AT,DimMat,DimMat);
SetLength(Jx,DimMat);
SetLength(Iy,DimMat);
SetLength(X2,DimMat);
SetLength(B2,DimMat);
for i:=0 to DimMat-1 do B2[i]:=B[i,0];

Result:=True;
if TrianguleO(A,AT,Jx,Iy,DimMat)
then begin ResouTrO(AT,X2,B2,Jx,Iy,DimMat);
SetLength(X,DimMat,1);
for i:=0 to DimMat-1 do X[i,0]:
...
...
```
...
Expand＞ ＜Close

## Want complete source code? Download it here

Point(s): 1

0 lines left, continue to read

## File list

Tips: You can preview the content of files by clicking file names^_^
Name Size Date
<数学小帮手>0.00 B11-01-05 09:31
AutoSaveDesk.dcr480.00 B18-09-01 11:19
AutoSaveDesk.pas14.46 kB18-06-02 21:17
DataMaker.pas44.69 kB18-06-02 21:25
HighLightedMemo.dcr488.00 B30-07-01 16:05
HighlightedMemo.pas8.69 kB18-06-02 11:24
HSBitMapSD.pas14.07 kB18-06-02 21:25
Matrix.pas11.63 kB18-06-02 21:26
PackageZXD.cfg286.00 B18-06-02 11:27
PackageZXD.dof1.03 kB18-06-02 11:27
PackageZXD.dpk786.00 B18-06-02 11:26
PackageZXD.res1.50 kB18-06-02 11:25
PScMain.cfg386.00 B18-06-02 21:28
PScMain.dof1.07 kB18-06-02 21:28
PScMain.dpr202.00 B21-03-02 16:30
PScMain.res876.00 B18-06-02 21:28
ScienceDraw.dcr480.00 B30-03-01 10:44
ScienceDraw.pas48.17 kB18-06-02 11:21
UnitMain.ddp51.00 B18-06-02 11:27
UnitMain.dfm189.78 kB18-06-02 11:27
UnitMain.pas98.32 kB18-06-02 11:27
ZXDDialogs.pas4.58 kB18-06-02 21:28
...
• Sent successfully!
• 1 point

### Matrix.pas (76.73 kB)

Need 1 point

Your Point isn't enough.

Get point immediately by PayPal

More(Debit card / Credit card / PayPal Credit / Online Banking)

Submit your source codes. Get more point

Don't have an account？ Register now
Need any help?
Mail to: support@codeforge.com

### 切换到中文版？

CodeForge Chinese Version
CodeForge English Version

### ^_^"Oops ...

Sorry!This guy is mysterious, its blog hasn't been opened, try another, please!

### Warm tip!

Favorite by Ctrl+D