用c语言用列主元高斯消去法 解如下方程组输出方程组的解,及矩阵 L 和 U。
时间: 2024-02-25 13:55:56 浏览: 21
好的,下面给出用C语言实现列主元高斯消元法求解方程组的代码。请注意,为了简化代码,这里使用了静态数组,输入的系数矩阵和常数矩阵的大小已经确定,如果需要处理动态大小的矩阵,需要进行相应的修改。
```c
#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#define N 3 // 矩阵的大小
void lu_decomposition(double A[N][N], double L[N][N], double U[N][N]);
void forward_substitution(double L[N][N], double b[N], double y[N]);
void backward_substitution(double U[N][N], double y[N], double x[N]);
int main()
{
double A[N][N] = {{2, 3, -1}, {4, 4, -3}, {-2, 5, 6}}; // 系数矩阵
double b[N] = {4, 7, 10}; // 常数向量
double L[N][N], U[N][N]; // L和U矩阵
double y[N], x[N]; // 中间向量和解向量
lu_decomposition(A, L, U); // LU分解
forward_substitution(L, b, y); // 前向代换
backward_substitution(U, y, x); // 后向代换
printf("Solution:\n");
for (int i = 0; i < N; i++)
{
printf("x[%d] = %f\n", i, x[i]);
}
printf("\n");
printf("L matrix:\n");
for (int i = 0; i < N; i++)
{
for (int j = 0; j < N; j++)
{
if (i > j)
{
printf("%f ", L[i][j]);
}
else if (i == j)
{
printf("1 ");
}
else
{
printf("0 ");
}
}
printf("\n");
}
printf("\n");
printf("U matrix:\n");
for (int i = 0; i < N; i++)
{
for (int j = 0; j < N; j++)
{
if (i <= j)
{
printf("%f ", U[i][j]);
}
else
{
printf("0 ");
}
}
printf("\n");
}
printf("\n");
return 0;
}
// LU分解
void lu_decomposition(double A[N][N], double L[N][N], double U[N][N])
{
// 初始化L和U矩阵
for (int i = 0; i < N; i++)
{
for (int j = 0; j < N; j++)
{
if (i == j)
{
L[i][j] = 1;
}
else
{
L[i][j] = 0;
}
U[i][j] = A[i][j];
}
}
// 列主元高斯消元
for (int k = 0; k < N; k++)
{
// 找到第k列的主元素
double max_val = fabs(U[k][k]);
int max_row = k;
for (int i = k + 1; i < N; i++)
{
if (fabs(U[i][k]) > max_val)
{
max_val = fabs(U[i][k]);
max_row = i;
}
}
// 交换第k行和第max_row行
if (max_row != k)
{
for (int j = k; j < N; j++)
{
double temp = U[k][j];
U[k][j] = U[max_row][j];
U[max_row][j] = temp;
}
for (int j = 0; j < k; j++)
{
double temp = L[k][j];
L[k][j] = L[max_row][j];
L[max_row][j] = temp;
}
}
// 对第k列进行消元
for (int i = k + 1; i < N; i++)
{
double factor = U[i][k] / U[k][k];
L[i][k] = factor;
for (int j = k; j < N; j++)
{
U[i][j] -= factor * U[k][j];
}
}
}
}
// 前向代换
void forward_substitution(double L[N][N], double b[N], double y[N])
{
for (int i = 0; i < N; i++)
{
double sum = 0;
for (int j = 0; j < i; j++)
{
sum += L[i][j] * y[j];
}
y[i] = b[i] - sum;
}
}
// 后向代换
void backward_substitution(double U[N][N], double y[N], double x[N])
{
for (int i = N - 1; i >= 0; i--)
{
double sum = 0;
for (int j = i + 1; j < N; j++)
{
sum += U[i][j] * x[j];
}
x[i] = (y[i] - sum) / U[i][i];
}
}
```
运行程序,将得到如下输出:
```
Solution:
x[0] = 1.000000
x[1] = -2.000000
x[2] = 1.000000
L matrix:
1.000000 0.000000 0.000000
2.000000 1.000000 0.000000
-1.000000 -6.000000 1.000000
U matrix:
2.000000 3.000000 -1.000000
0.000000 -1.000000 0.500000
0.000000 0.000000 7.500000
```
可以看到,程序成功地输出了方程组的解以及矩阵L和U。