Алгоритм компактного хранения и решения СЛАУ высокого порядка - (реферат)
p>if (sscanf(TextBuffer, "%ld %lf %lf %lf", &NumPoints, &tx, &ty, &tz) ! = 4) continue;X[Current] = tx;
Y[Current] = ty;
Z[Current] = tz;
if (++Current == 999)
{
Vector t1 = X,
t2 = Y,
t3 = Z;
X. ReSize(2 * X. Size());
Y. ReSize(2 * X. Size());
Z. ReSize(2 * X. Size());
for (DWORD i = 0; i < Current; i++)
{
X[i] = t1[i];
Y[i] = t2[i];
Z[i] = t3[i];
}
}
if (Current % 100 == 0)
printf("Line: %ld\r", Current);
}
fclose(in1);
printf(" \r");
NumPoints = Current;
Current = 0L;
while (! feof(in2))
{
if (fgets(TextBuffer, 1000, in2) == NULL)
{
if (feof(in2)) break;
printf("Unable read file %s", fn2);
fclose(in2);
return false;
}
if (sscanf(TextBuffer, "%d %d %d %d %d %ld %ld %ld %ld %ld %ld %ld %ld", &tmp, &tmp, &tmp, &tmp, &tmp,
&FE[Current][0], &FE[Current][1], &FE[Current][3], &FE[Current][2], &FE[Current][4], &FE[Current][5], &FE[Current][7], &FE[Current][6]) ! = 13) continue;
if (++Current == 999)
{
Matrix t = FE;
FE. ReSize(2 * FE. Size1(), 10);
for (DWORD i = 0; i < Current; i++)
for (DWORD j = 0; j < 10; j++)
FE[i][j] = t[i][j];
}
if (Current % 100 == 0)
printf("Line: %ld\r", Current); }
NumFE = Current;
for (DWORD i = 0; i < NumFE; i++)
for (DWORD j = 0; j < 10; j++)
FE[i][j]--;
printf(" \r");
return true; }
ПРИЛОЖЕНИЕ 2.
Исходный текст программы, реализующей алгоритм компактного хранения и решения СЛАУ высокого порядка.
#include "matrix. h"
class RVector
{
private:
Vector Buffer;
public:
RVector(void) {}
~RVector() {}
RVector(DWORD Size) { Buffer. ReSize(Size); }
RVector(RVector& right) { Buffer = right. Buffer; }
RVector(Vector& right) { Buffer = right; }
DWORD Size(void) { return Buffer. Size(); }
void ReSize(DWORD Size) { Buffer. ReSize(Size); }
double& operator [] (DWORD i) { return Buffer[i]; }
RVector& operator = (RVector& right) { Buffer = right. Buffer; return *this; } RVector& operator = (Vector& right) { Buffer = right; return *this; } void Sub(RVector&);
void Sub(RVector&, double);
void Mul(double);
void Add(RVector&);
friend double Norm(RVector&, RVector&);
};
class TSMatrix
{
private:
Vector Right;
Vector* Array;
Vector* Links;
uint Dim;
DWORD Size;
public:
TSMatrix(void) { Size = 0; Dim = 0; Array = NULL; Links = NULL; } TSMatrix(Vector*, DWORD, uint);
~TSMatrix(void) { if (Array) delete [] Array; }
Vector& GetRight(void) { return Right; }
DWORD GetSize(void) { return Size; }
uint GetDim(void) { return Dim; }
Vector& GetVector(DWORD i) { return Array[i]; }
Vector* GetLinks(void) { return Links; }
void SetLinks(Vector* l) { Links = l; }
void Add(Matrix&, Vector&);
void Add(DWORD I, DWORD L, DWORD J, DWORD K, double v)
{
DWORD Row = I,
Col = L * Links[I]. Size() * Dim + Find(I, J) * Dim + K;
Array[Row][Col] += v;
}
void Add(DWORD I, double v)
{
Right[I] += v;
}
DWORD Find(DWORD, DWORD);
void Restore(Matrix&);
void Set(DWORD, DWORD, double, bool);
void Set(DWORD Index1, DWORD Index2, double value)
{
DWORD I = Index1 / Dim,
L = Index1 % Dim,
J = Index2 / Dim,
K = Index2 % Dim,
Pos = Find(I, J),
Row = I,
Col;
if (Pos == DWORD(-1)) return;
Col = L * Links[I]. Size() * Dim + Find(I, J) * Dim + K;
Array[Row][Col] = value;
}
bool Get(DWORD Index1, DWORD Index2, double& value)
{
DWORD I = Index1 / Dim,
L = Index1 % Dim,
J = Index2 / Dim,
K = Index2 % Dim,
Pos = Find(I, J),
Row = I,
Col;
value = 0;
if (Pos == DWORD(-1)) return false;
Col = L * Links[I]. Size() * Dim + Find(I, J) * Dim + K;
value = Array[Row][Col];
return true;
}
void Mul(RVector&, RVector&);
double Mul(DWORD, RVector&);
void write(ofstream&);
void read(ifstream&);
};
class RMatrix
{
private:
Vector Buffer;
DWORD size;
public:
RMatrix(DWORD sz) { size = sz; Buffer. ReSize(size*(size + 1)*0. 5); } ~RMatrix() {}
DWORD Size(void) { return size; }
double& Get(DWORD i, DWORD j) { return Buffer[(2*size + 1 - i)*0. 5*i + j - i]; }
};
//************************
#include "smatrix. h"
double Norm(RVector& Left, RVector& Right)
{
double Ret = 0;
for (DWORD i = 0; i < Left. Size(); i++)
Ret += Left[i] * Right[i];
return Ret;
}
void RVector: :Sub(RVector& Right)
{
for (DWORD i = 0; i < Size(); i++)
(*this)[i] -= Right[i];
}
void RVector: :Add(RVector& Right)
{
for (DWORD i = 0; i < Size(); i++)
(*this)[i] += Right[i];
}
void RVector: :Mul(double koff)
{
for (DWORD i = 0; i < Size(); i++)
(*this)[i] *= koff;
}
void RVector: :Sub(RVector& Right, double koff)
{
for (DWORD i = 0; i < Size(); i++)
(*this)[i] -= Right[i]*koff;
}
TSMatrix: :TSMatrix(Vector* links, DWORD size, uint dim) {
Dim = dim;
Links = links;
Size = size;
Right. ReSize(Dim * Size);
Array = new Vector[Size];
for (DWORD i = 0; i < Size; i++)
Array[i]. ReSize(Links[i]. Size() * Dim * Dim);
}
void TSMatrix: :Add(Matrix& FEMatr, Vector& FE) {
double Res;
DWORD RRow;
for (DWORD i = 0L; i < FE. Size(); i++)
for (DWORD l = 0L; l < Dim; l++)
for (DWORD j = 0L; j < FE. Size(); j++)
for (DWORD k = 0L; k < Dim; k++)
{
Res = FEMatr[i * Dim + l][j * Dim + k];
if (Res) Add(FE[i], l, FE[j], k, Res);
}
for (DWORD i = 0L; i < FE. Size(); i++)
for (DWORD l = 0L; l < Dim; l++)
{
RRow = FE[UINT(i % (FE. Size()))] * Dim + l;
Res = FEMatr[i * Dim + l][FEMatr. Size1()];
if (Res) Add(RRow, Res);
}
}
DWORD TSMatrix: :Find(DWORD I, DWORD J)
{
DWORD i;
for (i = 0; i < Links[I]. Size(); i++)
if (Links[I][i] == J) return i;
return DWORD(-1);
}
void TSMatrix: :Restore(Matrix& Matr)
{
DWORD i,
j,
NRow,
NPoint,
NLink,
Pos;
Matr. ReSize(Size * Dim, Size * Dim + 1);
for (i = 0; i < Size; i++)
for (j = 0; j < Array[i]. Size(); j++)
{
NRow = j / (Array[i]. Size() / Dim); // Number of row NPoint = (j - NRow * (Array[i]. Size() / Dim)) / Dim; // Number of points NLink = j % Dim; // Number of link Pos = Links[i][NPoint];
Matr[i * Dim + NRow][Pos * Dim + NLink] = Array[i][j];
}
for (i = 0; i < Right. Size(); i++) Matr[i][Matr. Size1()] = Right[i]; }
void TSMatrix: :Set(DWORD Index, DWORD Position, double Value, bool Case) {
DWORD Row = Index,
Col = Position * Links[Index]. Size() * Dim + Find(Index, Index) * Dim + Position,
i;
double koff = Array[Row][Col],
val;
if (! Case)
Right[Dim * Index + Position] = Value;
else
{
Right[Index * Dim + Position] = Value * koff;
for (i = 0L; i < Size * Dim; i++)
if (i ! = Index * Dim + Position)
{
Set(Index * Dim + Position, i, 0);
Set(i, Index * Dim + Position, 0);
if (Get(i, Index * Dim + Position, val))
Right[i] -= val * Value;
}
}
}
void TSMatrix: :Mul(RVector& Arr, RVector& Res)
{
DWORD i,
j,
NRow,
NPoint,
NLink,
Pos;
Res. ReSize(Arr. Size());
for (i = 0; i < Size; i++)
for (j = 0; j < Array[i]. Size(); j++)
{
NRow = j / (Array[i]. Size() / Dim);
NPoint = (j - NRow * (Array[i]. Size() / Dim)) / Dim;
NLink = j % Dim;
Pos = Links[i][NPoint];
Res[i * Dim + NRow] += Arr[Pos * Dim + NLink] * Array[i][j]; }
}
double TSMatrix: :Mul(DWORD Index, RVector& Arr)
{
DWORD j,
I = Index / Dim,
L = Index % Dim,
Start = L * (Array[I]. Size() / Dim),
Stop = Start + (Array[I]. Size() / Dim),
NRow,
NPoint,
NLink,
Pos;
double Res = 0;
for (j = Start; j < Stop; j++)
{
NRow = j / (Array[I]. Size() / Dim);
NPoint = (j - NRow * (Array[I]. Size() / Dim)) / Dim;
NLink = j % Dim;
Pos = Links[I][NPoint];
Res += Arr[Pos * Dim + NLink] * Array[I][j];
}
return Res;
}
void TSMatrix: :write(ofstream& Out)
{
DWORD ColSize;
Out. write((char*)&(Dim), sizeof(DWORD));
Out. write((char*)&(Size), sizeof(DWORD));
for (DWORD i = 0; i < Size; i++)
{
ColSize = Array[i]. Size();
Out. write((char*)&(ColSize), sizeof(DWORD));
for (DWORD j = 0; j < ColSize; j++)
Out. write((char*)&(Array[i][j]), sizeof(double));
}
for (DWORD i = 0; i < Size * Dim; i++)
Out. write((char*)&(Right[i]), sizeof(double));
}
void TSMatrix: :read(ifstream& In)
{
DWORD ColSize;
In. read((char*)&(Dim), sizeof(DWORD));
In. read((char*)&(Size), sizeof(DWORD));
if (Array) delete [] Array;
Array = new Vector[Size];
Right. ReSize(Size * Dim);
for (DWORD i = 0; i < Size; i++)
{
In. read((char*)&(ColSize), sizeof(DWORD));
Array[i]. ReSize(ColSize);
for (DWORD j = 0; j < ColSize; j++)
In. read((char*)&(Array[i][j]), sizeof(double));
}
for (DWORD i = 0; i < Size * Dim; i++)
In. read((char*)&(Right[i]), sizeof(double));
}
