Программная реализация модального управления для линейных стационарных систем

procedure SetCell (ACol, ARow: Word; AValue: Float);

function GetItem (NumItem: LongInt): Float;

procedure SetItem (NumItem: LongInt; AValue: Float);

procedure SwitchRows (FirstRow, SecondRow: Word);

public

constructor Create (NCols, NRows: Word);

destructor Destroy; override;

procedure Assign (AMatrix: TMatrix);

procedure ReSize (NewCols, NewRows: Word);

procedure SetNull;

procedure SetSingle;

procedure SetNegative;

procedure AddConst (AConst: Float);

procedure AddMatrix (AMatrix: TMatrix);

procedure MultConst (MConst: Float);

procedure MultFromRight (MMatrix: TMatrix);

procedure MultFromLeft (MMatrix: TMatrix);

procedure NthPower (Power: Word);

procedure Transpose;

function Inverse: Boolean;

function Determinant: Float;

function Rang: Float;

property ColCount: Word read FCols;

property RowCount: Word read FRows;

property Cells [ACol, ARow: Word]: Float read GetCell write SetCell;

default;

property Items [NumItem: LongInt]: Float read GetItem write SetItem;

end;

implementation

uses Windows;

function IncPtr (p: Pointer; i: LongInt): Pointer;

asm

push EBX

mov EBX,EAX

add EBX,EDX

mov EAX,EBX

pop EBX

end;

function TMatrix.GetCell (ACol, ARow: Word): Float;

var

CellPtr: ^Float;

begin

CellPtr := IncPtr(DataPtr, (FRows * Pred(ACol) + Pred(ARow)) *

SizeOf(Float));

Result := CellPtr^;

end;

procedure TMatrix.SetCell (ACol, ARow: Word; AValue: Float);

var

CellPtr: ^Float;

begin

CellPtr := IncPtr(DataPtr, (FRows * Pred(ACol) + Pred(ARow)) *

SizeOf(Float));

CellPtr^ := AValue;

end;

function TMatrix.GetItem (NumItem: LongInt): Float;

var

CellPtr: ^Float;

begin

CellPtr := IncPtr(DataPtr, Pred(NumItem) * SizeOf(Float));

Result := CellPtr^;

end;

procedure TMatrix.SetItem (NumItem: LongInt; AValue: Float);

var

CellPtr: ^Float;

begin

CellPtr := IncPtr(DataPtr, Pred(NumItem) * SizeOf(Float));

CellPtr^ := AValue;

end;

procedure TMatrix.SwitchRows (FirstRow, SecondRow: Word);

var

i: Word;

Buffer: Float;

begin

for i := 1 to FCols do

begin

Buffer := GetCell(i, FirstRow);

SetCell(i, FirstRow, GetCell(i, SecondRow));

SetCell(i, SecondRow, Buffer);

end;

end;

constructor TMatrix.Create (NCols, NRows: Word);

begin

inherited Create;

FCols := NCols;

FRows := NRows;

DataPtr := AllocMem(FCols * FRows * SizeOf(Float));

end;

destructor TMatrix.Destroy;

begin

FreeMem(DataPtr);

inherited Destroy;

end;

procedure TMatrix.Assign (AMatrix: TMatrix);

var

NewMatrixSize: LongInt;

begin

NewMatrixSize := AMatrix.ColCount * AMatrix.RowCount * SizeOf(Float);

ReAllocMem(DataPtr, NewMatrixSize);

CopyMemory(DataPtr, AMatrix.DataPtr, NewMatrixSize);

FCols := AMatrix.ColCount;

FRows := AMatrix.RowCount

end;

procedure TMatrix.ReSize (NewCols, NewRows: Word);

var

NewMatrixSize: LongInt;

begin

NewMatrixSize := NewCols * NewRows * SizeOf(Float);

ReAllocMem(DataPtr, NewMatrixSize);

FCols := NewCols;

FRows := NewRows;

end;

procedure TMatrix.SetNull;

begin

ZeroMemory (DataPtr, FCols * FRows * SizeOf(Float));

end;

procedure TMatrix.SetSingle;

var

