最近在上计算方法这门课,要求是用MATLAB做练习题,但是我觉得C语言也很棒棒啊~
题目:
高斯消元法是线性方程组的直接解法,可能会造成很大的失真,尤其是高斯顺序消元法,对方法进行改进,使每次都选取绝对值最大的元素为主元,使其为乘数的分母,控制舍入误差的扩大,失真较小。代码都是上上个星期写的,暂时就不注释了……
使用VS2017,代码如下:
//使用全主元高斯消元法求解线性方程组
#include "stdafx.h"
#include<stdlib.h>
#include "math.h"
double **A, *b, *x;
unsigned int *x_number;
unsigned int RANK = 4;
unsigned int makematrix()
{
unsigned int r, c;
printf("请输入矩阵行列数,用空格隔开:");
scanf_s("%d %d", &r, &c);
A = (double**)malloc(sizeof(double*)*r);//创建一个指针数组,把指针数组的地址赋值给a ,*r是乘以r的意思
for (int i = 0; i < r; i++)
A[i] = (double*)malloc(sizeof(double)*c);//给第二维分配空间
for (int i = 0; i < r; i++) {
for (int j = 0; j < c; j++)
A[i][j] = 0.0;
}
b = (double*)malloc(sizeof(double)*r);
for (int i = 0; i < r; i++)
{
b[i] = 0.0;
}
x = (double*)malloc(sizeof(double)*c);
for (int i = 0; i < c; i++)
{
x[i] = 0.0;
}
x_number = (unsigned int*)malloc(sizeof(unsigned int)*c);
for (int i = 0; i < c; i++)
{
x_number[i] = i + 1;
}
return r;
}
void getmatrix(void)//输入矩阵并呈现
{
printf("请按行从左到右依次输入系数矩阵A,不同元素用空格隔开\n");
for (int i = 0; i < RANK; i++)
{
for (int j = 0; j<RANK; j++)
{
scanf_s("%lf", &A[i][j]);
}
}
printf("系数矩阵如下\n");
for (int i = 0; i < RANK; i++)
{
for (int j = 0; j<RANK; j++)
{
printf("%g\t", A[i][j]);
}
printf("\n");
}
printf("请按从上到下依次输入常数列b,不同元素用空格隔开\n");
for (int i = 0; i<RANK; i++)
{
scanf_s("%lf", &b[i]);
}
printf("常数列如下\n");
for (int i = 0; i<RANK; i++)
{
printf("%g\t", b[i]);
}printf("\n");
}
void exchange_matrix(unsigned int n,double **AA,double *bb, unsigned int *xx) //n为从第几行(列)开始,AA为系数矩阵,bb为常数列
{
double get_max = 0.0;
unsigned int get_max_i, get_max_j;
get_max_i = n - 1;
get_max_j = n - 1;
for (int i = n-1; i < RANK; i++)
{
for (int j = n-1; j<RANK; j++)
{
if (fabs(AA[i][j]) > fabs(get_max))
{
get_max = AA[i][j];
get_max_i = i;
get_max_j= j;
}
}
}
if (get_max_i != n - 1)//交换行
{
double *temp, temp2;
temp = AA[get_max_i];
AA[get_max_i] = AA[n - 1];
AA[n - 1] = temp;
temp2 = bb[get_max_i];
bb[get_max_i] = bb[n - 1];
bb[n - 1] = temp2;
}
if (get_max_j != n - 1)//交换列
{
double temp;
unsigned int temp2;
for (int i = 0; i < RANK; i++)//系数
{
temp = AA[i][get_max_j];
AA[i][get_max_j] = AA[i][n - 1];
AA[i][n - 1] = temp;
}
temp2 = xx[get_max_j];//x下标
xx[get_max_j] = xx[n - 1];
xx[n - 1] = temp2;
}
printf("第%d此换序后矩阵如下\n",n);
for (int i = 0; i < RANK; i++)
{
for (int j = 0; j<RANK; j++)
{
printf("%g\t", AA[i][j]);
}
printf(" %g", bb[i]);
printf("\n");
}
}
void Gauss_calculation(void)//Gauss全主元消去法解线性方程组
{
double get_A = 0.0;
printf("利用以上A与b组成的增广阵进行全主元高斯消去法计算方程组\n");
for (int i = 1; i < RANK; i++)
{
exchange_matrix(i, A, b,x_number);//换序
for (int j = i; j<RANK; j++)
{
get_A = A[j][i - 1] / A[i - 1][i - 1];
b[j] = b[j] - get_A * b[i - 1];
for (int k = i - 1; k < RANK; k++)
{
A[j][k] = A[j][k] - get_A * A[i - 1][k];
}
}
}
printf("消元后的上三角系数增广矩阵如下\n");
for (int i = 0; i < RANK; i++)
{
for (int j = 0; j<RANK; j++)
{
printf("%g\t", A[i][j]);
}
printf(" %g", b[i]);
printf("\n");
}
printf("利用回代法求解上三角方程组,解得:\n");
for (int i = 0; i < RANK; i++)
{
double get_x = 0.0;
for (int j = 0; j < RANK; j++)
{
get_x = get_x + A[RANK - 1 - i][j] * x[j];//把左边全部加起来了,下面需要多减去一次Xn*Ann
}
x[RANK - 1 - i] = (b[RANK - 1 - i] - get_x + A[RANK - 1 - i][RANK - 1 - i] * x[RANK - 1 - i]) / A[RANK - 1 - i][RANK - 1 - i];
}
for (int i = 0; i < RANK; i++)
{
printf("x%d = %g\n", x_number[i], x[i]);
}
printf("计算完成,按回车退出程序或按1重新输入矩阵\n");
}
int main()
{
_again:
RANK = makematrix();
getmatrix();
Gauss_calculation();
getchar();
if ('1' == getchar())
goto _again;
return 0;
}
按设计的提示老老实实 输入题目的系数矩阵和常数向量后,得到运行结果:
一般来说直接求法是很快的,但是一般来说只用在不大的方程组上,因为失真会不断被放大,那就很扎心了。