谁能帮我设计一个C语言程序,求解方程组的,急啊,哪位高手帮下咯,万分感谢

Question

如题所述

Answer 1

假定你要的是线性方程组，下面的
float *GauseSeidel(float *a,int n)
是高斯赛德尔法求解线性方程组的通用子程序。
N 是 迭代次数极限。

main()里写了调用的例子。

#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define N 100
float *GauseSeidel(float *a,int n)
{
int i,j,nu=0;
float *x,dx;
x=(float *)malloc(n*sizeof(float));
for(i=0;i<=n-1;i++)
x[i]=0.0;
do {
for(i=0;i<=n-1;i++) {
float d=0.0;
for(j=0;j<=n-1;j++)
d+=*(a+i*(n+1)+j)*x[j];
dx=(*(a+i*(n+1)+n)-d)/(*(a+i*(n+1)+i));
x[i]+=dx;
}
if(nu>=N)
{
printf("iter divergence\n");
exit(1);
}
nu++;
}
while(fabs(dx)>1e-6);
return x;
}
void main()
{
int i;
float *x;
float c[12]={5,2,1,8,2,8,-3,21,1,-3,-6,1};
float *GauseSeidel(float *,int);
x=GauseSeidel(c,3);
for(i=0;i<=2;i++)
printf("x[%d]=%f\n",i,x[i]);
getch();
}

Answer 2

解方程组有好几种方法：
一：高斯消去法，
程序段：
#define N 3
#include <iostream.h>
#include <math.h>
void main()
{
double a[N][N+1] = {0.101,2.304,3.555,1.183,
-1.347,3.712,4.623,2.137,
-2.835,1.072,5.643,3.035}; //增广矩阵
double x[N]; //存储解向量
double m; // 中间变量
int i, j, k;

//消元过程
for(k=0;k<N-1;k++){ // k<N-1 ？
if(a[k][k]==0) {
cout << "求解失败!";
break;
}
for(i=k+1;i<N;i++){
m = a[i][k] / a[k][k];
a[i][k] = 0;
for(j=k+1;j<=N;j++) // j<=N？
a[i][j] = a[i][j] - a[k][j] * m;
}
}
// 回代过程
x[N-1] = a[N-1][N] / a[N-1][N-1];
for(k=N-2;k>=0;k--){
for(j=k+1;j<N;j++){
a[k][N] = a[k][N] - a[k][j] * x[j];
}
x[k] = a[k][N] / a[k][k];
}
// 输出结果
for(k=0;k<N;k++)
cout << x[k]<<endl;
}
二：列主元高斯消去算法，
程序段：
#define N 3
#include <iostream.h>
#include <math.h>
void main()
{
double a[N][N+1] = {0.101,2.304,3.555,1.183,
-1.347,3.712,4.623,2.137,
-2.835,1.072,5.643,3.035}; //增广矩阵
double x[N]; //存储解向量
double m; // 中间变量
int i, j, k,s;

//消元过程
for(k=0;k<N-1;k++)
{ // k<N-1 ？
s=k;
for(i=k+1;i<N;i++)
{
if(fabs(a[s][k])<fabs(a[i][k]))
s=i;
}
if(s!=k)
{
for(j=0;j<=N;j++)
{
m=a[k][j];
a[k][j]=a[s][j];
a[s][j]=m;
}
}
if(a[k][k]==0)
{
cout << "求解失败!";
break;
}
for(i=k+1;i<N;i++)
{
m = a[i][k] / a[k][k];
a[i][k] = 0;
for(j=k+1;j<=N;j++) // j<=N？
a[i][j] = a[i][j] - a[k][j] * m;
}
}
// 回代过程
x[N-1] = a[N-1][N] / a[N-1][N-1];
for(k=N-2;k>=0;k--)
{
for(j=k+1;j<N;j++)
{
a[k][N] = a[k][N] - a[k][j] * x[j];
}
x[k] = a[k][N] / a[k][k];
}
// 输出结果
for(k=0;k<N;k++)
cout << x[k]<<endl;
}
三：顺序高斯消去算法相对应的LU分解方法：
程序段：
#include <iostream.h>
#include <math.h>
#define N 4