i: Word;

begin

if FCols <> FRows then

Raise EMatrixOperatingError.Create ('Единичная матрица должна быть '+

'квадратной')

else

begin

SetNull;

for i := 1 to FCols do SetCell (i, i, 1);

end;

end;

procedure TMatrix.SetNegative;

var

i: LongInt;

begin

for i := 1 to FCols * FRows do SetItem(i, - GetItem(i));

end;

procedure TMatrix.AddConst (AConst: Float);

var

i: LongInt;

begin

for i := 1 to FCols * FRows do SetItem (i, GetItem(i) + AConst);

end;

procedure TMatrix.AddMatrix (AMatrix: TMatrix);

var

i: LongInt;

begin

for i := 1 to FCols * FRows do SetItem (i, GetItem(i) + AMatrix.Items

[i]);

end;

procedure TMatrix.MultConst (MConst: Float);

var

i: LongInt;

begin

for i := 1 to FCols * FRows do SetItem (i, GetItem(i) * MConst);

end;

procedure TMatrix.MultFromRight (MMatrix: TMatrix);

var

j, i, k: Word;

DummyRes: Float;

DummyMatrix: TMatrix;

begin

DummyMatrix := TMatrix.Create (MMatrix.ColCount, FRows);

if FCols <> MMatrix.RowCount then

Raise EMatrixOperatingError.Create ('Перемножаемые матрицы должны

быть '+

'соответствующей размерности')

else

for i := 1 to FRows do

for j := 1 to MMatrix.ColCount do

begin

DummyRes := 0;

for k := 1 to FCols do

DummyRes := DummyRes + Cells[k, i] * MMatrix[j, k];

DummyMatrix[j, i] := DummyRes;

end;

Assign(DummyMatrix);

DummyMatrix.Free;

end;

procedure TMatrix.MultFromLeft (MMatrix: TMatrix);

var

j, i, k: Word;

DummyRes: Float;

DummyMatrix: TMatrix;

begin

DummyMatrix := TMatrix.Create (FCols, MMatrix.RowCount);

if MMatrix.ColCount <> FRows then

Raise EMatrixOperatingError.Create ('Перемножаемые матрицы должны

быть '+

'соответствующей размерности')

else

for i := 1 to MMatrix.ColCount do

for j := 1 to FCols do

begin

DummyRes := 0;

for k := 1 to MMatrix.ColCount do

DummyRes := DummyRes + MMatrix[k, i] * Cells[j, k];

DummyMatrix[j, i] := DummyRes;

end;

Assign(DummyMatrix);

DummyMatrix.Free;

end;

procedure TMatrix.NthPower (Power: Word);

var

i: Word;

DummyMatrix: TMatrix;

begin

DummyMatrix := TMatrix.Create (FCols, FRows);

DummyMatrix.Assign (Self);

if FCols <> FRows then

Raise EMatrixOperatingError.Create ('Возводимая в степень матрица

должна '+

'быть квадратной')

else

case Power of

0 : SetSingle;

1 : begin end;

else

for i := 2 to Power do MultFromRight (DummyMatrix);

end;

DummyMatrix.Free;

end;

procedure TMatrix.Transpose;

var

i, j: Word;

Dummy: Float;

begin

if FCols <> FRows then

Raise EMatrixOperatingError.Create ('Транспонируемая матрица должна

быть '+

'квадратной')

else

for i := 1 to FCols do

for j := 1 to FRows do

if j > i then

begin

Dummy := GetCell(j, i);

SetCell(j, i, GetCell(i, j));

SetCell(i, j, Dummy);

end

end;

function TMatrix.Inverse: Boolean;

var

DummyMatrix: TMatrix;

Divisor, Multiplier: Float;

Row, RefRow, NewRow, Term: Word;

Singular: Boolean;

begin

Singular := False;

DummyMatrix := TMatrix.Create (FCols, FRows);

if (FCols <> FRows) or (FCols = 0) then

