# [计算方法]列主元消去法 LU分解法求线性方程解

/*计算方法实验二
*WYC
*计科1604
*/
#include <cstdio>
#include <cmath>
#include <iostream>
using namespace std;

const int MAXSIZE = 20;
void columnPivoting(double a[MAXSIZE][MAXSIZE],int n, double b[]);   //声明函数columnPivoting

int main(){
double matrix[MAXSIZE][MAXSIZE];//系数
double b[MAXSIZE]; //b
int sizeM;
cin>> sizeM;
for(int i=0;i<sizeM;i++){
for(int j=0;j<sizeM;j++)
cin>> matrix[i][j];
cin>>b[i];
}
columnPivoting(matrix,sizeM,b);
return 0;
}

//列主元消去法
void columnPivoting(double a[MAXSIZE][MAXSIZE], int n, double b[]){
int i, j, k;
int col = 0, row = 0;
for (k = 0; k < n - 1; k++){
double tmp = 0;
//找出消元列中最大的那个元素所在的位置
for (i = k; i < n ; i++)
if (fabs(a[i][k]) > tmp){
tmp = fabs(a[i][k]);
//cout << "MAX " << tmp << endl;//
row = i;
col = k;
}

//如果该对角线元素是0，同样不能用高斯消元法来求解
if (a[row][row] == 0){
cout << "can't solve" << endl;
return;
}
//将找出的行进行交换
if (k != row){
for (i = 0; i < n; i++){
swap(a[row][i], a[k][i]);
swap(b[k], b[row]);
}
}
//消元过程
double c[n];
for (j = k + 1; j < n; j++){
c[j] = a[j][k] / a[k][k];
//cout << c[j] << endl;
}
for (i = k + 1; i < n; i++){
for (j = 1; j < n; j++){
a[i][j] = a[i][j] - c[i] * a[k][j];
}
b[i] = b[i] - c[i] * b[k];
}

}

double x[n];
x[n - 1] = b[n - 1] / a[n - 1][n - 1];
for (i = n - 2; i >= 0; i--){
double sum = 0;
for(j = i + 1; j < n; j++)
sum += a[i][j] * x[j];
x[i] = (b[i] - sum)/a[i][i];
}
for(int i=0;i<n;i++)
printf("x[%d] = %.4f\n",i+1,x[i]);  //输出结果
printf("\n");
}

/*
4
1 1 0 3 4
2 1 -1 1 1
3 -1 -1 3 -3
-1 2 3 -1 4
3
1 1 -1 1
1 2 -2 0
-2 1 1 1
4
3 1 0 0 1
2 3 1 0 0
0 2 3 1 1
0 0 1 3 0

*/
/*计算方法实验二 LU Decomposition
*WYC
*计科1604
*/
#include<iostream>
using namespace std;

const int MAXSIZE = 20;
void LUDecomposition(double a[MAXSIZE][MAXSIZE], double b[MAXSIZE],int n);

int main()
{
double a[MAXSIZE][MAXSIZE];
double b[MAXSIZE];
int n;
cin>> n;
for(int i=0;i<n;i++){
for(int j=0;j<n;j++)
cin>> a[i][j];
cin>>b[i];
}
LUDecomposition(a, b, n);
return 0;
}

//矩阵的ALU分解
void LUDecomposition(double a[MAXSIZE][MAXSIZE], double b[MAXSIZE],int n)
{
double l[MAXSIZE][MAXSIZE] = { 0 };
double u[MAXSIZE][MAXSIZE] = { 0 };
int i, r, k;
//进行U的第一行的赋值
for (i = 0; i<n; i++)
{
u[0][i] = a[0][i];
}

//进行L的第一列的赋值
for (i = 1; i<n; i++)
{
l[i][0] = a[i][0] / u[0][0];
}

//计算U的剩下的行数和L的剩下的列数
for (r = 1; r<n; r++)
{
for (i = r; i <n; i++)
{
double sum1 = 0;
for (k = 0; k < r; k++)
{
sum1 += l[r][k] * u[k][i];
//cout << "" << r << "" << sum1 << endl;
}
u[r][i] = a[r][i] - sum1;
}

if(r!=n)
for(i=r+1;i<n;i++)
{
double sum2 = 0;
for (k = 0; k<r; k++)
{
sum2 += l[i][k] * u[k][r];
}
l[i][r] = (a[i][r] - sum2) / u[r][r];
}

}

double y[MAXSIZE] = { 0 };
y[0] = b[0];
for (i = 1; i<n; i++)
{
double sum3 = 0;
for (k = 0; k<i; k++)
sum3 += l[i][k] * y[k];
y[i] = b[i] - sum3;
}

double x[MAXSIZE] = { 0 };
x[n - 1] = y[n - 1] / u[n - 1][n - 1];
for (i = n - 2; i >= 0; i--)
{
double sum4 = 0;
for (k = i + 1; k<n; k++)
sum4 += u[i][k] * x[k];
x[i] = (y[i] - sum4) / u[i][i];
}
for (i = 0; i<n; i++)
cout << "x[" << i + 1 << "]=" << x[i] << endl;

}



• 广告
• 抄袭
• 版权
• 政治
• 色情
• 无意义
• 其他

120