Приложение 2. Код программы.
//---------------------------------------------------------------------------
#include <vcl.h>
#include <time.h>
#include <math.h>
#include <stdio.h>
#pragma hdrstop
#include "Unit1.h"
//---------------------------------------------------------------------------
#pragma package(smart_init)
#pragma resource "*.dfm"
TForm1 *Form1;
int n,m;
float epsilon; //Точность
double**a; //Матрица
double *f; //Столбец правых частей
double *x; //Решение
double *ax;
double**b;
bool ch=0;
double *tm;
//---------------------------------------------------------------------------
__fastcall TForm1::TForm1(TComponent* Owner)
: TForm(Owner)
{
sg[0]=StringGrid2;
sg[1]=StringGrid3;
sg[2]=StringGrid4;
sg[3]=StringGrid5;
sg[4]=StringGrid6;
sg[5]=StringGrid7;
sg[6]=StringGrid8;
}
//---------------------------------------------------------------------------
//Матрично-векторное умножение
void mv_mult(double **a,double *x, double *ax)
{
for (int i=0;i<n;i++)
{
ax[i]=0;
for (int j=0;j<n;j++)
{
ax[i]+=a[i][j]*x[j];
}
}
}
//---------------------------------------------------------------------------
//Векторное вычитание
void vv_substr(double *ax, double *f, double *r)
{
for (int i=0;i<n;i++)
{
r[i]=ax[i]-f[i];
}
}
//---------------------------------------------------------------------------
//Векторное сложение
void vv_addit(double *ax, double *f, double *r)
{
for (int i=0;i<n;i++)
{
r[i]=ax[i]+f[i];
}
}
//---------------------------------------------------------------------------
//Норма вектора
double norma(double * vec)
{
double nr=0;
for (int i=0;i<n;i++)
{
nr+=vec[i]*vec[i];
}
return (sqrt(nr));
}
//---------------------------------------------------------------------------
//Скалярное умножение векторов
double vv_mult(double *a,double *x)
{
double ax=0;
for (int i=0;i<n;i++)
{
ax+=a[i]*x[i];
}
return(ax);
}
//---------------------------------------------------------------------------
//Проверка на заданную точность нормы
bool check(double *a)
{
double mod;
mod=norma(a);
if (mod<epsilon) return(1);
else return(0);
}
//---------------------------------------------------------------------------
//Построение матрицы B для метода минимальных поправок
void mkB()
{
ch=1;
b=new double* [n];
for (int i=0; i<n; i++)
{
b[i]=new double[n];
for (int j=0; j<n; j++)
if (i!=j) b[i][j]=0;
else b[i][j]=a[i][j]-300;
}
}
//---------------------------------------------------------------------------
//Заполнение матрицы
void __fastcall TForm1::Button5Click(TObject *Sender)
{
Button1->Enabled=true;
Button2->Enabled=true;
Button3->Enabled=true;
Button4->Enabled=true;
CheckBox2->Enabled=true;
time_t t;
time(&t); srand((unsigned)t);
n=StrToInt(Edit1->Text);
m=StrToInt(Edit3->Text);
epsilon=StrToFloat(Edit2->Text);
StringGrid1->ColCount=n;
StringGrid1->RowCount=n;
for (int i=0;i<6;i++)
sg[i]->RowCount=n;
f=new double [n];
if (CheckBox2->Checked==false){
a=new double* [n];
tm=new double [n];
x=new double [n];
for (int i=0;i<n;i++)
a[i]=new double [n];
for (int i=0;i<n;i++)
{
double b=rand()%200-100;
double c=rand()%100+1;
x[i]=b/c;
sg[1]->Cells[0][i]=FloatToStr(x[i]);
for (int j=0;j<=i;j++)
{
double b=rand()%200-100;
double c=rand()%100+1;
a[i][j]=b/c;
a[j][i]=a[i][j];
if ((i==j)&&(CheckBox1->Checked==true)) a[i][j]+=m;
if (i==j) tm[i]=a[i][j];
StringGrid1->Cells[j][i]=FloatToStr(a[i][j]);
StringGrid1->Cells[i][j]=FloatToStr(a[j][i]);
}
}}
else
for (int i=0;i<n;i++)
{
a[i][i]=tm[i]+(double)m;
tm[i]+=m;
StringGrid1->Cells[i][i]=FloatToStr(a[i][i]);
StringGrid1->Cells[i][i]=FloatToStr(a[i][i]);
}
mv_mult(a,x,f);
for(int i=0;i<n;i++)
sg[0]->Cells[0][i]=FloatToStr(f[i]);
}
//---------------------------------------------------------------------------
void __fastcall TForm1::FormClose(TObject *Sender, TCloseAction &Action)
{
for (int i=0;i<n;i++)
{
delete [] a[i];
if (ch==1) delete [] b[i];
}
}
//---------------------------------------------------------------------------
//Метод минимальных невязок
void __fastcall TForm1::Button1Click(TObject *Sender)
{
double *rk; //невязка
double *t; //вспомогательный вектор
double *xk;
double *xk1;
double tau,sm,nrm;
ax=new double [n];
xk=new double [n];
xk1=new double [n];
rk=new double [n];
t=new double [n];
for (int i=0;i<n;i++)//Задаём начальное приближение
xk[i]=0;
int k=0;
mv_mult(a,xk,ax); //ax=A*xk
vv_substr(ax,f,rk); //Вычисляем невязку: rk=A*xk-f
bool c=check(rk);
while (c==0)
{
mv_mult(a,rk,ax); //A*rk
sm=vv_mult(ax,rk); //Считаем скалярное произведение: (A*rk,rk)
nrm=norma(ax);
tau=sm/(nrm*nrm);
for (int i=0;i<n;i++)
{
t[i]=rk[i]*tau;
ax[i]=ax[i]*tau;
}
vv_substr(xk,t,xk1);
for (int i=0;i<n;i++)
{
xk[i]=xk1[i];
}
vv_substr(rk,ax,t);
for (int i=0;i<n;i++)
{
rk[i]=t[i];
}
c=check(rk);k++;
}
for(int i=0;i<n;i++)
sg[2]->Cells[0][i]=FloatToStr(xk[i]);
sg[6]->Cells[0][0]="Метод мин.невязок";
sg[6]->Cells[1][0]=IntToStr(k);;
delete []rk; delete [] xk; delete []xk1; delete [] t;
}
//---------------------------------------------------------------------------
//Метод минимальных поправок
void __fastcall TForm1::Button2Click(TObject *Sender)
{
double *rk; //невязка
double *wk; //поправка
double *vk;
double *xk;
double *xk1;
double *t;
double tau,sm,smd;
ax=new double [n];
xk=new double [n];
xk1=new double [n];
rk=new double [n];
wk=new double [n];
vk=new double [n];
t=new double [n];
for (int i=0;i<n;i++)//Задаём начальное приближение
xk[i]=0;
bool c=0;
mkB();int k=0;
mv_mult(a,xk,ax); //ax=A*xk
vv_substr(ax,f,rk); //Вычисляем невязку: rk=A*xk-f
mv_mult(b,rk,wk); //Вычисляем поправку wk=B*rk
while (c==0)
{
mv_mult(a,wk,ax); //ax=A*wk
mv_mult(b,ax,vk); //vk=B*A*wk
sm=vv_mult(ax,wk); //Считаем скалярное произведение: (A*wk,wk)
smd=vv_mult(vk,ax);
tau=sm/smd;
for (int i=0;i<n;i++)
t[i]=wk[i]*tau;
vv_substr(xk,t,xk1);
for (int i=0;i<n;i++)
{
xk[i]=xk1[i];
vk[i]=vk[i]*tau;
}
vv_substr(wk,vk,t);
for (int i=0;i<n;i++)
{
wk[i]=t[i];
}
c=check(wk);k++;
}
for(int i=0;i<n;i++)
sg[3]->Cells[0][i]=FloatToStr(xk[i]);
sg[6]->Cells[0][1]="Метод мин.поправок";
sg[6]->Cells[1][1]=IntToStr(k);
delete []rk; delete [] xk; delete []xk1;delete []wk; delete [] vk; delete [] t;
}
//---------------------------------------------------------------------------
//Метод скорейшего спуска
void __fastcall TForm1::Button3Click(TObject *Sender)
{
double *rk; //невязка
double *xk;
double *xk1;
double *t;
double tau,sm,smd;
ax=new double [n];
xk=new double [n];
xk1=new double [n];
rk=new double [n];
t=new double [n];
for (int i=0;i<n;i++)//Задаём начальное приближение
xk[i]=0;
bool c=0;
int k=0;
mv_mult(a,xk,ax); //ax=A*xk
vv_substr(ax,f,rk); //Вычисляем невязку: rk=A*xk-f
while (c==0)
{
mv_mult(a,rk,ax); //A*rk
sm=vv_mult(rk,rk); //Считаем скалярное произведение: (rk,rk)
smd=vv_mult(ax,rk); //Считаем скалярное произведение: (A*rk,rk)
tau=sm/smd;
for (int i=0;i<n;i++)
{
t[i]=rk[i]*tau;
ax[i]=ax[i]*tau;
}
vv_substr(xk,t,xk1);
for (int i=0;i<n;i++)
xk[i]=xk1[i];
vv_substr(rk,ax,t);
for (int i=0;i<n;i++)
{
rk[i]=t[i];
}
c=check(rk);k++;
}
for(int i=0;i<n;i++)
sg[4]->Cells[0][i]=FloatToStr(xk[i]);
sg[6]->Cells[0][2]="Метод наиск.спуска";
sg[6]->Cells[1][2]=IntToStr(k);
delete []rk; delete [] xk; delete []xk1;
}
//---------------------------------------------------------------------------
//Метод сопряженных градиентов
void __fastcall TForm1::Button4Click(TObject *Sender)
{
double *rk; //невязка
double *rk1;
double *pk;
double *xk;
double *xk1;
double *t;
double alpha,beta,sm,smd;
ax=new double [n];
xk=new double [n];
xk1=new double [n];
rk=new double [n];
rk1=new double [n];
pk=new double [n];
t=new double [n];
for (int i=0;i<n;i++)//Задаём начальное приближение
xk[i]=0;
bool c=0;
int k=0;
mv_mult(a,xk,ax); //ax=A*xk
vv_substr(ax,f,rk); //Вычисляем невязку: r0=A*x0-f
vv_substr(ax,f,pk); //p0=r0
while (c==0)
{
mv_mult(a,pk,ax); //A*pk
sm=vv_mult(rk,rk); //Считаем скалярное произведение: (rk,rk)
smd=vv_mult(ax,pk); //Считаем скалярное произведение: (A*pk,pk)
alpha=sm/smd;
for (int i=0;i<n;i++)
{
t[i]=pk[i]*alpha;
ax[i]=ax[i]*alpha;
}
vv_addit(xk,t,xk1);
for (int i=0;i<n;i++)
xk[i]=xk1[i];
vv_substr(rk,ax,rk1);
sm=vv_mult(rk1,rk1);
smd=vv_mult(rk,rk);
beta=sm/smd;
c=check(rk);
for (int i=0;i<n;i++)
pk[i]=pk[i]*beta;
vv_addit(rk1,pk,pk);
for (int i=0;i<n;i++)
rk[i]=rk1[i];
k++;
}
for(int i=0;i<n;i++)
{
xk[i]=-1*xk[i];
sg[5]->Cells[0][i]=FloatToStr(xk[i]);
}
sg[6]->Cells[0][3]="Метод сопр.градиентов";
sg[6]->Cells[1][3]=IntToStr(k);
delete []rk; delete [] xk; delete []xk1;delete []rk1; delete [] pk;
}
//---------------------------------------------------------------------------
void __fastcall TForm1::Edit1Change(TObject *Sender)
{
CheckBox2->Enabled=false;
CheckBox2->Checked=false;
}
//---------------------------------------------------------------------------
//Чтение из файла
void __fastcall TForm1::Button6Click(TObject *Sender)
{
FILE *matr;
FILE *pch;
String s=ComboBox1->Text,p=s+"p";
Button1->Enabled=true;Button2->Enabled=true;Button3->Enabled=true;
Button4->Enabled=true;CheckBox2->Enabled=true;
epsilon=StrToFloat(Edit2->Text);
matr=fopen(s.c_str(),"rt");
fscanf(matr,"%d\n",&n);
pch=fopen(p.c_str(),"rt");
int tn=StrToInt(Edit1->Text);
if (tn>n) MessageBox(0,"Размерность матрицы в файле меньше заданной","Ошибка",MB_OK);
else {
n=tn;
f=new double [n];
a=new double* [n];
StringGrid1->ColCount=n;
StringGrid1->RowCount=n;
for (int i=0;i<6;i++)
{
sg[i]->RowCount=n;
sg[1]->Cells[0][i]="";
}
for (int i=0;i<n;i++)
a[i]=new double [n];
for (int i=0;i<n;i++)
{
float r;
fscanf(pch,"%f",&r);
f[i]=r;
sg[0]->Cells[0][i]=FloatToStr(f[i]);
for (int j=0;j<n;j++)
{
fscanf(matr,"%f",&r);
a[i][j]=r;
StringGrid1->Cells[j][i]=FloatToStr(a[i][j]);
}
}}
fclose(matr);fclose(pch);
}
//---------------------------------------------------------------------------
//Запись в файл
void __fastcall TForm1::Button7Click(TObject *Sender)
{
FILE *matr;
FILE *pch;
String s=ComboBox2->Text,p=s+"p";
matr=fopen(s.c_str(),"w");
fprintf(matr,"%d\n",n);
pch=fopen(p.c_str(),"w");
for (int i=0;i<n;i++)
{
fprintf(pch,"%f\n",f[i]);
for (int j=0;j<n;j++)
fprintf(matr,"%f ",a[i][j]);
fprintf(matr,"\n");
}
fclose(matr);fclose(pch);
}