void main()
{
double a[N][N]={0.3e-15,59.14,3,1,
5.291,-6.13,-1,2,
11.2,9,5,2,
1,2,1,1 }; //存储矩阵A
double x[N],b[N]={59.17,46.78,1,2};
int i,j,k;
//****** 开始LU分解 *************
for(k=0;k<N-1;k++)
{
if(a[k][k]==0)
{
cout << "主元素为零！" << endl;
break;
}
for(i=k+1;i<N;i++)
{
a[i][k] = a[i][k] / a[k][k]; // 矩阵A的严格下三角部分存储L矩阵的严格下三角部分
for(j=k+1;j<N;j++)
{
a[i][j] = a[i][j] - a[i][k] * a[k][j];
} //计算U矩阵
}
}
//********* LU 分解完成 **************
for(i=0;i<N;i++)
{
for(j=0;j<N;j++)
{
cout << a[i][j] << " ";
}
cout << endl;
} // 输出 LU 分解
x[0]=b[0];
for(k=1;k<N;k++)
{
for(j=0;j<k;j++)
{
b[k]=b[k]-a[k][j]*x[j];
}
x[k]=b[k];
}
x[N-1]=x[N-1]/a[N-1][N-1];
for(k=N-2;k>=0;k--)
{
for(j=k+1;j<N;j++)
{
x[k] = x[k] - a[k][j] * x[j];
}
x[k] = x[k] / a[k][k];
}
for(k=0;k<N;k++)
cout << x[k]<<endl;

}
四：列主元高斯消去算法所对应的LU分解方法：
程序段：
#include <iostream.h>
#include <math.h>
#define N 4

