Нашел программу для решения систем алгебраических уравнений.
Работает, но на маленьких матрицах. Запустил на матрице порядка 3000 - не работает. в чем ошибка не пойму
| Код |
#include <stdio.h>
#include <process.h>
#include <sstream>
#include <fstream>
#include <iostream>
using namespace std;
float **a, *b, *x;
int N;
char filename[256];
FILE* InFile=NULL;
void count_num_lines(){
//count number of lines in input file - number of equations
int nelf=0; //non empty line flag
do{
nelf = 0;
while(fgetc(InFile)!='\n' && !feof(InFile)) nelf=1;
if(nelf) N++;
}while(!feof(InFile));
}
void freematrix(){
//free memory for matrixes
int i;
for(i=0; i<N; i++){
delete [] a[i];
}
delete [] a;
delete [] b;
delete [] x;
}
void allocmatrix(){
//allocate memory for matrixes
int i,j;
x = new float[N]();
b = new float[N]();
a = new float*[N]();
if(x==NULL || b==NULL || a==NULL){
printf("\nNot enough memory to allocate for %d equations.\n", N);
exit(-1);
}
for(i=0; i<N; i++){
a[i] = new float[N];
if(a[i]==NULL){
printf("\nNot enough memory to allocate for %d equations.\n", N);
}
}
void readmatrix(){
int i=0,j=0;
//read matrixes a and b from input file
for(i=0; i<N; i++){
for(j=0; j<N; j++){
fscanf(InFile, "%f", &a[i][j]);
}
fscanf(InFile, "%f", &b[i]);
}
}
//void printmatrix(){
// //print matrix "a"
// int i=0,j=0;
// printf("\n");
// for(i=0; i<N; i++){
// for(j=0; j<N; j++){
// printf(" %+f*X%d", a[i][j], j);
// }
// printf(" =%f\n", b[i]);
// }
//}
//void testsolve(){
// //test that ax=b
// int i=0,j=0;
// printf("\n");
// for(i=0; i<N; i++){
// float s = 0;
// for(j=0; j<N; j++){
// s += a[i][j]*x[j];
// }
// printf("%f\t%f\n", s, b[i]);
// }
//}
void printresult(){
int i=0;
ofstream ofsx("x.txt");
printf("\n");
printf("Result\n");
for(i=0; i<N; i++){
ofsx << x[i] << endl;
//printf("X%d = %f\n", i, x[i]);
}
}
void diagonal(){
int i, j, k;
float temp=0;
for(i=0; i<N; i++){
if(a[i][i]==0){
for(j=0; j<N; j++){
if(j==i) continue;
if(a[j][i] !=0 && a[i][j]!=0){
for(k=0; k<N; k++){
temp = a[j][k];
a[j][k] = a[i][k];
a[i][k] = temp;
}
temp = b[j];
b[j] = b[i];
b[i] = temp;
break;
}
}
}
}
}
void cls(){
for(int i=0; i<25; i++) printf("\n");
}
void main(){
int i=0,j=0, k=0;
cls();
do{
printf("\nInput filename: ");
scanf("%s", filename);
InFile = fopen("3.txt", "rt");
}while(InFile==NULL);
count_num_lines();
allocmatrix();
rewind(InFile);
//read data from file
readmatrix();
//check if there are 0 on main diagonal and exchange rows in that case
diagonal();
fclose(InFile);
//printmatrix();
//process rows
for(k=0; k<N; k++){
for(i=k+1; i<N; i++){
if(a[k][k]==0){
printf("\nSolution is not exist.\n");
return;
}
float M = a[i][k] / a[k][k];
for(j=k; j<N; j++){
a[i][j] -= M * a[k][j];
}
b[i] -= M*b[k];
}
}
//printmatrix();
for(i=N-1; i>=0; i--){
float s = 0;
for(j = i; j<N; j++){
s = s+a[i][j]*x[j];
}
x[i] = (b[i] - s) / a[i][i];
}
InFile = fopen(filename, "rt");
readmatrix();
fclose(InFile);
//printmatrix();
//testsolve();
printresult();
freematrix();
}
|
программа считывает из файла матрицу и столбец свободных членов.
Матрица и столбец, которые я должен решить, прикрепил в топику. |