Raise EMatrixOperatingError.Create ('Инвертируемая матрица должна

быть '+

'квадратной и ненулевого

размера');

if FCols = 1 then

if ABS(GetItem(1)) < NearlyZero then Singular := True

else DummyMatrix.Items[1] := 1 / GetItem(1);

if FCols > 1 then

begin

DummyMatrix.SetSingle;

RefRow := 0;

repeat

Inc(RefRow);

if ABS(Cells[RefRow, RefRow]) < NearlyZero then

begin

Singular := TRUE;

NewRow := RefRow;

repeat

Inc(NewRow);

if ABS(Cells[RefRow, NewRow]) > NearlyZero then

begin

SwitchRows(NewRow, RefRow);

DummyMatrix.SwitchRows(NewRow, RefRow);

Singular := False;

end;

until (not Singular) or (NewRow >= FCols);

end;

if not Singular then

begin

Divisor := Cells[RefRow, RefRow];

for Term := 1 to FCols do

begin

SetCell(Term, RefRow, GetCell(Term, RefRow)/Divisor);

DummyMatrix[Term, RefRow] := DummyMatrix[Term,

RefRow]/Divisor;

end;

for Row := 1 to FCols do

if (Row <> RefRow) and (ABS(Cells[RefRow, Row]) >

NearlyZero) then

begin

Multiplier := - Cells[RefRow, Row] / Cells[RefRow,

RefRow];

for Term := 1 to FCols do

begin

SetCell(Term, Row, GetCell(Term, Row) +

Multiplier * GetCell(Term,

RefRow));

DummyMatrix[Term, Row] := DummyMatrix[Term, Row] +

Multiplier * DummyMatrix[Term,

RefRow];

end

end;

end;

until Singular or (RefRow >= FCols);

end;

Assign(DummyMatrix);

DummyMatrix.Free;

if not Singular then Result := True

else Result := False;

end;

function TMatrix.Determinant: Float;

begin

Result := 0;

end;

function TMatrix.Rang: Float;

begin

Result := 0;

end;

end.

unit Operates;

interface

uses Matrix, Grids, SysUtils;

const

MaxArraySize = 30;

type

Float = Extended;

TOrder = 1..MaxArraySize;

ESingularMatrix = class (Exception);

type

TComplex = record

Re, Im : Float;

end;

TComplexVector = record

Data : array [1..MaxArraySize] of TComplex;

Dim : TOrder;

end;

function SymmetricalFunction (Roots: TComplexVector; K: byte): Float;

procedure DiffSystemSolve (matrixA,

matrixB: TMatrix;

LowerLimit,

UpperLimit: Float;

InitialValues: TMatrix;

NumReturn,

NumIntervals: Word;

SolutionValues: TMatrix);

implementation

function SymmetricalFunction (Roots: TComplexVector; K: byte): Float;

var

Z: TComplex;

function SummComplex (FirstNC, SecondNC: TComplex): TComplex;

begin

Result.Re := FirstNC.Re + SecondNC.Re;

Result.Im := FirstNC.Im + SecondNC.Im;

end;

function MultComplex (FirstNC, SecondNC: TComplex): TComplex;

begin

Result.Re := FirstNC.Re * SecondNC.Re - FirstNC.Im * SecondNC.Im;

Result.Im := FirstNC.Re * SecondNC.Im + FirstNC.Im * SecondNC.Re;

end;

function DivComplex (FirstNC, SecondNC: TComplex): TComplex;

var

Z: Float;

begin

Z := Sqr(SecondNC.Re) + Sqr(SecondNC.Im);

Result.Re := (FirstNC.Re * SecondNC.Re + FirstNC.Im * SecondNC.Im) / Z;

Result.Im := (FirstNC.Im * SecondNC.Re - FirstNC.Re * SecondNC.Im) / Z;

end;

function CombinationSumm (LowLimit, HighLimit, K: byte): TComplex;

var i: byte;

begin

Result.Re := 0;

Result.Im := 0;

if LowLimit = HighLimit then Result := Roots.Data[LowLimit]

else

for i := LowLimit to HighLimit - K + 1 do

if K = 1 then Result := SummComplex(Result, Roots.Data [i])

else Result := SummComplex(Result,

MultComplex(Roots.Data [i],

CombinationSumm(i + 1,

HighLimit, K-

1)));

end;

begin

Z := CombinationSumm(1, Roots.Dim, K);

Result := Z.Re;

end;

procedure DiffSystemSolve (matrixA, matrixB: TMatrix;

LowerLimit, UpperLimit: Float;

InitialValues: TMatrix;

NumReturn, NumIntervals: Word;

SolutionValues: TMatrix);

type

Ptr = ^Data;

Data = record

Values: TMatrix;

Next: Ptr;

end;

var

ValuesStack: Ptr;

Spacing, HalfSpacing: Float;

Index, Term: Word;

F1, F2, F3, F4,

CurrentValues,

TempValues: TMatrix;

NumEquations, NumTimeCol: Word;

function TargetALL (matrixA, mayrixB: TMatrix; Values: TMatrix; KRow:

Word): Float;

var

j: Word;

begin

try

Result := matrixB.Items[KRow];

for j := 1 to NumEquations do

Result := Result + matrixA[j, KRow] * Values.Items[j];

except

on EO: EOverflow do EO.Message := 'Не буду считать !!!'#10 +

'С уменьшите разброс

коэффициентов в матрице А'#10 +

'либо измените опции (уменьшите

их pls.)';