void main()
{
double a[N][N]={0.3e-15,59.14,3,1,
5.291,-6.13,-1,2,
11.2,9,5,2,
1,2,1,1 }; //存储矩阵A
double x[N],b[N]={59.17,46.78,1,2};
int p[N]; // 记录分解过程中的行交换!
int i,j,k;

for(k=0;k<N;k++)
p[k] = k+1; // 初始化，相当于单位阵
//****** 开始LU分解 *************
for(k=0;k<N-1;k++)
{
int s,t; //s存储当前列绝对值最大的元素所在的行号,t为中间变量
s=k;
for(i=k+1;i<N;i++)
{
if(fabs(a[i][k])>fabs(a[s][k]))
s=i;
} // 选主元素
if(s!=k)
{
double m;
for(j=0;j<N;j++)// 注意这里要交换k、s两行所有的元素
{
m = a[k][j];
a[k][j] = a[s][j];
a[s][j] = m;
} // 交换行
t= p[s];
p[s] = p[k];
p[k] = t; //记录当前行交换
}

if(a[k][k]==0)
{
cout << "存在为零的主元素！" << endl;
continue;
}// 判定主元素为零，矩阵奇异，方程组解不唯一

for(i=k+1;i<N;i++)
{
a[i][k] = a[i][k] / a[k][k]; // 矩阵A的严格下三角部分存储L矩阵的严格下三角部分
for(j=k+1;j<N;j++)
{
a[i][j] = a[i][j] - a[i][k] * a[k][j];
} //计算U矩阵
}
}
//********* LU 分解完成 **************
for(i=0;i<N;i++)
{
for(j=0;j<N;j++)
{
cout << a[i][j] << " ";
}
cout << endl;
} // 输出 LU 分解
x[0]=b[0];
for(k=1;k<N;k++)
{
for(j=0;j<k;j++)
{
b[k]=b[k]-a[k][j]*x[j];
}
x[k]=b[k];
}
x[N-1]=x[N-1]/a[N-1][N-1];
for(k=N-2;k>=0;k--)
{
for(j=k+1;j<N;j++)
{
x[k] = x[k] - a[k][j] * x[j];
}
x[k] = x[k] / a[k][k];
}
for(k=0;k<N;k++)
cout << x[k]<<endl;
}
还有一些迭代法；
一：雅可比迭代
#define N 4 // 线性方程组的阶数
#include <iostream.h>
#include <math.h>
void main()
{
double a[N][N]={-4,1,1,1,1,-4,1,1,1,1,-4,1,1,1,1,-4}, //系数矩阵
b[N]={1,1,1,1}; //右端常数向量
double x0[N]={1,1,1,1},x[N]; // 迭代初始向量和迭代向量
double e=1e-5; // 精度要求
int M=5000; //最大迭代次数
int i,j,c_M=0;
double sum,current_e;

do{
current_e=0;
for(i=0;i<N;i++)
{
sum = 0;
for(j=0;j<N;j++)
{
if(j!=i)
{
sum = sum + a[i][j] * x0[j];
}
}
x[i] = (b[i] - sum) / a[i][i];
}// 更新迭代向量

c_M++; //迭代次数加1
for(i=0;i<N;i++)
{
if(fabs(x[i]-x0[i])>current_e)
current_e = fabs(x[i]-x0[i]);
} //计算当前误差
for(i=0;i<N;i++)
x0[i] = x[i]; //更新初始向量
}while(current_e>e&&c_M<M); //判断是否仍未达到精度要求且未达到最大迭代次数

for(i=0;i<N;i++)
cout << x[i] << endl; //输出结果
cout << c_M << endl; //输出迭代次数
}
二：高斯-塞德尔(Gauss-Seidel)迭代算法
程序
#define N 3 // 线性方程组的阶数
#include <iostream.h>
#include <math.h>
void main()
{
double a[N][N]={5,2,1,
2,8,-3,
1,-3,-6}, //系数矩阵
b[N]={8,21,1}; //右端常数向量
double x0[N]={1,1,1},x[N]; // 迭代初始向量和迭代向量
double e=1e-5; // 精度要求
int M=5000; //最大迭代次数
int i,j,c_M=0;
double sum,current_e;

do{
current_e=0;
for(i=0;i<N;i++)
{
sum = 0;
for(j=0;j<N;j++)
{
if(j<i)
{
sum = sum + a[i][j] * x[j];
}
if(j>i)
{
sum = sum + a[i][j] * x0[j];
}
}
x[i] = (b[i] - sum) / a[i][i];
}// 更新迭代向量

c_M++; //迭代次数加1
for(i=0;i<N;i++)
{
if(fabs(x[i]-x0[i])>current_e)
current_e = fabs(x[i]-x0[i]);
} //计算当前误差
for(i=0;i<N;i++)
x0[i] = x[i]; //更新初始向量
}while(current_e>e&&c_M<M); //判断是否仍未达到精度要求且未达到最大迭代次数

for(i=0;i<N;i++)
cout << x[i] << endl; //输出结果
cout << c_M << endl; //输出迭代次数
}
三：SOR迭代算法
程序
#define N 4 // 线性方程组的阶数
#include <iostream.h>
#include <math.h>
void main()
{
double a[N][N]={-4,1,1,1,
1,-4,1,1,
1,1,-4,1,
1,1,1,-4}, //系数矩阵
b[N]={1,1,1,1}; //右端常数向量
double x0[N]={1,1,1,1},x[N]; // 迭代初始向量和迭代向量
double e=1e-5,w=0.5; // 精度要求
int M=5000; //最大迭代次数
int i,j,c_M=0;
double sum,current_e;

do{
current_e=0;
for(i=0;i<N;i++)
{
sum = 0;
for(j=0;j<N;j++)
{
if(j<i)
{
sum = sum + a[i][j] * x[j];
}
if(j>i)
{
sum = sum + a[i][j] * x0[j];
}

}
x[i] = (1-w)*x0[i] + w*(b[i] - sum) / a[i][i];
}// 更新迭代向量

c_M++; //迭代次数加1
for(i=0;i<N;i++)
{
if(fabs(x[i]-x0[i])>current_e)
current_e = fabs(x[i]-x0[i]);
} //计算当前误差
for(i=0;i<N;i++)
x0[i] = x[i]; //更新初始向量
}while(current_e>e&&c_M<M); //判断是否仍未达到精度要求且未达到最大迭代次数

for(i=0;i<N;i++)
cout << x[i] << endl; //输出结果
cout << c_M << endl; //输出迭代次数
}