| Цитата(Амортизатор) | | Есть ли в стандартной библиотеке средства для решения сиситемы лин. уравнений, заданной матрицей коэффициентов и свободных членов. |
В math.h? Там только калькуляторные функции.
| Цитата(Амортизатор) | | Приведите, п-ста, функцию с опис. параметров, если таковая имеется. |
Существуют специальные алгоритмы. Например:
Метод квадратного корня
| Код |
#include <stdio.h>
#include <math.h>
#define N 4 // size of matrix
#define N1 N+1
float matrix[N][N1]=
{{10.9,1.2,2.1,0.9,23.2},
{1.2,11.2,1.5,2.5,38.1},
{2.7,1.5,9.8,1.3,40.3},
{0.9,2.5,1.3,12.1,58.2}
};
void ShowMatrix(void)
{
for (int i=0;i<N;i++)
{
for (int j=0;j<N;j++)
printf("%+f*x%d",matrix[i][j],i+1);
printf("=%f\n",matrix[i][N]);
}
}
void main(void)
{
// Variables declaration
register char i,j,k;
float s[N][N],x[N],y[N],d[N],sum;
// Printing given matrix
ShowMatrix();
// Building matrixes S and D
for (j=0;j<N;j++)
for (i=0;i<=j;i++)
{
sum=matrix[i][j];
for (k=0;k<i;k++)
sum-=s[k][i]*d[k]*s[k][j];
if (i==j)
{
d[i]=(sum>0?1:-1);
s[i][j]=sqrt(fabs(sum));
}
else s[i][j]=sum/(d[i]*s[i][i]);
}
// Solving equation Gy=b (G: g[I][j]=s[j][I]*d[j])
for (i=0;i<N;i++)
{
y[i]=matrix[i][N];
for (j=0;j<i;j++) y[i]-=s[j][i]*d[j]*y[j];
y[i]/=s[i][i]*d[i];
}
// Solving equation Sx=y
for (i=N-1;i>=0;i--)
{
x[i]=y[i];
for (j=i+1;j<N;j++) x[i]-=s[i][j]*x[j];
x[i]/=s[i][i];
}
// Printing solution
printf("\nSolution:\n");
for (i=0;i<N;i++)
printf("x%d=%f\n",i+1,x[i]);
}
|
Метод Некрасова
| Код |
#include <stdio.h>
#include <math.h>
#define N 4 // size of matrix
#define N1 N+1
float matrix[N][N1]=
{{10.9,1.2,2.1,0.9,23.2},
{1.2,11.2,1.5,2.5,38.1},
{2.7,1.5,9.8,1.3,40.3},
{0.9,2.5,1.3,12.1,58.2}
};
void ShowMatrix(void)
{
for (int i=0;i<N;i++)
{
for (int j=0;j<N;j++)
printf("%+f*x%d",matrix[i][j],i+1);
printf("=%f\n",matrix[i][N]);
}
}
// this function calculates ||x||
float norma(float *x)
{
float sum=0;
for (unsigned i=0;i<N;i++) sum+=fabs(*(x+i));
return sum;
}
void main(void)
{
// Variables declaration
register char i,j;
float x_current[N],x_next[N],sum1,sum2,epsilon;
epsilon=1e-7;
// Printing given matrix
ShowMatrix();
// Solving
do
{
for (i=0;i<N;i++) x_current[i]=x_next[i];
for (i=0;i<N;i++)
{
sum1=0;sum2=0;
for(j=0;j<i;j++) sum1+=matrix[i][j]*x_next[j];
for(j=i+1;j<N;j++) sum2+=matrix[i][j]*x_current[j];
x_next[i]=(matrix[i][N]-sum1-sum2)/matrix[i][i];
}
} while (fabs(norma(x_current)-norma(x_next))>epsilon);
// Printing solution
printf("\nSolution:\n");
for (i=0;i<N;i++)
printf("x%d=%f\n",i+1,x_current[i]);
}
|
Метод Гаусса
| Код |
#include <stdio.h>
#include <math.h>
#define N 4 // size of matrix
#define N1 N+1
float matrix[N][N1]=
{{2,5,4,1,28},
{1,3,2,1,17},
{2,10,9,7,77},
{3,8,9,2,54}
};
void ShowMatrix(void)
{
for (int i=0;i<N;i++)
{
for (int j=0;j<N;j++)
printf("%+f*x%d",matrix[i][j],i+1);
printf("=%f\n",matrix[i][N]);
}
}
void main(void)
{
// Variables declaration
float tmp,x[N1];
register short int i,j,k;
// Printing given matrix
ShowMatrix();
// Main loop
for (i=0;i<N;i++)
{
// Excluding variable x[i] from equations
tmp=matrix[i][i];
for (j=N;j>=i;j--) matrix[i][j]/=tmp;
for (j=i+1;j<N;j++)
{
tmp=matrix[j][i];
for (k=N;k>=i;k--)
matrix[j][k]-=tmp*matrix[i][k];
}
}
// Calculating vector x
x[N-1]=matrix[N-1][N];
for (i=N-2;i>=0;i--)
{
x[i]=matrix[i][N];
for (j=i+1;j<N;j++) x[i]-=matrix[i][j]*x[j];
}
// Printing solution
printf("\nSolution:\n");
for (i=0;i<N;i++)
printf("x%d=%f\n",i+1,x[i]);
}
|
Метод итераций
| Код |
#include <stdio.h>
#include <math.h>
#define N 3 // size of matrix
#define N1 N+1
float matrix[N][N1]=
{{14.38,-2.41,1.39,5.86},
{1.84,25.36,-3.31,-2.28},
{2.46,-3.49,16.37,4.47}
};
int maxiterations=10; // maximum number of iterations
float epsilon=0.0001; // required accuracy
void ShowMatrix(void)
{
for (int i=0;i<N;i++)
{
for (int j=0;j<N;j++)
printf("%+f*x%d",matrix[i][j],i+1);
printf("=%f\n",matrix[i][N]);
}
}
void main(void)
{
// Variables declaration
float x[N],y[N],t;
register short int i,j,k;
int iterations=0;
// Printing given matrix
ShowMatrix();
// setting first iteration of vector X
for (i=0;i<N;i++) x[i]=matrix[i][N];
do
{
for (i=0;i<N;i++)
{
t=-matrix[i][N];
for (j=0;j<N;j++)
t+=matrix[i][j]*x[j];
y[i]=(-t+matrix[i][i]*x[i])/matrix[i][i];
}
iterations++;
k=0;
// checking solution
while (fabs(x[k]-y[k])<epsilon && k<N) k++;
// new iteration becomes old
for(i=0;i<N;i++) x[i]=y[i];
} while (k!=N && iterations<maxiterations);
if (iterations==maxiterations)
{
printf("Iterations are very slow...");
} else
{
// Printing solution
printf("\nSolution:\n");
for (i=0;i<N;i++)
printf("x%d=%f\n",i+1,x[i]);
printf("%d iterations were made",iterations);
}
}
|
Метод Монте-Карло
| Код |
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#define N 3 // size of matrix
#define N1 N+1
#define M 1000 // number of tries
float matrix[N][N1]=
{{3,1,-1,-2},
{1,4,-2,-3},
{2,-1,5,-15}
};
void ShowMatrix(void)
{
for (int i=0;i<N;i++)
{
printf("x%d=",i+1);
for (int j=0;j<N;j++)
printf("%+f*x%d",matrix[i][j],j+1);
printf("%+f\n",matrix[i][N]);
}
}
void main(void)
{
// Variables declaration
float b[N][N],w[N],y,c,x[N],l;
register short int i,j,v;
unsigned char finish;
int row, // currently analysed row
tries; // number of tries to calculate true root by generating randoms
// Transforming matrix
for (i=0;i<N;i++)
{
l=0;
for (j=0;j<N1;j++) l+=fabs(matrix[i][j]);
v=(matrix[i][i]>0?-1:1);
for (j=0;j<N1;j++) matrix[i][j]=v*(matrix[i][j])/l+(i==j?1:0);
}
// Printing transformed matrix
ShowMatrix();
// filling matrix B
for (i=0;i<N;i++)
{
b[i][0]=fabs(matrix[i][0]);
for(j=1;j<N;j++)
b[i][j]=b[i][j-1]+fabs(matrix[i][j]);
}
// filling matrix w
for(i=0;i<N;i++) w[i]=matrix[i][N]/(1-b[i][N-1]);
// for each equation in the system
for(i=0;i<N;i++)
{
tries=0;
row=i;
y=0;
v=1;
while (tries<M)
{
// generating random number
c=random(32767)/32767.0;
j=N-1;
finish=0;
while(j>=0 && !finish)
{
if (c>b[row][j])
{
if (j==N-1)
{
tries++;
y+=v*w[row];
row=i;v=1;
}
else
{
v*=(matrix[row][j+1]>=0?1:-1);
row=j+1;
}
finish=1;
}
j--;
}
if (!finish)
{
v*=(matrix[row][j+1]>=0?1:-1);
row=0;
}
}
// calculating root
x[i]=y/M;
}
// Printing solution
printf("\nSolution:\n");
for (i=0;i<N;i++)
printf("x%d=%f\n",i+1,x[i]);
}
|
Метод ортогонализации
| Код |
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#define N 3 // size of matrix
#define N1 N+1
double matrix[N1][N1]=
{{3,1,-1,-2},
{1,4,-2,-3},
{2,-1,5,-15}
};
void ShowMatrix(void)
{
for (int i=0;i<N;i++)
{
for (int j=0;j<N;j++)
printf("%+f*x%d",matrix[i][j],i+1);
printf("%+f=0\n",matrix[i][N]);
}
}
void main(void)
{
// Variables declaration
double c[N],x[N],s;
int i,j,k,l,m;
ShowMatrix();
// expanding system by vector (0,0,0,...,0,1)
for (i=0;i<N;i++) matrix[N][i]=0;
matrix[N][N]=1;
// for all equations
for (i=0;i<N1;i++)
{
if (i>0)
{
k=0;
// make some iterations to increase accuracy of calculations
while (k<=3)
{
for (l=0;l<i;l++)
{
c[l]=0;
for(m=0;m<N1;m++) c[l]+=matrix[l][m]*matrix[i][m];
}
for (j=0;j<N1;j++)
{
s=0;
for(l=0;l<i;l++) s+=c[l]*matrix[l][j];
matrix[i][j]-=s;
}
k++;
}
}
// normalizing vector
s=0;
for (k=0;k<N1;k++) s+=matrix[i][k]*matrix[i][k];
s=sqrt(s);
for (j=0;j<N1;j++) matrix[i][j]/=s;
}
s=matrix[N][N];
for (j=0;j<N;j++) x[j]=matrix[N][j]/s;
// Printing solution
printf("\nSolution:\n");
for (i=0;i<N;i++)
printf("x%d=%f\n",i+1,x[i]);
}
|
|