end;

end;

procedure Push (var ValuesStack: Ptr;

CurrentValues: TMatrix);

var

NewNode : Ptr;

begin

New(NewNode);

NewNode^.Values := TMatrix.Create(NumTimeCol, 1);

NewNode^.Values.Assign(CurrentValues);

NewNode^.Next := ValuesStack;

ValuesStack := NewNode;

end; { procedure Push }

procedure Pop (var ValuesStack: Ptr;

CurrentValues: TMatrix);

var

OldNode : Ptr;

begin

OldNode := ValuesStack;

ValuesStack := OldNode^.Next;

CurrentValues.Assign(OldNode^.Values);

OldNode^.Values.Free;

Dispose(OldNode);

end; { procedure Pop }

procedure GetValues(NumReturn, NumIntervals: Word; var ValuesStack: Ptr;

SolutionValues: TMatrix);

var

Index, Term: Integer;

j: Word;

CurrValues: TMatrix;

begin

SolutionValues.ReSize(NumTimeCol, Succ(NumReturn));

CurrValues := TMatrix.Create(NumTimeCol, 1);

Term := NumIntervals;

for Index := NumReturn downto 0 do

begin

Pop(ValuesStack, CurrValues);

Dec(Term);

while (Term / NumIntervals >= Index / NumReturn) and (Term >= 0) do

begin

Pop(ValuesStack, CurrValues);

Dec(Term);

end;

for j := 1 to NumTimeCol do

SolutionValues[j, Succ(Index)] := CurrValues.Items[j];

end;

CurrValues.Free;

end; { procedure GetValues }

procedure Step(Spacing: Float; CurrentValues: TMatrix; F: TMatrix);

var

i : byte;

begin

for i := 1 to NumEquations do

F.Items[i] := Spacing * TargetALL (matrixA, matrixB, CurrentValues,

i);

end; { procedure Step }

begin

NumEquations := matrixA.RowCount;

NumTimeCol := Succ(NumEquations);

ValuesStack := nil;

Spacing := (UpperLimit - LowerLimit) / NumIntervals;

CurrentValues := TMatrix.Create(1, 1);

CurrentValues.Assign(InitialValues);

CurrentValues.ReSize(NumTimeCol, 1);

CurrentValues.Items[NumTimeCol] := LowerLimit;

TempValues := TMatrix.Create(NumTimeCol, 1);

F1 := TMatrix.Create(NumTimeCol, 1);

F2 := TMatrix.Create(NumTimeCol, 1);

F3 := TMatrix.Create(NumTimeCol, 1);

F4 := TMatrix.Create(NumTimeCol, 1);

Push(ValuesStack, CurrentValues);

HalfSpacing := Spacing / 2;

for Index := 1 to NumIntervals do

begin

{ First step - calculate F1 }

Step(Spacing, CurrentValues, F1);

TempValues.Items[NumTimeCol] := CurrentValues.Items[NumTimeCol] +

HalfSpacing;

for Term := 1 to NumEquations do

TempValues.Items[Term] := CurrentValues.Items[Term] + 0.5 *

F1.Items[Term];

{ 2nd step - calculate F2 }

Step(Spacing, TempValues, F2);

for Term := 1 to NumEquations do

TempValues.Items[Term] := CurrentValues.Items[Term] + 0.5 *

F2.Items[Term];

{ Third step - calculate F3 }

Step(Spacing, TempValues, F3);

TempValues.Items[NumTimeCol] := CurrentValues.Items[NumTimeCol] +

Spacing;

for Term := 1 to NumEquations do

TempValues.Items[Term] := CurrentValues.Items[Term] +

F3.Items[Term];

{ Fourth step - calculate F4[1]; first equation }

Step(Spacing, TempValues, F4);

{ Combine F1, F2, F3, and F4 to get }

{ the solution at this mesh point }

CurrentValues.Items[NumTimeCol] := CurrentValues.Items[NumTimeCol]

+ Spacing;

for Term := 1 to NumEquations do

CurrentValues.Items[Term] := CurrentValues.Items[Term] +

(F1.Items[Term] + 2 * F2.Items[Term]

+

2 * F3.Items[Term] +

F4.Items[Term]) /6;

Push(ValuesStack, CurrentValues);

end;

GetValues(NumReturn, NumIntervals, ValuesStack, SolutionValues);

F1.Free;

F2.Free;

F3.Free;

F4.Free;

CurrentValues.Free;

TempValues.Free;

end;

end.

unit HelpUnit;

interface

uses

Windows, Messages, SysUtils, Classes, Graphics, Controls, Forms,

Dialogs,

StdCtrls, ExtCtrls, Buttons;

type

TForm_Help = class(TForm)

BitBtn1: TBitBtn;

Bevel1: TBevel;

Label1: TLabel;

Label2: TLabel;

Label3: TLabel;

Label4: TLabel;

Label5: TLabel;

Label6: TLabel;

Label7: TLabel;

private

{ Private declarations }

public

{ Public declarations }

end;

var

Form_Help: TForm_Help;

implementation

{$R *.DFM}

end.

unit OptsUnit;

interface

uses

Windows, Messages, SysUtils, Classes, Graphics, Controls, Forms,

Dialogs,

StdCtrls, Spin, Buttons, ExtCtrls;

type

TForm_Options = class(TForm)

CheckBox_Link: TCheckBox;

SpinEdit1: TSpinEdit;

SpinEdit2: TSpinEdit;

SpinEdit3: TSpinEdit;

Label_UpLimit: TLabel;

Label_PointsNumber: TLabel;

Label_Intervals: TLabel;

Label_1: TLabel;

BitBtn_Ok: TBitBtn;

BitBtn_Cancel: TBitBtn;

SpinEdit0: TSpinEdit;

Label1: TLabel;

Bevel1: TBevel;

procedure SpinEdit0Change(Sender: TObject);

procedure SpinEdit2Change(Sender: TObject);

procedure CheckBox_LinkClick(Sender: TObject);

private

{ Private declarations }

public

{ Public declarations }

end;

var

Form_Options: TForm_Options;

implementation

uses MainUnit, SubUnit;

{$R *.DFM}

procedure TForm_Options.SpinEdit0Change(Sender: TObject);

begin

SpinEdit1.MinValue := Succ(SpinEdit0.Value);

if SpinEdit1.Value < SpinEdit1.MinValue then SpinEdit1.Value:=

SpinEdit1.MinValue;

end;

procedure TForm_Options.SpinEdit2Change(Sender: TObject);

begin

SpinEdit3.MinValue := SpinEdit2.Value;

if SpinEdit3.Value < SpinEdit3.MinValue then SpinEdit3.Value:=

SpinEdit3.MinValue;

end;

procedure TForm_Options.CheckBox_LinkClick(Sender: TObject);

begin

if CheckBox_Link.State = cbChecked then Form_Main.BindGrids

else Form_Main.UnBindGrids

end;

end.

Страницы: 1, 2, 3