<手搓有限元 干翻Ansys> 【1】C++实现矩阵Matrix类 实现基本运算
2025年2月5日 19:57浏览：1328 评论：4 收藏：7
本系列文章致力于实现“手搓有限元，干翻Ansys的目标”，基本框架为前端显示使用QT实现交互，后端计算采用Visual Studio C++。
Matrix类
矩阵基本类，用于有限元矩阵计算。
1、public function
1.1、构造函数与析构函数
构造函数用来初始化矩阵，析构函数用来释放内存。
Matrix.h声明文件：
//******************************构造函数与析构函数********************************// /* 函数名称: 无参构造函数 */ Matrix(); /* 函数名称: 矩阵有参构造函数，初始化为row行、col列的0矩阵 row: 矩阵行数 col: 矩阵列数 */ Matrix(int row, int col); /* 函数名称: 矩阵有参构造函数，初始化为row行、col列、数值为mat的矩阵 row: 矩阵行数 col: 矩阵列数 *mat: 矩阵数值一维数组 */ Matrix(int row, int col, double* mat); /* 函数名称: 深拷贝构造函数 mat: 需要复制的矩阵 */ Matrix(const Matrix& mat); /* 函数名称: 析构函数 */ ~Matrix();
Matrix.cpp函数实现文件：
Matrix::Matrix()
{
}
//初始化矩阵 默认值为0
Matrix::Matrix(int row, int col)
{
this->m_Row = row;
this->m_Col = col;
//开辟内存
this->m_Matrix = new double* [row];
for (int i = 0; i < row; i++)
{
this->m_Matrix[i] = new double[col] {0.0};
}
}
//初始化矩阵 设定数值
Matrix::Matrix(int row, int col, double *mat)
{
this->m_Row = row;
this->m_Col = col;
//开辟内存
this->m_Matrix = new double* [row];
for (int i = 0; i < row; i++)
{
this->m_Matrix[i] = new double[col] {0.0};
}
//矩阵赋值
for(int i = 0; i<row; i++)
{
for (int j = 0; j < col; j++)
{
this->m_Matrix[i][j] = mat[i * col + j];
}
}
}
//深拷贝
Matrix::Matrix(const Matrix& mat)
{
//行列传递
this->m_Row = mat.m_Row;
this->m_Col = mat.m_Col;
//矩阵深拷贝
this->m_Matrix = new double* [this->m_Row];
for (int i = 0; i < this->m_Row; i++)
{
this->m_Matrix[i] = new double[this->m_Col];
memcpy(this->m_Matrix[i], mat.m_Matrix[i], sizeof(double) * this->m_Col);
}
}
//析构函数
Matrix::~Matrix()
{
//释放矩阵每一行
for (int i = 0; i < this->m_Row; i++)
{
if (this->m_Matrix[i] != NULL)
{
delete[]this->m_Matrix[i];
this->m_Matrix[i] = NULL;
}
}
//释放矩阵顶点
if (this->m_Matrix != NULL)
{
delete[]this->m_Matrix;
this->m_Matrix = NULL;
}
}
1.2、获取矩阵数值
可以获取矩阵指定位置数值、打印矩阵。
Matrix.h声明文件：
//*******************获取矩阵*****************// /* 函数名称: 获取矩阵的第row行、第col列元素数值 row: 矩阵行数 col: 矩阵列数 */ double GetMatrixEle(int row, int col);
Matrix.cpp函数实现文件：
//获取矩阵某个元素 某行某列
double Matrix::GetMatrixEle(int row, int col)
{
if (row >= this->m_Row)
{
std::cout << "Error: <GetMatrixEle> Input row >= m_Row" << std::endl;
return 0.0;
}
else if (col >= this->m_Col)
{
std::cout << "Error: <GetMatrixEle> Input col >= m_Col" << std::endl;
return 0.0;
}
else
{
return this->m_Matrix[row][col];
}
}
//矩阵输出
void Matrix::PrintMat()
{
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < this->m_Col; j++)
{
std::cout.setf(std::ios::scientific); //科学计数法表示
std::cout << this->m_Matrix[i][j] << "\t";
}
std::cout << std::endl;
}
std::cout << std::endl;
}
1.3、设置矩阵
可进行设置矩阵指定位置数值，以及深拷贝矩阵。
Matrix.h声明文件：
//*******************设置矩阵*****************// /* 函数名称: 设置矩阵第row行、第col列数值 row: 矩阵行数 col: 矩阵列数 value: 设置的矩阵数值 */ void SetMatrixEle(int row, int col, double value); /* 函数名称: 深拷贝矩阵 mat: 需要复制的矩阵 */ Matrix CopyMat(const Matrix mat);
Matrix.cpp函数实现文件：
//*******************设置矩阵*****************//
void Matrix::SetMatrixEle(int row, int col, double value)
{
if (row >= this->m_Row)
{
std::cout << "Error: <SetMatrixEle> Input row >= m_Row" << std::endl;
return;
}
else if (col >= this->m_Col)
{
std::cout << "Error: <SetMatrixEle> Input col >= m_Col" << std::endl;
return;
}
else
{
this->m_Matrix[row][col] = value;
return;
}
}
//深拷贝矩阵
Matrix Matrix::CopyMat(const Matrix mat)
{
//行列传递
this->m_Row = mat.m_Row;
this->m_Col = mat.m_Col;
//矩阵深拷贝
this->m_Matrix = new double* [this->m_Row];
for (int i = 0; i < this->m_Row; i++)
{
this->m_Matrix[i] = new double[this->m_Col];
memcpy(this->m_Matrix[i], mat.m_Matrix[i], sizeof(double) * this->m_Col);
}
return *this;
}
1.4、矩阵转置、单位化
可进行矩阵转置，单位化，注意返回值类型为自身的引用，可实现链式编程。
Matrix.h声明文件：
//*****************矩阵基本操作***************// /* 函数名称: 矩阵转置，返回的是自身引用，可链式调用 */ Matrix& Transpose(); /* 函数名称: 等维度的单位矩阵，前提是方阵 */ Matrix& Uint();
Matrix.cpp函数实现文件：
//矩阵转置
Matrix& Matrix::Transpose()
{
Matrix* resMat = new Matrix(this->m_Col, this->m_Row);
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < this->m_Col; j++)
{
resMat->m_Matrix[j][i] = this->m_Matrix[i][j];
}
}
return *resMat;
}
//求等长度单位矩阵
Matrix& Matrix::Uint()
{
//矩阵是否为方阵
if (this->m_Col != this->m_Row)
{
std::cout << "Error: <Uint> Row != Col" << std::endl;
Matrix* resMat = new Matrix(this->m_Row, this->m_Row);
return *resMat;
}
else
{
//单位矩阵初始化
Matrix* resMat = new Matrix(this->m_Row, this->m_Col);
//单位矩阵生成
for (int i = 0; i < this->m_Row; i++)
{
resMat->m_Matrix[i][i] = 1.0;
}
return *resMat;
}
}
1.5、矩阵的删除与替换
可进行矩阵指定行、列的删除与替换，注意返回值类型为自身的引用，可实现链式编程。
Matrix.h声明文件：
//****************矩阵保留与剔除**************// /* 函数名称: 剔除矩阵中以index为行标和列标的行和列，num代表index的大小 *index: 矩阵中的行号与列号一维数组 num: index动态数组长度 */ Matrix& DeleteMat(int *index, int num); /* 函数名称: 剔除矩阵中以index为行标和列标的行和列，num代表index的大小 *index: 矩阵中的行号与列号一维动态数组 num: index动态数组长度 */ Matrix& DeleteMat(std::vector<int> index, int num); /* 函数名称: 剔除矩阵中以index为行标的行，num代表index的大小 *index: 矩阵中的行号一维数组 num: index动态数组长度 */ Matrix& DeleteRow(int* index, int num); /* 函数名称: 剔除矩阵中以index为行标的行，num代表index的大小 *index: 矩阵中的行号一维动态数组 num: index动态数组长度 */ Matrix& DeleteRow(std::vector<int> index, int num); /* 函数名称: 剔除矩阵中以index为列标的列，num代表index的大小 *index: 矩阵中的列号一维数组 num: index动态数组长度 */ Matrix& DeleteCol(int* index, int num); /* 函数名称: 剔除矩阵中以index为列标的列，num代表index的大小 *index: 矩阵中的列号一维动态数组 num: index动态数组长度 */ Matrix& DeleteCol(std::vector<int> index, int num); //******************矩阵的替换****************// /* 函数名称: 替换矩阵中行标和列标为 index中的行与列，num代表index的大小, mat是需要替换 *index: 矩阵中的行标和列标的一维数组 num: index动态数组长度 mat: 需要替换的矩阵 */ Matrix& ReplaceMat(int* index, int num, Matrix& mat); /* 函数名称: 替换矩阵中行标和列标为 index中的行与列，num代表index的大小, mat是需要替换 *index: 矩阵中的行标和列标的一维动态数组 num: index动态数组长度 mat: 需要替换的矩阵 */ Matrix& ReplaceMat(std::vector<int> index, int num, Matrix& mat); /* 函数名称: 替换矩阵中行标为 index中的行，num代表index的大小, mat是需要替换的矩阵 *index: 矩阵中的行标的一维数组 num: index动态数组长度 mat: 需要替换的矩阵 */ Matrix& ReplaceRow(int* index, int num, Matrix& mat); /* 函数名称: 替换矩阵中行标为 index中的行，num代表index的大小, mat是需要替换的矩阵 *index: 矩阵中的行标的一动态维数组 num: index动态数组长度 mat: 需要替换的矩阵 */ Matrix& ReplaceRow(std::vector<int> index, int num, Matrix& mat); /* 函数名称: 替换矩阵中列标为 index中的列，num代表index的大小, mat是需要替换的矩阵 *index: 矩阵中的列标的一维数组 num: index动态数组长度 mat: 需要替换的矩阵 */ Matrix& ReplaceCol(int* index, int num, Matrix& mat); /* 函数名称: 替换矩阵中列标为 index中的列，num代表index的大小, mat是需要替换的矩阵 *index: 矩阵中的列标的一维动态数组 num: index动态数组长度 mat: 需要替换的矩阵 */ Matrix& ReplaceCol(std::vector<int> index, int num, Matrix& mat);
Matrix.cpp函数实现文件：
//****************矩阵保留与剔除**************//
//剔除矩阵的 index中的行与列，num代表index的大小
Matrix& Matrix::DeleteMat(int* index, int num)
{
//结果矩阵
Matrix* resMat = new Matrix(this->m_Row - num, this->m_Col - num);
int recIndex[MAX_COUNT];
int currIndex = 0;
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
else if (index[i] >= this->m_Col)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Col" << std::endl;
return *this;
}
}
//筛选出剔除后行数
for (int iRow = 0; iRow < this->m_Row; iRow++)
{
for (int iNum = 0; iNum < num; iNum++)
{
if (iRow == index[iNum])
{
break;
}
if (iNum == num-1)
{
recIndex[currIndex++] = iRow;
}
}
}
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[iRow][iCol] = this->m_Matrix[recIndex[iRow]][recIndex[iCol]];
}
}
return *resMat;
}
Matrix& Matrix::DeleteMat(std::vector<int> index, int num)
{
//结果矩阵
Matrix* resMat = new Matrix(this->m_Row - num, this->m_Col - num);
int recIndex[MAX_COUNT];
int currIndex = 0;
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
else if (index[i] >= this->m_Col)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Col" << std::endl;
return *this;
}
}
//筛选出剔除后行数
for (int iRow = 0; iRow < this->m_Row; iRow++)
{
for (int iNum = 0; iNum < num; iNum++)
{
if (iRow == index[iNum])
{
break;
}
if (iNum == num - 1)
{
recIndex[currIndex++] = iRow;
}
}
}
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[iRow][iCol] = this->m_Matrix[recIndex[iRow]][recIndex[iCol]];
}
}
return *resMat;
}
//剔除矩阵的 index中的行，num代表index的大小
Matrix& Matrix::DeleteRow(int* index, int num)
{
//结果矩阵
Matrix* resMat = new Matrix(this->m_Row - num, this->m_Col);
int recIndex[MAX_COUNT];
int currIndex = 0;
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
}
//筛选出剔除后行数
for (int iRow = 0; iRow < this->m_Row; iRow++)
{
for (int iNum = 0; iNum < num; iNum++)
{
if (iRow == index[iNum])
{
break;
}
if (iNum == num - 1)
{
recIndex[currIndex++] = iRow;
}
}
}
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[iRow][iCol] = this->m_Matrix[recIndex[iRow]][iCol];
}
}
return *resMat;
}
Matrix& Matrix::DeleteRow(std::vector<int> index, int num)
{
//结果矩阵
Matrix* resMat = new Matrix(this->m_Row - num, this->m_Col);
int recIndex[MAX_COUNT];
int currIndex = 0;
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
}
//筛选出剔除后行数
for (int iRow = 0; iRow < this->m_Row; iRow++)
{
for (int iNum = 0; iNum < num; iNum++)
{
if (iRow == index[iNum])
{
break;
}
if (iNum == num - 1)
{
recIndex[currIndex++] = iRow;
}
}
}
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[iRow][iCol] = this->m_Matrix[recIndex[iRow]][iCol];
}
}
return *resMat;
}
Matrix& Matrix::DeleteCol(int* index, int num)
{
//结果矩阵
Matrix* resMat = new Matrix(this->m_Row, this->m_Col - num);
int recIndex[MAX_COUNT];
int currIndex = 0;
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
}
//筛选出剔除后行数
for (int iRow = 0; iRow < this->m_Row; iRow++)
{
for (int iNum = 0; iNum < num; iNum++)
{
if (iRow == index[iNum])
{
break;
}
if (iNum == num - 1)
{
recIndex[currIndex++] = iRow;
}
}
}
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[iRow][iCol] = this->m_Matrix[iRow][recIndex[iCol]];
}
}
return *resMat;
}
Matrix& Matrix::DeleteCol(std::vector<int> index, int num)
{
//结果矩阵
Matrix* resMat = new Matrix(this->m_Row, this->m_Col - num);
int recIndex[MAX_COUNT];
int currIndex = 0;
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
}
//筛选出剔除后行数
for (int iRow = 0; iRow < this->m_Row; iRow++)
{
for (int iNum = 0; iNum < num; iNum++)
{
if (iRow == index[iNum])
{
break;
}
if (iNum == num - 1)
{
recIndex[currIndex++] = iRow;
}
}
}
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[iRow][iCol] = this->m_Matrix[iRow][recIndex[iCol]];
}
}
return *resMat;
}
//******************矩阵的替换****************//
//替换矩阵中的行和列 index中的行与列，num代表index的大小
Matrix& Matrix::ReplaceMat(int* index, int num, Matrix& mat)
{
//错误判定 方阵
if (this->m_Row != this->m_Col)
{
std::cout << "Error: <ReplaceMat> this m_Col != m_Row" << std::endl;
return *this;
}
//检验插入矩阵为方阵
if (mat.m_Row != mat.m_Col)
{
std::cout << "Error: <ReplaceMat> mat m_Col != m_Row" << std::endl;
return *this;
}
//检验插入矩阵大小与num保持一致
if (mat.m_Col != num)
{
std::cout << "Error: <ReplaceMat> num != mat.m_Col" << std::endl;
return *this;
}
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <ReplaceMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
else if (index[i] >= this->m_Col)
{
std::cout << "Error: <ReplaceMat> Input index[" << i << "] = " << index[i] << " >= m_Col" << std::endl;
return *this;
}
}
//结果矩阵
Matrix* resMat = new Matrix(*this);
//加入元素
for (int iRow = 0; iRow < num; iRow++)
{
for (int iCol = 0; iCol < num; iCol++)
{
resMat->m_Matrix[index[iRow]][index[iCol]] = mat.m_Matrix[iRow][iCol];
}
}
return *resMat;
}
Matrix& Matrix::ReplaceMat(std::vector<int> index, int num, Matrix& mat)
{
//错误判定 方阵
if (this->m_Row != this->m_Col)
{
std::cout << "Error: <ReplaceMat> this m_Col != m_Row" << std::endl;
return *this;
}
//检验插入矩阵为方阵
if (mat.m_Row != mat.m_Col)
{
std::cout << "Error: <ReplaceMat> mat m_Col != m_Row" << std::endl;
return *this;
}
//检验插入矩阵大小与num保持一致
if (mat.m_Col != num)
{
std::cout << "Error: <ReplaceMat> num != mat.m_Col" << std::endl;
return *this;
}
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <ReplaceMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
else if (index[i] >= this->m_Col)
{
std::cout << "Error: <ReplaceMat> Input index[" << i << "] = " << index[i] << " >= m_Col" << std::endl;
return *this;
}
}
//结果矩阵
Matrix* resMat = new Matrix(*this);
//加入元素
for (int iRow = 0; iRow < num; iRow++)
{
for (int iCol = 0; iCol < num; iCol++)
{
resMat->m_Matrix[index[iRow]][index[iCol]] = mat.m_Matrix[iRow][iCol];
}
}
return *resMat;
}
//替换矩阵中的行 index中的行，num代表index的大小, mat是需要替换的矩阵
Matrix& Matrix::ReplaceRow(int* index, int num, Matrix& mat)
{
//检验插入矩阵大小与num保持一致
if (mat.m_Row != num)
{
std::cout << "Error: <ReplaceRow> num != mat.m_Row" << std::endl;
return *this;
}
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <ReplaceRow> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
}
//当前矩阵列数应与mat列数一致
if (this->m_Col != mat.m_Col)
{
std::cout << "Error: <ReplaceRow> this->m_Col != mat.m_Col" << std::endl;
return *this;
}
//结果矩阵
Matrix* resMat = new Matrix(*this);
//加入元素
for (int iRow = 0; iRow < num; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[index[iRow]][iCol] = mat.m_Matrix[iRow][iCol];
}
}
return *resMat;
}
Matrix& Matrix::ReplaceRow(std::vector<int> index, int num, Matrix& mat)
{
//检验插入矩阵大小与num保持一致
if (mat.m_Row != num)
{
std::cout << "Error: <ReplaceRow> num != mat.m_Row" << std::endl;
return *this;
}
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <ReplaceRow> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
}
//当前矩阵列数应与mat列数一致
if (this->m_Col != mat.m_Col)
{
std::cout << "Error: <ReplaceRow> this->m_Col != mat.m_Col" << std::endl;
return *this;
}
//结果矩阵
Matrix* resMat = new Matrix(*this);
//加入元素
for (int iRow = 0; iRow < num; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[index[iRow]][iCol] = mat.m_Matrix[iRow][iCol];
}
}
return *resMat;
}
//替换矩阵中的列 index中的列，num代表index的大小, mat是需要替换的矩阵
Matrix& Matrix::ReplaceCol(int* index, int num, Matrix& mat)
{
//检验插入矩阵大小与num保持一致
if (mat.m_Col != num)
{
std::cout << "Error: <ReplaceCol> mat.m_Col != num" << std::endl;
return *this;
}
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Col)
{
std::cout << "Error: <ReplaceCol> Input index[" << i << "] = " << index[i] << " >= m_Col" << std::endl;
return *this;
}
}
//当前矩阵行数应与mat行数一致
if (this->m_Row != mat.m_Row)
{
std::cout << "Error: <ReplaceCol> this->m_Row != mat.m_Row" << std::endl;
return *this;
}
//结果矩阵
Matrix* resMat = new Matrix(*this);
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < num; iCol++)
{
resMat->m_Matrix[iRow][index[iCol]] = mat.m_Matrix[iRow][iCol];
}
}
return *resMat;
}
Matrix& Matrix::ReplaceCol(std::vector<int> index, int num, Matrix& mat)
{
//检验插入矩阵大小与num保持一致
if (mat.m_Col != num)
{
std::cout << "Error: <ReplaceCol> mat.m_Col != num" << std::endl;
return *this;
}
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Col)
{
std::cout << "Error: <ReplaceCol> Input index[" << i << "] = " << index[i] << " >= m_Col" << std::endl;
return *this;
}
}
//当前矩阵行数应与mat行数一致
if (this->m_Row != mat.m_Row)
{
std::cout << "Error: <ReplaceCol> this->m_Row != mat.m_Row" << std::endl;
return *this;
}
//结果矩阵
Matrix* resMat = new Matrix(*this);
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < num; iCol++)
{
resMat->m_Matrix[iRow][index[iCol]] = mat.m_Matrix[iRow][iCol];
}
}
return *resMat;
}
1.6、矩阵初等变换
可实现矩阵的初等变化，注意返回值类型为自身的引用，可实现链式编程。
Matrix.h声明文件：
//*****************矩阵初等变化***************// /* 函数名称: 交换矩阵中行标为row0与row1的元素 row0: 矩阵行标0 row1: 矩阵行标1 */ Matrix& SwapRow(int row0, int row1); /* 函数名称: 交换矩阵中列标为col0与col1的元素 col0: 矩阵列标0 col1: 矩阵列标1 */ Matrix& SwapCol(int col0, int col1); /* 函数名称: 矩阵行加法 rowLocal = rowLocal + rate *rowAdd rowLocal: 矩阵行标，被加数 rowAdd: 矩阵行标，加数 rate: 加数前倍数 */ Matrix& AddRow(int rowLocal, int rowAdd, double rate = 1.0); //矩阵加法 某列 + 倍数*某列 /* 函数名称: 矩阵列加法 colLocal = colLocal + rate * colAdd colLocal: 矩阵列标，被加数 colAdd: 矩阵列标，加数 rate: 加数前倍数 */ Matrix& AddCol(int colLocal, int colAdd, double rate = 1.0);
Matrix.cpp函数实现文件：
//*****************矩阵初等变化***************//
Matrix& Matrix::SwapRow(int row0, int row1)
{
//错误判定 越界
if ((this->m_Row <= row0) || (this->m_Col <= row1))
{
std::cout << "Error: <SwapRow> Input row0 Or row1 More Than m_Row" << std::endl;
return *this;
}
else if ((0 > row0) || (0 > row1))
{
std::cout << "Error: <SwapRow> Input row0 Or row1 Less 0" << std::endl;
return *this;
}
else
{
//结果矩阵初始化
Matrix* resMat = new Matrix(*this);
//中转临时变量
double temp = 0.0;
for (int j = 0; j < resMat->m_Col; j++)
{
temp = resMat->m_Matrix[row0][j];
resMat->m_Matrix[row0][j] = resMat->m_Matrix[row1][j];
resMat->m_Matrix[row1][j] = temp;
}
return*resMat;
}
}
Matrix& Matrix::SwapCol(int col0, int col1)
{
//错误判定 越界
if ((this->m_Col <= col0) || (this->m_Col <= col1))
{
std::cout << "Error: <SwapCol> Input col0 Or col1 More Than m_Col" << std::endl;
return *this;
}
else if ((0 > col0) || (0 > col1))
{
std::cout << "Error: <SwapCol> Input col0 Or col1 Less 0" << std::endl;
return *this;
}
else
{
//结果矩阵初始化
Matrix* resMat = new Matrix(*this);
//中转临时变量
double temp = 0.0;
for (int i = 0; i < resMat->m_Row; i++)
{
temp = resMat->m_Matrix[i][col0];
resMat->m_Matrix[i][col0] = resMat->m_Matrix[i][col1];
resMat->m_Matrix[i][col1] = temp;
}
return*resMat;
}
}
//矩阵加法 某行 + 倍数*某行
Matrix& Matrix::AddRow(int rowLocal, int rowAdd, double rate)
{
if ((this->m_Row <= rowLocal) || (this->m_Row <= rowAdd))
{
std::cout << "Error: <AddRow> Input rowLocal Or rowAdd More Than m_Row" << std::endl;
return *this;
}
else if ((0 > rowLocal) || (0 > rowAdd))
{
std::cout << "Error: <AddRow> Input rowLocal Or rowAdd Less 0" << std::endl;
return *this;
}
else
{
//结果矩阵初始化
Matrix* resMat = new Matrix(*this);
//指定行相加
for (int j = 0; j < resMat->m_Col; j++)
{
resMat->m_Matrix[rowLocal][j] += rate * resMat->m_Matrix[rowAdd][j];
}
return *resMat;
}
}
//矩阵加法 某列 + 倍数*某列
Matrix& Matrix::AddCol(int colLocal, int colAdd, double rate)
{
if ((this->m_Col <= colLocal) || (this->m_Col <= colAdd))
{
std::cout << "Error: <AddCol> Input colLocal Or colAdd More Than m_Col" << std::endl;
return *this;
}
else if ((0 > colLocal) || (0 > colAdd))
{
std::cout << "Error: <AddCol> Input colLocal Or colAdd Less 0" << std::endl;
return *this;
}
else
{
//结果矩阵初始化
Matrix* resMat = new Matrix(*this);
//指定列相加
for (int i = 0; i < resMat->m_Row; i++)
{
resMat->m_Matrix[i][colLocal] += rate * resMat->m_Matrix[i][colAdd];
}
return *resMat;
}
}
1.7、矩阵加法
实现矩阵基本加法，注意返回值类型为自身的引用，可实现链式编程。
Matrix.h声明文件：
//*******************矩阵加法*****************// /* 函数名称: 矩阵加法 本矩阵 = 本矩阵 + mat 前提是两个矩阵维度一致 mat: 加数矩阵 */ Matrix& AddMat(Matrix& mat);
Matrix.cpp函数实现文件：
//*******************矩阵加法*****************//
Matrix& Matrix::AddMat(Matrix& mat)
{
Matrix* ResMat = new Matrix(*this);
for (int i = 0; i < ResMat->m_Row; i++)
{
for (int j = 0; j < ResMat->m_Col; j++)
{
ResMat->m_Matrix[i][j] += mat.m_Matrix[i][j];
}
}
return *ResMat;
}
1.8、矩阵乘法
实现矩阵基本乘法，注意返回值类型为自身的引用，可实现链式编程。
Matrix.h声明文件：
//*******************矩阵乘法*****************// /* 函数名称: 矩阵乘法 本矩阵 = 本矩阵*num num: 矩阵乘数 */ Matrix& MultNum(double num); /* 函数名称: 矩阵乘法（运算符重载） 本矩阵 = 本矩阵*num num: 矩阵乘数 */ Matrix& operator * (double num); /* 函数名称: 矩阵某行乘数值row = row*num num: 矩阵某列乘数 row: 矩阵行标 */ Matrix& MultRow(double num, int row); /* 函数名称: 矩阵某列乘数值col = col *num num: 矩阵某列乘数 col: 矩阵列标 */ Matrix& MultCol(double num, int col); /* 函数名称: 矩阵乘法，按照矩阵相乘规则 inputMat: 乘数矩阵 */ Matrix& MultMat(Matrix& inputMat);
Matrix.cpp函数实现文件：
//*******************矩阵乘法*****************//
//矩阵数乘
Matrix& Matrix::MultNum(double num)
{
//结果矩阵初始化
Matrix* resMat = new Matrix(this->m_Row, this->m_Col);
//乘后矩阵生成
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < this->m_Col; j++)
{
resMat->m_Matrix[i][j] = num * this->m_Matrix[i][j];
}
}
return *resMat;
}
//运算符重载 矩阵数乘
Matrix& Matrix::operator*(double num)
{
//结果矩阵初始化
Matrix* resMat = new Matrix(this->m_Row, this->m_Col);
//乘后矩阵生成
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < this->m_Col; j++)
{
resMat->m_Matrix[i][j] = num * this->m_Matrix[i][j];
}
}
return *resMat;
}
//矩阵某行乘数值 行标从0开始计数
Matrix& Matrix::MultRow(double num, int row)
{
if (this->m_Row <= row)
{
std::cout << "Error: <MultRow> Input row More Than m_Row" << std::endl;
return *this;
}
else if (0 > row)
{
std::cout << "Error: <MultRow> Input row Less 0" << std::endl;
return *this;
}
else
{
//结果矩阵初始化
Matrix* resMat = new Matrix(*this);
//乘后矩阵生成
for (int j = 0; j < this->m_Col; j++)
{
resMat->m_Matrix[row][j] = num * this->m_Matrix[row][j];
}
return *resMat;
}
}
//矩阵某列乘数值 列标从0开始计数
Matrix& Matrix::MultCol(double num, int col)
{
if (this->m_Col <= col)
{
std::cout << "Error: <MultCol> Input col More Than m_Row" << std::endl;
return *this;
}
else if (0 > col)
{
std::cout << "Error: <MultCol> Input col Less 0" << std::endl;
return *this;
}
else
{
//结果矩阵初始化
Matrix* resMat = new Matrix(*this);
//乘后矩阵生成
for (int i = 0; i < this->m_Row; i++)
{
resMat->m_Matrix[i][col] = num * this->m_Matrix[i][col];
}
return *resMat;
}
}
//矩阵相乘
Matrix& Matrix::MultMat(Matrix& inputMat)
{
Matrix *resMat = new Matrix(this->m_Row, inputMat.m_Col);
if (this->m_Col != inputMat.m_Row)
{
std::cout << "Matrix Mult Error!" << std::endl;
return *resMat;
}
else
{
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < inputMat.m_Col; j++)
{
for (int k = 0; k < this->m_Col; k++)
{
resMat->m_Matrix[i][j] += this->m_Matrix[i][k] * inputMat.m_Matrix[k][j];
}
}
}
return *resMat;
}
}
1.9、行列式相关操作
实现行列式计算相关操作。
Matrix.h声明文件：
//******************行列式相关操作***********************// /* 函数名称: 求解矩阵对应行列式数值，前提为方阵，按照定义求解，时间复杂度为O(n!*n)，一般不用此方法求解 */ double Det(); /* 函数名称: 求解矩阵对应行列式的顺序主子式，前提为方阵，按照定义求解，时间复杂度为O(n!*n)，一般不用此方法求解 order: 阶数 */ double Det(int order); /* 函数名称: 矩阵行标为row、列标为col的余子式 row: 矩阵行标 col: 矩阵列标 */ Matrix& ChildMatrix(int row, int col); /* 函数名称: 通过高斯列主消元求解矩阵行列式数值，最为常用 */ double DetRow();
Matrix.cpp函数实现文件：
//矩阵的行列式数值
double Matrix::Det()
{
double res = 0.0;
int sign = 1;
if (this->m_Row != this->m_Col)
{
//错误判定
std::cout << "Error: <Det> Matrix Col != Row" << std::endl;
return 0;
}
else if (this->m_Row <= 1)
{
//程序终止出口
return this->m_Matrix[0][0];
}
else
{
for (int i = 0; i < this->m_Col; i++)
{
Matrix* temp = &(this->ChildMatrix(0, i));
res += sign * this->m_Matrix[0][i] * (temp->Det());
sign = -1*sign;
delete temp;
}
}
}
//矩阵行列式顺序主子式 order阶数
double Matrix::Det(int order)
{
if (this->m_Row != this->m_Col)
{
//错误判定
std::cout << "Error: <Det> Matrix Col != Row" << std::endl;
return 0;
}
else if (order < 0)
{
std::cout << "Error: <Det> Input Order Less 0" << std::endl;
return 0;
}
else if (order >= this->m_Row)
{
std::cout << "Error: <Det> Input Order More Than Row" << std::endl;
return 0;
}
else
{
Matrix tempMat(order + 1, order + 1);
for (int i = 0; i < tempMat.m_Col; i++)
{
for (int j = 0; j < tempMat.m_Row; j++)
{
tempMat.m_Matrix[i][j] = this->m_Matrix[i][j];
}
}
return tempMat.Det();
}
}
//求解余子式
Matrix& Matrix::ChildMatrix(int row, int col)
{
if (this->m_Row != this->m_Col)
{
std::cout << "Error: <ChildMatrix> Matrix row != col" << std::endl;
return *this;
}
else if (this->m_Row <= 1)
{
std::cout << "Error: <ChildMatrix> Matrix Row Less 1 " << std::endl;
return *this;
}
else if ((row > this->m_Row) || (col > this->m_Col))
{
std::cout << "Error: <ChildMatrix> Input Row Or Col More Than Matix Max Row Or Col" << std::endl;
return* this;
}
else
{
Matrix* resMat = new Matrix(this->m_Row-1, this->m_Col-1);
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < this->m_Col; j++)
{
if ((i < row) && (j < col))
resMat->m_Matrix[i][j] = this->m_Matrix[i][j];
else if((i > row) && (j < col))
resMat->m_Matrix[i-1][j] = this->m_Matrix[i][j];
else if((i < row) && (j > col))
resMat->m_Matrix[i][j - 1] = this->m_Matrix[i][j];
else if((i > row) && (j > col))
resMat->m_Matrix[i - 1][j - 1] = this->m_Matrix[i][j];
}
}
return *resMat;
}
}
//列主消元处理为上三角矩阵
double Matrix::DetRow()
{
//交换标志位 1代表偶数次交换 -1代表奇数次交换
int flagShift = 1;
//本矩阵
Matrix *localMat = new Matrix(*this);
//行列式数值
double resDet = 1.0;
//*******************通过交换 num1*i + num2*j 实现下三角为0***************//
for (int i = 0; i < localMat->m_Row - 1; i++)
{
//记录最大行所在行标
int tempMaxRow = i;
for (int i1 = i + 1; i1 < localMat->m_Row; i1++)
{
if (abs(localMat->m_Matrix[i1][i]) > abs(localMat->m_Matrix[tempMaxRow][i]))
{
tempMaxRow = i1;
}
}
if (tempMaxRow != i)
{
//std::cout << i << " 行交换" << tempMaxRow << " 行" << std::endl;
//进行交换 将当前第i行与第tempMaxRow行进行互换 初等行变换
*localMat = localMat->SwapRow(i, tempMaxRow);
//记录交换次数
flagShift = -flagShift;
//localMat->PrintMat();
}
//此对角线以下的元素通过初等变化为0
for (int i2 = i + 1; i2 < localMat->m_Row; i2++)
{
if (localMat->m_Matrix[i2][i] != 0)
{
//std::cout << "<" << localMat->m_Matrix[i][i] << "> *" << i2 << " 行 + <" << -1.0 * (localMat->m_Matrix[i2][i]) << "> *" << i << " 行" << std::endl;
*localMat = localMat->AddRow(i2, i, -1.0 * (localMat->m_Matrix[i2][i]) / localMat->m_Matrix[i][i]);
//localMat->PrintMat();
}
}
}
//计算行列式数值 对角线相乘
for (int i = 0; i < localMat->m_Row; i++)
{
resDet = resDet * localMat->m_Matrix[i][i];
}
//矩阵交换一次就会变号
resDet = flagShift * resDet;
//清理localMatrix
delete localMat;
return resDet;
}
1.10、矩阵求逆
实现矩阵求逆相关操作
Matrix.h声明文件：
//*********************矩阵求逆********************// /* 函数名称: 矩阵求逆，按照定义求解，1/|A|*(A*)，时间复杂度为O(n!*n)，一般不用此方法 */ Matrix& Inverse(); /* 函数名称: 矩阵求逆，通过行初等变化，高斯列主消元法求解 */ Matrix& InverseRow(); /* 函数名称: 矩阵求逆，只针对于下三角矩阵进行求解 */ Matrix& InverseDownTriangle(); /* 函数名称: 矩阵求逆，只针对于上三角矩阵进行求解 */ Matrix& InverseUpTriangle(); //矩阵LU分解 /* 函数名称: 矩阵LU分解 LMat: 矩阵分解后的L矩阵 UMat: 矩阵分解后的U矩阵 */ void ResolveLU(Matrix& LMat, Matrix& UMat); /* 函数名称: 矩阵的LUP分解 P*A = L*U 添加了列主消元功能 LMat: 矩阵分解后的L矩阵 UMat: 矩阵分解后的U矩阵 PMat: 矩阵分解后的P矩阵 */ void ResolveLUP(Matrix& LMat, Matrix& UMat, Matrix& PMat);
Matrix.cpp函数实现文件：
//矩阵求逆
Matrix& Matrix::Inverse()
{
if (abs(this->DetRow()) < MIN_DET)
{
std::cout << "Error: <Inverse> Matrix Det Near 0" << std::endl;
return *this;
}
else
{
Matrix* resMat = new Matrix(this->m_Row, this->m_Col);
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < this->m_Col; j++)
{
Matrix* temp = &(this->ChildMatrix(j, i));
resMat->m_Matrix[i][j] = pow(-1.0, (i + j)) / this->DetRow() * (temp->DetRow());
delete temp;
}
}
return *resMat;
}
}
//矩阵求逆 行初等变化
Matrix& Matrix::InverseRow()
{
//错误判断
if (abs(this->DetRow()) < MIN_DET)
{
std::cout << "Error: <InverseRow> Matrix Det Near 0" << std::endl;
return *this;
}
else if (this->m_Row <= 1)
{
std::cout << "Error: <InverseRow> Size Less 2" << std::endl;
return *this;
}
else
{
//单位矩阵 与带转换矩阵维度相同的
Matrix uint = this->Uint();
//结果矩阵 逆矩阵 初始状态与本矩阵相同 为不使本矩阵发生改变
Matrix temp(this->m_Row, this->m_Col);
Matrix* resMat = new Matrix(temp.Uint());
//本矩阵
Matrix localMat(*this);
//*******************通过交换 num1*i + num2*j 实现下三角为0***************//
for (int i = 0; i < localMat.m_Row - 1; i++)
{
//记录最大行所在行标
int tempMaxRow = i;
for (int i1 = i + 1; i1 < localMat.m_Row; i1++)
{
if (abs(localMat.m_Matrix[i1][i]) > abs(localMat.m_Matrix[tempMaxRow][i]))
{
tempMaxRow = i1;
}
}
if (tempMaxRow != i)
{
//std::cout << i << " 行交换" << tempMaxRow << " 行" << std::endl;
//进行交换 将当前第i行与第tempMaxRow行进行互换 初等行变换
localMat = localMat.SwapRow(i, tempMaxRow);
*resMat = resMat->SwapRow(i, tempMaxRow);
//localMat.PrintMat();
}
//此对角线以下的元素通过初等变化为0
for (int i2 = i + 1; i2 < localMat.m_Row; i2++)
{
if (localMat.m_Matrix[i2][i] != 0)
{
//std::cout << "<" << localMat.m_Matrix[i][i] << "> *" << i2 << " 行 + <" << -1.0 * (localMat.m_Matrix[i2][i]) << "> *" << i << " 行" << std::endl;
*resMat = resMat->AddRow(i2, i, -1.0 * (localMat.m_Matrix[i2][i]) / localMat.m_Matrix[i][i]);
localMat = localMat.AddRow(i2, i, -1.0 * (localMat.m_Matrix[i2][i]) / localMat.m_Matrix[i][i]);
//localMat.PrintMat();
}
}
}
//错误判断
if (localMat.m_Matrix[localMat.m_Row - 1][localMat.m_Col - 1] == 0)
{
std::cout << "Error: <InverseRow> marix[" << localMat.m_Row - 1 << "][" << localMat.m_Col - 1 <<"] == 0" << std::endl;
return *this;
}
//*******************通过 num1*i + num2*j 实现上三角为0***************//
for (int i = localMat.m_Row - 1; i > 0; i--)
{
for (int i2 = i - 1; i2 >= 0; i2--)
{
if (localMat.m_Matrix[i2][i] != 0)
{
//std::cout << "<" << localMat.m_Matrix[i][i] << "> *" << i2 << " 行 + <" << -1.0 * (localMat.m_Matrix[i2][i]) << "> *" << i << " 行" << std::endl;
*resMat = resMat->AddRow(i2, i, -1.0 * (localMat.m_Matrix[i2][i]) / localMat.m_Matrix[i][i]);
localMat = localMat.AddRow(i2, i, -1.0 * (localMat.m_Matrix[i2][i]) / localMat.m_Matrix[i][i]);
//localMat.PrintMat();
}
}
}
//*******************通过 i*num 实现矩阵为单位矩阵***************//
for (int i = 0; i < localMat.m_Row; i++)
{
if (localMat.m_Matrix[i][i] == 0)
{
std::cout << "Error: <InverseRow> matrix[" << i << "]" << "[" << i << "] == 0" << std::endl;
return *this;
}
else
{
//std::cout << "<" << 1 / localMat.m_Matrix[i][i] << "> *" << i << " 行" << std::endl;
*resMat = resMat->MultRow(1 / localMat.m_Matrix[i][i], i);
localMat = localMat.MultRow(1 / localMat.m_Matrix[i][i], i);
//localMat.PrintMat();
}
}
return *resMat;
}
}
//矩阵求逆 下三角矩阵
Matrix& Matrix::InverseDownTriangle()
{
//错误判断 方阵检测
if (this->m_Row != this->m_Col)
{
std::cout << "Error: <InverseDownTriangle> Matrix Col != Row" << std::endl;
return *this;
}
//下三角求逆
Matrix* resMat = new Matrix(*this);
for (int i = 0; i < resMat->m_Row; i++)
{
for (int j = 0; j <= i; j++)
{
//分段求解 对角线为倒数
if (i == j)
{
resMat->m_Matrix[i][j] = 1 / resMat->m_Matrix[i][j];
}
else
{
//分段求解 非对角线元素
double tempSum = 0.0;
for (int k = j; k <= i - 1; k++)
{
tempSum += resMat->m_Matrix[i][k] * resMat->m_Matrix[k][j];
}
resMat->m_Matrix[i][j] = -1.0*tempSum / resMat->m_Matrix[i][i];
}
}
}
return *resMat;
}
//矩阵求逆 上三角矩阵
Matrix& Matrix::InverseUpTriangle()
{
//错误判断 方阵检测
if (this->m_Row != this->m_Col)
{
std::cout << "Error: <InverseUpTriangle> Matrix Col != Row" << std::endl;
return *this;
}
//上三角求逆
Matrix* resMat = new Matrix(*this);
for (int j = resMat->m_Col-1; j >=0; j--)
{
for (int i = j; i >=0; i--)
{
//分段求解 对角线为倒数
if (i == j)
{
resMat->m_Matrix[i][j] = 1 / resMat->m_Matrix[i][j];
}
else
{
//分段求解 非对角线元素
double tempSum = 0.0;
for (int k = j; k >= i+1; k--)
{
tempSum += resMat->m_Matrix[i][k] * resMat->m_Matrix[k][j];
}
resMat->m_Matrix[i][j] = -1.0 * tempSum / resMat->m_Matrix[i][i];
}
}
}
return *resMat;
}
//矩阵LU分解 顺序分解 对于病态矩阵可能存在精度问题
void Matrix::ResolveLU(Matrix& LMat, Matrix& UMat)
{
if (this->m_Col != this->m_Row)
{
std::cout << "Error: <ResolveLU> Is Not Square Matrix" << std::endl;
return;
}
//存在性判定 顺序主子式不为0
for (int i = 0; i < this->m_Row; i++)
{
if (this->Det(i) == 0)
{
std::cout << "Error: <ResolveLU> order Det = 0" << std::endl;
return;
}
}
//LU 分解
//L矩阵为单位矩阵
LMat = this->Uint();
//U矩阵初始化为空矩阵
Matrix temp(this->m_Row, this->m_Col);
UMat = temp;
for (int i = 0; i < this->m_Row; i++)
{
//计算U
for (int j1 = i; j1 < this->m_Col; j1++)
{
double tempSum1 = 0.0;
if (i != 0)
{
for (int j2 = 0; j2 <= i - 1; j2++)
{
tempSum1 += LMat.m_Matrix[i][j2] * UMat.m_Matrix[j2][j1];
}
}
UMat.m_Matrix[i][j1] = this->m_Matrix[i][j1] - tempSum1;
}
//计算L
for (int i1 = i; i1 < this->m_Row; i1++)
{
double tempSum2 = 0.0;
if (i != 0)
{
for (int j2 = 0; j2 <= i - 1; j2++)
{
tempSum2 += LMat.m_Matrix[i1][j2] * UMat.m_Matrix[j2][i];
}
}
LMat.m_Matrix[i1][i] = (this->m_Matrix[i1][i] - tempSum2)/UMat.m_Matrix[i][i];
}
}
}
//矩阵的LUP分解 P*A = L*U 添加了列主消元功能
//L为主对角线元素为1的下三角矩阵 U为上二角矩阵 P为行交换矩阵 P*A=L*U
void Matrix::ResolveLUP(Matrix& LMat, Matrix& UMat, Matrix& PMat)
{
//条件判断 矩阵行列式不为0
if (this->Det() == 0)
{
std::cout << "Error: <ResolveLUP> Can't Resolve Matrix To L U P" << std::endl;
return;
}
//初始化 L U P
LMat = this->Uint();
PMat = this->Uint();
UMat = *this;
//进行分解计算
for (int i = 0; i < UMat.m_Row - 1; i++)
{
//记录最大行所在行标
int tempMaxRow = i;
for (int i1 = i + 1; i1 < UMat.m_Row; i1++)
{
if (abs(UMat.m_Matrix[i1][i]) > abs(UMat.m_Matrix[tempMaxRow][i]))
{
tempMaxRow = i1;
}
}
//进行交换 将当前第i行与第tempMaxRow行进行互换 初等行变换
UMat = UMat.SwapRow(i, tempMaxRow);
//L矩阵做出对应交换 先交换<itempMaxRow>列再交换<itempMaxRow>行
LMat = LMat.SwapCol(i, tempMaxRow);
LMat = LMat.SwapRow(i, tempMaxRow);
//P矩阵做出对应变换 交换<itempMaxRow>行
PMat = PMat.SwapRow(i, tempMaxRow);
//高斯消元 V矩阵消除下三角区域，L矩阵添加下三角区域
for (int i1 = i + 1; i1 < UMat.m_Row; i1++)
{
//记录消元系数
double deleteVar = UMat.m_Matrix[i1][i] / UMat.m_Matrix[i][i];
//L矩阵列填充
LMat.m_Matrix[i1][i] = deleteVar;
//U矩阵列消除
UMat = UMat.MultRow(UMat.m_Matrix[i][i], i1).AddRow(i1, i, -1.0 * UMat.m_Matrix[i1][i]).MultRow(1 / UMat.m_Matrix[i][i], i1);
}
}
return;
}
2、private variable
double** m_Matrix; //矩阵 int m_Row; //矩阵行数 int m_Col; //矩阵列数
3、全部源码
为了方便大家复制应用，这里直接贴出源码。
Matrix.h声明文件：
#ifndef _MATRIX_H_
#define _MATRIX_H_
#include <iostream>
#include <math.h>
#include <vector>
//矩阵最大容量
#define MAX_COUNT 500
#define MIN_DET 1e-12 //行列式最小数值
class Matrix
{
public:
//******************************构造函数与析构函数********************************//
/*
函数名称: 无参构造函数
*/
Matrix();
/*
函数名称: 矩阵有参构造函数，初始化为row行、col列的0矩阵
row: 矩阵行数
col: 矩阵列数
*/
Matrix(int row, int col);
/*
函数名称: 矩阵有参构造函数，初始化为row行、col列、数值为mat的矩阵
row: 矩阵行数
col: 矩阵列数
*mat: 矩阵数值一维数组
*/
Matrix(int row, int col, double* mat);
/*
函数名称: 深拷贝构造函数
mat: 需要复制的矩阵
*/
Matrix(const Matrix& mat);
/*
函数名称: 析构函数
*/
~Matrix();
//*******************获取矩阵*****************//
/*
函数名称: 获取矩阵的第row行、第col列元素数值
row: 矩阵行数
col: 矩阵列数
*/
double GetMatrixEle(int row, int col);
//*******************设置矩阵*****************//
/*
函数名称: 设置矩阵第row行、第col列数值
row: 矩阵行数
col: 矩阵列数
value: 设置的矩阵数值
*/
void SetMatrixEle(int row, int col, double value);
/*
函数名称: 深拷贝矩阵
mat: 需要复制的矩阵
*/
Matrix CopyMat(const Matrix mat);
//********************************矩阵的相关计算**********************************//
//*******************打印矩阵*****************//
/*
函数名称: 打印矩阵
*/
void PrintMat();
//*****************矩阵基本操作***************//
/*
函数名称: 矩阵转置，返回的是自身引用，可链式调用
*/
Matrix& Transpose();
/*
函数名称: 等维度的单位矩阵，前提是方阵
*/
Matrix& Uint();
//****************矩阵保留与剔除**************//
/*
函数名称: 剔除矩阵中以index为行标和列标的行和列，num代表index的大小
*index: 矩阵中的行号与列号一维数组
num: index动态数组长度
*/
Matrix& DeleteMat(int *index, int num);
/*
函数名称: 剔除矩阵中以index为行标和列标的行和列，num代表index的大小
*index: 矩阵中的行号与列号一维动态数组
num: index动态数组长度
*/
Matrix& DeleteMat(std::vector<int> index, int num);
/*
函数名称: 剔除矩阵中以index为行标的行，num代表index的大小
*index: 矩阵中的行号一维数组
num: index动态数组长度
*/
Matrix& DeleteRow(int* index, int num);
/*
函数名称: 剔除矩阵中以index为行标的行，num代表index的大小
*index: 矩阵中的行号一维动态数组
num: index动态数组长度
*/
Matrix& DeleteRow(std::vector<int> index, int num);
/*
函数名称: 剔除矩阵中以index为列标的列，num代表index的大小
*index: 矩阵中的列号一维数组
num: index动态数组长度
*/
Matrix& DeleteCol(int* index, int num);
/*
函数名称: 剔除矩阵中以index为列标的列，num代表index的大小
*index: 矩阵中的列号一维动态数组
num: index动态数组长度
*/
Matrix& DeleteCol(std::vector<int> index, int num);
//******************矩阵的替换****************//
/*
函数名称: 替换矩阵中行标和列标为 index中的行与列，num代表index的大小, mat是需要替换的矩阵
*index: 矩阵中的行标和列标的一维数组
num: index动态数组长度
mat: 需要替换的矩阵
*/
Matrix& ReplaceMat(int* index, int num, Matrix& mat);
/*
函数名称: 替换矩阵中行标和列标为 index中的行与列，num代表index的大小, mat是需要替换的矩阵
*index: 矩阵中的行标和列标的一维动态数组
num: index动态数组长度
mat: 需要替换的矩阵
*/
Matrix& ReplaceMat(std::vector<int> index, int num, Matrix& mat);
/*
函数名称: 替换矩阵中行标为 index中的行，num代表index的大小, mat是需要替换的矩阵
*index: 矩阵中的行标的一维数组
num: index动态数组长度
mat: 需要替换的矩阵
*/
Matrix& ReplaceRow(int* index, int num, Matrix& mat);
/*
函数名称: 替换矩阵中行标为 index中的行，num代表index的大小, mat是需要替换的矩阵
*index: 矩阵中的行标的一动态维数组
num: index动态数组长度
mat: 需要替换的矩阵
*/
Matrix& ReplaceRow(std::vector<int> index, int num, Matrix& mat);
/*
函数名称: 替换矩阵中列标为 index中的列，num代表index的大小, mat是需要替换的矩阵
*index: 矩阵中的列标的一维数组
num: index动态数组长度
mat: 需要替换的矩阵
*/
Matrix& ReplaceCol(int* index, int num, Matrix& mat);
/*
函数名称: 替换矩阵中列标为 index中的列，num代表index的大小, mat是需要替换的矩阵
*index: 矩阵中的列标的一维动态数组
num: index动态数组长度
mat: 需要替换的矩阵
*/
Matrix& ReplaceCol(std::vector<int> index, int num, Matrix& mat);
//*****************矩阵初等变化***************//
/*
函数名称: 交换矩阵中行标为row0与row1的元素
row0: 矩阵行标0
row1: 矩阵行标1
*/
Matrix& SwapRow(int row0, int row1);
/*
函数名称: 交换矩阵中列标为col0与col1的元素
col0: 矩阵列标0
col1: 矩阵列标1
*/
Matrix& SwapCol(int col0, int col1);
/*
函数名称: 矩阵行加法 rowLocal = rowLocal + rate *rowAdd
rowLocal: 矩阵行标，被加数
rowAdd: 矩阵行标，加数
rate: 加数前倍数
*/
Matrix& AddRow(int rowLocal, int rowAdd, double rate = 1.0);
//矩阵加法 某列 + 倍数*某列
/*
函数名称: 矩阵列加法 colLocal = colLocal + rate * colAdd
colLocal: 矩阵列标，被加数
colAdd: 矩阵列标，加数
rate: 加数前倍数
*/
Matrix& AddCol(int colLocal, int colAdd, double rate = 1.0);
//*******************矩阵加法*****************//
/*
函数名称: 矩阵加法 本矩阵 = 本矩阵 + mat 前提是两个矩阵维度一致
mat: 加数矩阵
*/
Matrix& AddMat(Matrix& mat);
//*******************矩阵乘法*****************//
/*
函数名称: 矩阵乘法 本矩阵 = 本矩阵*num
num: 矩阵乘数
*/
Matrix& MultNum(double num);
/*
函数名称: 矩阵乘法（运算符重载） 本矩阵 = 本矩阵*num
num: 矩阵乘数
*/
Matrix& operator * (double num);
/*
函数名称: 矩阵某行乘数值row = row*num
num: 矩阵某列乘数
row: 矩阵行标
*/
Matrix& MultRow(double num, int row);
/*
函数名称: 矩阵某列乘数值col = col *num
num: 矩阵某列乘数
col: 矩阵列标
*/
Matrix& MultCol(double num, int col);
/*
函数名称: 矩阵乘法，按照矩阵相乘规则
inputMat: 乘数矩阵
*/
Matrix& MultMat(Matrix& inputMat);
//******************行列式相关操作***********************//
/*
函数名称: 求解矩阵对应行列式数值，前提为方阵，按照定义求解，时间复杂度为O(n!*n)，一般不用此方法求解
*/
double Det();
/*
函数名称: 求解矩阵对应行列式的顺序主子式，前提为方阵，按照定义求解，时间复杂度为O(n!*n)，一般不用此方法求解
order: 阶数
*/
double Det(int order);
/*
函数名称: 矩阵行标为row、列标为col的余子式
row: 矩阵行标
col: 矩阵列标
*/
Matrix& ChildMatrix(int row, int col);
/*
函数名称: 通过高斯列主消元求解矩阵行列式数值，最为常用
*/
double DetRow();
//*********************矩阵求逆********************//
/*
函数名称: 矩阵求逆，按照定义求解，1/|A|*(A*)，时间复杂度为O(n!*n)，一般不用此方法
*/
Matrix& Inverse();
/*
函数名称: 矩阵求逆，通过行初等变化，高斯列主消元法求解
*/
Matrix& InverseRow();
/*
函数名称: 矩阵求逆，只针对于下三角矩阵进行求解
*/
Matrix& InverseDownTriangle();
/*
函数名称: 矩阵求逆，只针对于上三角矩阵进行求解
*/
Matrix& InverseUpTriangle();
//矩阵LU分解
/*
函数名称: 矩阵LU分解
LMat: 矩阵分解后的L矩阵
UMat: 矩阵分解后的U矩阵
*/
void ResolveLU(Matrix& LMat, Matrix& UMat);
/*
函数名称: 矩阵的LUP分解 P*A = L*U 添加了列主消元功能
LMat: 矩阵分解后的L矩阵
UMat: 矩阵分解后的U矩阵
PMat: 矩阵分解后的P矩阵
*/
void ResolveLUP(Matrix& LMat, Matrix& UMat, Matrix& PMat);
private:
double** m_Matrix; //矩阵
int m_Row; //矩阵行数
int m_Col; //矩阵列数
};
#endif
Matrix.cpp函数实现文件：
#include "Matrix.h"
//******************************构造函数与析构函数********************************//
Matrix::Matrix()
{
}
//初始化矩阵 默认值为0
Matrix::Matrix(int row, int col)
{
this->m_Row = row;
this->m_Col = col;
//开辟内存
this->m_Matrix = new double* [row];
for (int i = 0; i < row; i++)
{
this->m_Matrix[i] = new double[col] {0.0};
}
}
//初始化矩阵 设定数值
Matrix::Matrix(int row, int col, double *mat)
{
this->m_Row = row;
this->m_Col = col;
//开辟内存
this->m_Matrix = new double* [row];
for (int i = 0; i < row; i++)
{
this->m_Matrix[i] = new double[col] {0.0};
}
//矩阵赋值
for(int i = 0; i<row; i++)
{
for (int j = 0; j < col; j++)
{
this->m_Matrix[i][j] = mat[i * col + j];
}
}
}
//深拷贝
Matrix::Matrix(const Matrix& mat)
{
//行列传递
this->m_Row = mat.m_Row;
this->m_Col = mat.m_Col;
//矩阵深拷贝
this->m_Matrix = new double* [this->m_Row];
for (int i = 0; i < this->m_Row; i++)
{
this->m_Matrix[i] = new double[this->m_Col];
memcpy(this->m_Matrix[i], mat.m_Matrix[i], sizeof(double) * this->m_Col);
}
}
Matrix::~Matrix()
{
//释放矩阵每一行
for (int i = 0; i < this->m_Row; i++)
{
if (this->m_Matrix[i] != NULL)
{
delete[]this->m_Matrix[i];
this->m_Matrix[i] = NULL;
}
}
//释放矩阵顶点
if (this->m_Matrix != NULL)
{
delete[]this->m_Matrix;
this->m_Matrix = NULL;
}
}
//获取矩阵某个元素 某行某列
double Matrix::GetMatrixEle(int row, int col)
{
if (row >= this->m_Row)
{
std::cout << "Error: <GetMatrixEle> Input row >= m_Row" << std::endl;
return 0.0;
}
else if (col >= this->m_Col)
{
std::cout << "Error: <GetMatrixEle> Input col >= m_Col" << std::endl;
return 0.0;
}
else
{
return this->m_Matrix[row][col];
}
}
//*******************设置矩阵*****************//
void Matrix::SetMatrixEle(int row, int col, double value)
{
if (row >= this->m_Row)
{
std::cout << "Error: <SetMatrixEle> Input row >= m_Row" << std::endl;
return;
}
else if (col >= this->m_Col)
{
std::cout << "Error: <SetMatrixEle> Input col >= m_Col" << std::endl;
return;
}
else
{
this->m_Matrix[row][col] = value;
return;
}
}
Matrix Matrix::CopyMat(const Matrix mat)
{
//行列传递
this->m_Row = mat.m_Row;
this->m_Col = mat.m_Col;
//矩阵深拷贝
this->m_Matrix = new double* [this->m_Row];
for (int i = 0; i < this->m_Row; i++)
{
this->m_Matrix[i] = new double[this->m_Col];
memcpy(this->m_Matrix[i], mat.m_Matrix[i], sizeof(double) * this->m_Col);
}
return *this;
}
//*******************打印矩阵*****************//
//矩阵输出
void Matrix::PrintMat()
{
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < this->m_Col; j++)
{
std::cout.setf(std::ios::scientific); //科学计数法表示
std::cout << this->m_Matrix[i][j] << "\t";
}
std::cout << std::endl;
}
std::cout << std::endl;
}
//*****************矩阵基本操作***************//
//矩阵转置
Matrix& Matrix::Transpose()
{
Matrix* resMat = new Matrix(this->m_Col, this->m_Row);
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < this->m_Col; j++)
{
resMat->m_Matrix[j][i] = this->m_Matrix[i][j];
}
}
return *resMat;
}
//求等长度单位矩阵
Matrix& Matrix::Uint()
{
//矩阵是否为方阵
if (this->m_Col != this->m_Row)
{
std::cout << "Error: <Uint> Row != Col" << std::endl;
Matrix* resMat = new Matrix(this->m_Row, this->m_Row);
return *resMat;
}
else
{
//单位矩阵初始化
Matrix* resMat = new Matrix(this->m_Row, this->m_Col);
//单位矩阵生成
for (int i = 0; i < this->m_Row; i++)
{
resMat->m_Matrix[i][i] = 1.0;
}
return *resMat;
}
}
//****************矩阵保留与剔除**************//
//剔除矩阵的 index中的行与列，num代表index的大小
Matrix& Matrix::DeleteMat(int* index, int num)
{
//结果矩阵
Matrix* resMat = new Matrix(this->m_Row - num, this->m_Col - num);
int recIndex[MAX_COUNT];
int currIndex = 0;
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
else if (index[i] >= this->m_Col)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Col" << std::endl;
return *this;
}
}
//筛选出剔除后行数
for (int iRow = 0; iRow < this->m_Row; iRow++)
{
for (int iNum = 0; iNum < num; iNum++)
{
if (iRow == index[iNum])
{
break;
}
if (iNum == num-1)
{
recIndex[currIndex++] = iRow;
}
}
}
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[iRow][iCol] = this->m_Matrix[recIndex[iRow]][recIndex[iCol]];
}
}
return *resMat;
}
Matrix& Matrix::DeleteMat(std::vector<int> index, int num)
{
//结果矩阵
Matrix* resMat = new Matrix(this->m_Row - num, this->m_Col - num);
int recIndex[MAX_COUNT];
int currIndex = 0;
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
else if (index[i] >= this->m_Col)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Col" << std::endl;
return *this;
}
}
//筛选出剔除后行数
for (int iRow = 0; iRow < this->m_Row; iRow++)
{
for (int iNum = 0; iNum < num; iNum++)
{
if (iRow == index[iNum])
{
break;
}
if (iNum == num - 1)
{
recIndex[currIndex++] = iRow;
}
}
}
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[iRow][iCol] = this->m_Matrix[recIndex[iRow]][recIndex[iCol]];
}
}
return *resMat;
}
//剔除矩阵的 index中的行，num代表index的大小
Matrix& Matrix::DeleteRow(int* index, int num)
{
//结果矩阵
Matrix* resMat = new Matrix(this->m_Row - num, this->m_Col);
int recIndex[MAX_COUNT];
int currIndex = 0;
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
}
//筛选出剔除后行数
for (int iRow = 0; iRow < this->m_Row; iRow++)
{
for (int iNum = 0; iNum < num; iNum++)
{
if (iRow == index[iNum])
{
break;
}
if (iNum == num - 1)
{
recIndex[currIndex++] = iRow;
}
}
}
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[iRow][iCol] = this->m_Matrix[recIndex[iRow]][iCol];
}
}
return *resMat;
}
Matrix& Matrix::DeleteRow(std::vector<int> index, int num)
{
//结果矩阵
Matrix* resMat = new Matrix(this->m_Row - num, this->m_Col);
int recIndex[MAX_COUNT];
int currIndex = 0;
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
}
//筛选出剔除后行数
for (int iRow = 0; iRow < this->m_Row; iRow++)
{
for (int iNum = 0; iNum < num; iNum++)
{
if (iRow == index[iNum])
{
break;
}
if (iNum == num - 1)
{
recIndex[currIndex++] = iRow;
}
}
}
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[iRow][iCol] = this->m_Matrix[recIndex[iRow]][iCol];
}
}
return *resMat;
}
Matrix& Matrix::DeleteCol(int* index, int num)
{
//结果矩阵
Matrix* resMat = new Matrix(this->m_Row, this->m_Col - num);
int recIndex[MAX_COUNT];
int currIndex = 0;
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
}
//筛选出剔除后行数
for (int iRow = 0; iRow < this->m_Row; iRow++)
{
for (int iNum = 0; iNum < num; iNum++)
{
if (iRow == index[iNum])
{
break;
}
if (iNum == num - 1)
{
recIndex[currIndex++] = iRow;
}
}
}
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[iRow][iCol] = this->m_Matrix[iRow][recIndex[iCol]];
}
}
return *resMat;
}
Matrix& Matrix::DeleteCol(std::vector<int> index, int num)
{
//结果矩阵
Matrix* resMat = new Matrix(this->m_Row, this->m_Col - num);
int recIndex[MAX_COUNT];
int currIndex = 0;
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <DeleteMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
}
//筛选出剔除后行数
for (int iRow = 0; iRow < this->m_Row; iRow++)
{
for (int iNum = 0; iNum < num; iNum++)
{
if (iRow == index[iNum])
{
break;
}
if (iNum == num - 1)
{
recIndex[currIndex++] = iRow;
}
}
}
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[iRow][iCol] = this->m_Matrix[iRow][recIndex[iCol]];
}
}
return *resMat;
}
//******************矩阵的替换****************//
//替换矩阵中的行和列 index中的行与列，num代表index的大小
Matrix& Matrix::ReplaceMat(int* index, int num, Matrix& mat)
{
//错误判定 方阵
if (this->m_Row != this->m_Col)
{
std::cout << "Error: <ReplaceMat> this m_Col != m_Row" << std::endl;
return *this;
}
//检验插入矩阵为方阵
if (mat.m_Row != mat.m_Col)
{
std::cout << "Error: <ReplaceMat> mat m_Col != m_Row" << std::endl;
return *this;
}
//检验插入矩阵大小与num保持一致
if (mat.m_Col != num)
{
std::cout << "Error: <ReplaceMat> num != mat.m_Col" << std::endl;
return *this;
}
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <ReplaceMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
else if (index[i] >= this->m_Col)
{
std::cout << "Error: <ReplaceMat> Input index[" << i << "] = " << index[i] << " >= m_Col" << std::endl;
return *this;
}
}
//结果矩阵
Matrix* resMat = new Matrix(*this);
//加入元素
for (int iRow = 0; iRow < num; iRow++)
{
for (int iCol = 0; iCol < num; iCol++)
{
resMat->m_Matrix[index[iRow]][index[iCol]] = mat.m_Matrix[iRow][iCol];
}
}
return *resMat;
}
Matrix& Matrix::ReplaceMat(std::vector<int> index, int num, Matrix& mat)
{
//错误判定 方阵
if (this->m_Row != this->m_Col)
{
std::cout << "Error: <ReplaceMat> this m_Col != m_Row" << std::endl;
return *this;
}
//检验插入矩阵为方阵
if (mat.m_Row != mat.m_Col)
{
std::cout << "Error: <ReplaceMat> mat m_Col != m_Row" << std::endl;
return *this;
}
//检验插入矩阵大小与num保持一致
if (mat.m_Col != num)
{
std::cout << "Error: <ReplaceMat> num != mat.m_Col" << std::endl;
return *this;
}
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <ReplaceMat> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
else if (index[i] >= this->m_Col)
{
std::cout << "Error: <ReplaceMat> Input index[" << i << "] = " << index[i] << " >= m_Col" << std::endl;
return *this;
}
}
//结果矩阵
Matrix* resMat = new Matrix(*this);
//加入元素
for (int iRow = 0; iRow < num; iRow++)
{
for (int iCol = 0; iCol < num; iCol++)
{
resMat->m_Matrix[index[iRow]][index[iCol]] = mat.m_Matrix[iRow][iCol];
}
}
return *resMat;
}
//替换矩阵中的行 index中的行，num代表index的大小, mat是需要替换的矩阵
Matrix& Matrix::ReplaceRow(int* index, int num, Matrix& mat)
{
//检验插入矩阵大小与num保持一致
if (mat.m_Row != num)
{
std::cout << "Error: <ReplaceRow> num != mat.m_Row" << std::endl;
return *this;
}
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <ReplaceRow> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
}
//当前矩阵列数应与mat列数一致
if (this->m_Col != mat.m_Col)
{
std::cout << "Error: <ReplaceRow> this->m_Col != mat.m_Col" << std::endl;
return *this;
}
//结果矩阵
Matrix* resMat = new Matrix(*this);
//加入元素
for (int iRow = 0; iRow < num; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[index[iRow]][iCol] = mat.m_Matrix[iRow][iCol];
}
}
return *resMat;
}
Matrix& Matrix::ReplaceRow(std::vector<int> index, int num, Matrix& mat)
{
//检验插入矩阵大小与num保持一致
if (mat.m_Row != num)
{
std::cout << "Error: <ReplaceRow> num != mat.m_Row" << std::endl;
return *this;
}
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Row)
{
std::cout << "Error: <ReplaceRow> Input index[" << i << "] = " << index[i] << " >= m_Row" << std::endl;
return *this;
}
}
//当前矩阵列数应与mat列数一致
if (this->m_Col != mat.m_Col)
{
std::cout << "Error: <ReplaceRow> this->m_Col != mat.m_Col" << std::endl;
return *this;
}
//结果矩阵
Matrix* resMat = new Matrix(*this);
//加入元素
for (int iRow = 0; iRow < num; iRow++)
{
for (int iCol = 0; iCol < resMat->m_Col; iCol++)
{
resMat->m_Matrix[index[iRow]][iCol] = mat.m_Matrix[iRow][iCol];
}
}
return *resMat;
}
//替换矩阵中的列 index中的列，num代表index的大小, mat是需要替换的矩阵
Matrix& Matrix::ReplaceCol(int* index, int num, Matrix& mat)
{
//检验插入矩阵大小与num保持一致
if (mat.m_Col != num)
{
std::cout << "Error: <ReplaceCol> mat.m_Col != num" << std::endl;
return *this;
}
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Col)
{
std::cout << "Error: <ReplaceCol> Input index[" << i << "] = " << index[i] << " >= m_Col" << std::endl;
return *this;
}
}
//当前矩阵行数应与mat行数一致
if (this->m_Row != mat.m_Row)
{
std::cout << "Error: <ReplaceCol> this->m_Row != mat.m_Row" << std::endl;
return *this;
}
//结果矩阵
Matrix* resMat = new Matrix(*this);
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < num; iCol++)
{
resMat->m_Matrix[iRow][index[iCol]] = mat.m_Matrix[iRow][iCol];
}
}
return *resMat;
}
Matrix& Matrix::ReplaceCol(std::vector<int> index, int num, Matrix& mat)
{
//检验插入矩阵大小与num保持一致
if (mat.m_Col != num)
{
std::cout << "Error: <ReplaceCol> mat.m_Col != num" << std::endl;
return *this;
}
//检验数据有效性
for (int i = 0; i < num; i++)
{
//越界判定
if (index[i] >= this->m_Col)
{
std::cout << "Error: <ReplaceCol> Input index[" << i << "] = " << index[i] << " >= m_Col" << std::endl;
return *this;
}
}
//当前矩阵行数应与mat行数一致
if (this->m_Row != mat.m_Row)
{
std::cout << "Error: <ReplaceCol> this->m_Row != mat.m_Row" << std::endl;
return *this;
}
//结果矩阵
Matrix* resMat = new Matrix(*this);
//加入元素
for (int iRow = 0; iRow < resMat->m_Row; iRow++)
{
for (int iCol = 0; iCol < num; iCol++)
{
resMat->m_Matrix[iRow][index[iCol]] = mat.m_Matrix[iRow][iCol];
}
}
return *resMat;
}
//*****************矩阵初等变化***************//
Matrix& Matrix::SwapRow(int row0, int row1)
{
//错误判定 越界
if ((this->m_Row <= row0) || (this->m_Col <= row1))
{
std::cout << "Error: <SwapRow> Input row0 Or row1 More Than m_Row" << std::endl;
return *this;
}
else if ((0 > row0) || (0 > row1))
{
std::cout << "Error: <SwapRow> Input row0 Or row1 Less 0" << std::endl;
return *this;
}
else
{
//结果矩阵初始化
Matrix* resMat = new Matrix(*this);
//中转临时变量
double temp = 0.0;
for (int j = 0; j < resMat->m_Col; j++)
{
temp = resMat->m_Matrix[row0][j];
resMat->m_Matrix[row0][j] = resMat->m_Matrix[row1][j];
resMat->m_Matrix[row1][j] = temp;
}
return*resMat;
}
}
Matrix& Matrix::SwapCol(int col0, int col1)
{
//错误判定 越界
if ((this->m_Col <= col0) || (this->m_Col <= col1))
{
std::cout << "Error: <SwapCol> Input col0 Or col1 More Than m_Col" << std::endl;
return *this;
}
else if ((0 > col0) || (0 > col1))
{
std::cout << "Error: <SwapCol> Input col0 Or col1 Less 0" << std::endl;
return *this;
}
else
{
//结果矩阵初始化
Matrix* resMat = new Matrix(*this);
//中转临时变量
double temp = 0.0;
for (int i = 0; i < resMat->m_Row; i++)
{
temp = resMat->m_Matrix[i][col0];
resMat->m_Matrix[i][col0] = resMat->m_Matrix[i][col1];
resMat->m_Matrix[i][col1] = temp;
}
return*resMat;
}
}
//矩阵加法 某行 + 倍数*某行
Matrix& Matrix::AddRow(int rowLocal, int rowAdd, double rate)
{
if ((this->m_Row <= rowLocal) || (this->m_Row <= rowAdd))
{
std::cout << "Error: <AddRow> Input rowLocal Or rowAdd More Than m_Row" << std::endl;
return *this;
}
else if ((0 > rowLocal) || (0 > rowAdd))
{
std::cout << "Error: <AddRow> Input rowLocal Or rowAdd Less 0" << std::endl;
return *this;
}
else
{
//结果矩阵初始化
Matrix* resMat = new Matrix(*this);
//指定行相加
for (int j = 0; j < resMat->m_Col; j++)
{
resMat->m_Matrix[rowLocal][j] += rate * resMat->m_Matrix[rowAdd][j];
}
return *resMat;
}
}
//矩阵加法 某列 + 倍数*某列
Matrix& Matrix::AddCol(int colLocal, int colAdd, double rate)
{
if ((this->m_Col <= colLocal) || (this->m_Col <= colAdd))
{
std::cout << "Error: <AddCol> Input colLocal Or colAdd More Than m_Col" << std::endl;
return *this;
}
else if ((0 > colLocal) || (0 > colAdd))
{
std::cout << "Error: <AddCol> Input colLocal Or colAdd Less 0" << std::endl;
return *this;
}
else
{
//结果矩阵初始化
Matrix* resMat = new Matrix(*this);
//指定列相加
for (int i = 0; i < resMat->m_Row; i++)
{
resMat->m_Matrix[i][colLocal] += rate * resMat->m_Matrix[i][colAdd];
}
return *resMat;
}
}
//*******************矩阵加法*****************//
Matrix& Matrix::AddMat(Matrix& mat)
{
Matrix* ResMat = new Matrix(*this);
for (int i = 0; i < ResMat->m_Row; i++)
{
for (int j = 0; j < ResMat->m_Col; j++)
{
ResMat->m_Matrix[i][j] += mat.m_Matrix[i][j];
}
}
return *ResMat;
}
//*******************矩阵乘法*****************//
//矩阵数乘
Matrix& Matrix::MultNum(double num)
{
//结果矩阵初始化
Matrix* resMat = new Matrix(this->m_Row, this->m_Col);
//乘后矩阵生成
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < this->m_Col; j++)
{
resMat->m_Matrix[i][j] = num * this->m_Matrix[i][j];
}
}
return *resMat;
}
//运算符重载 矩阵数乘
Matrix& Matrix::operator*(double num)
{
//结果矩阵初始化
Matrix* resMat = new Matrix(this->m_Row, this->m_Col);
//乘后矩阵生成
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < this->m_Col; j++)
{
resMat->m_Matrix[i][j] = num * this->m_Matrix[i][j];
}
}
return *resMat;
}
//矩阵某行乘数值 行标从0开始计数
Matrix& Matrix::MultRow(double num, int row)
{
if (this->m_Row <= row)
{
std::cout << "Error: <MultRow> Input row More Than m_Row" << std::endl;
return *this;
}
else if (0 > row)
{
std::cout << "Error: <MultRow> Input row Less 0" << std::endl;
return *this;
}
else
{
//结果矩阵初始化
Matrix* resMat = new Matrix(*this);
//乘后矩阵生成
for (int j = 0; j < this->m_Col; j++)
{
resMat->m_Matrix[row][j] = num * this->m_Matrix[row][j];
}
return *resMat;
}
}
//矩阵某列乘数值 列标从0开始计数
Matrix& Matrix::MultCol(double num, int col)
{
if (this->m_Col <= col)
{
std::cout << "Error: <MultCol> Input col More Than m_Row" << std::endl;
return *this;
}
else if (0 > col)
{
std::cout << "Error: <MultCol> Input col Less 0" << std::endl;
return *this;
}
else
{
//结果矩阵初始化
Matrix* resMat = new Matrix(*this);
//乘后矩阵生成
for (int i = 0; i < this->m_Row; i++)
{
resMat->m_Matrix[i][col] = num * this->m_Matrix[i][col];
}
return *resMat;
}
}
//矩阵相乘
Matrix& Matrix::MultMat(Matrix& inputMat)
{
Matrix *resMat = new Matrix(this->m_Row, inputMat.m_Col);
if (this->m_Col != inputMat.m_Row)
{
std::cout << "Matrix Mult Error!" << std::endl;
return *resMat;
}
else
{
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < inputMat.m_Col; j++)
{
for (int k = 0; k < this->m_Col; k++)
{
resMat->m_Matrix[i][j] += this->m_Matrix[i][k] * inputMat.m_Matrix[k][j];
}
}
}
return *resMat;
}
}
//矩阵的行列式数值
double Matrix::Det()
{
double res = 0.0;
int sign = 1;
if (this->m_Row != this->m_Col)
{
//错误判定
std::cout << "Error: <Det> Matrix Col != Row" << std::endl;
return 0;
}
else if (this->m_Row <= 1)
{
//程序终止出口
return this->m_Matrix[0][0];
}
else
{
for (int i = 0; i < this->m_Col; i++)
{
Matrix* temp = &(this->ChildMatrix(0, i));
res += sign * this->m_Matrix[0][i] * (temp->Det());
sign = -1*sign;
delete temp;
}
}
}
//矩阵行列式顺序主子式 order阶数
double Matrix::Det(int order)
{
if (this->m_Row != this->m_Col)
{
//错误判定
std::cout << "Error: <Det> Matrix Col != Row" << std::endl;
return 0;
}
else if (order < 0)
{
std::cout << "Error: <Det> Input Order Less 0" << std::endl;
return 0;
}
else if (order >= this->m_Row)
{
std::cout << "Error: <Det> Input Order More Than Row" << std::endl;
return 0;
}
else
{
Matrix tempMat(order + 1, order + 1);
for (int i = 0; i < tempMat.m_Col; i++)
{
for (int j = 0; j < tempMat.m_Row; j++)
{
tempMat.m_Matrix[i][j] = this->m_Matrix[i][j];
}
}
return tempMat.Det();
}
}
//求解余子式
Matrix& Matrix::ChildMatrix(int row, int col)
{
if (this->m_Row != this->m_Col)
{
std::cout << "Error: <ChildMatrix> Matrix row != col" << std::endl;
return *this;
}
else if (this->m_Row <= 1)
{
std::cout << "Error: <ChildMatrix> Matrix Row Less 1 " << std::endl;
return *this;
}
else if ((row > this->m_Row) || (col > this->m_Col))
{
std::cout << "Error: <ChildMatrix> Input Row Or Col More Than Matix Max Row Or Col" << std::endl;
return* this;
}
else
{
Matrix* resMat = new Matrix(this->m_Row-1, this->m_Col-1);
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < this->m_Col; j++)
{
if ((i < row) && (j < col))
resMat->m_Matrix[i][j] = this->m_Matrix[i][j];
else if((i > row) && (j < col))
resMat->m_Matrix[i-1][j] = this->m_Matrix[i][j];
else if((i < row) && (j > col))
resMat->m_Matrix[i][j - 1] = this->m_Matrix[i][j];
else if((i > row) && (j > col))
resMat->m_Matrix[i - 1][j - 1] = this->m_Matrix[i][j];
}
}
return *resMat;
}
}
//列主消元处理为上三角矩阵
double Matrix::DetRow()
{
//交换标志位 1代表偶数次交换 -1代表奇数次交换
int flagShift = 1;
//本矩阵
Matrix *localMat = new Matrix(*this);
//行列式数值
double resDet = 1.0;
//*******************通过交换 num1*i + num2*j 实现下三角为0***************//
for (int i = 0; i < localMat->m_Row - 1; i++)
{
//记录最大行所在行标
int tempMaxRow = i;
for (int i1 = i + 1; i1 < localMat->m_Row; i1++)
{
if (abs(localMat->m_Matrix[i1][i]) > abs(localMat->m_Matrix[tempMaxRow][i]))
{
tempMaxRow = i1;
}
}
if (tempMaxRow != i)
{
//std::cout << i << " 行交换" << tempMaxRow << " 行" << std::endl;
//进行交换 将当前第i行与第tempMaxRow行进行互换 初等行变换
*localMat = localMat->SwapRow(i, tempMaxRow);
//记录交换次数
flagShift = -flagShift;
//localMat->PrintMat();
}
//此对角线以下的元素通过初等变化为0
for (int i2 = i + 1; i2 < localMat->m_Row; i2++)
{
if (localMat->m_Matrix[i2][i] != 0)
{
//std::cout << "<" << localMat->m_Matrix[i][i] << "> *" << i2 << " 行 + <" << -1.0 * (localMat->m_Matrix[i2][i]) << "> *" << i << " 行" << std::endl;
*localMat = localMat->AddRow(i2, i, -1.0 * (localMat->m_Matrix[i2][i]) / localMat->m_Matrix[i][i]);
//localMat->PrintMat();
}
}
}
//计算行列式数值 对角线相乘
for (int i = 0; i < localMat->m_Row; i++)
{
resDet = resDet * localMat->m_Matrix[i][i];
}
//矩阵交换一次就会变号
resDet = flagShift * resDet;
//清理localMatrix
delete localMat;
return resDet;
}
//矩阵求逆
Matrix& Matrix::Inverse()
{
if (abs(this->DetRow()) < MIN_DET)
{
std::cout << "Error: <Inverse> Matrix Det Near 0" << std::endl;
return *this;
}
else
{
Matrix* resMat = new Matrix(this->m_Row, this->m_Col);
for (int i = 0; i < this->m_Row; i++)
{
for (int j = 0; j < this->m_Col; j++)
{
Matrix* temp = &(this->ChildMatrix(j, i));
resMat->m_Matrix[i][j] = pow(-1.0, (i + j)) / this->DetRow() * (temp->DetRow());
delete temp;
}
}
return *resMat;
}
}
//矩阵求逆 行初等变化
Matrix& Matrix::InverseRow()
{
//错误判断
if (abs(this->DetRow()) < MIN_DET)
{
std::cout << "Error: <InverseRow> Matrix Det Near 0" << std::endl;
return *this;
}
else if (this->m_Row <= 1)
{
std::cout << "Error: <InverseRow> Size Less 2" << std::endl;
return *this;
}
else
{
//单位矩阵 与带转换矩阵维度相同的
Matrix uint = this->Uint();
//结果矩阵 逆矩阵 初始状态与本矩阵相同 为不使本矩阵发生改变
Matrix temp(this->m_Row, this->m_Col);
Matrix* resMat = new Matrix(temp.Uint());
//本矩阵
Matrix localMat(*this);
//*******************通过交换 num1*i + num2*j 实现下三角为0***************//
for (int i = 0; i < localMat.m_Row - 1; i++)
{
//记录最大行所在行标
int tempMaxRow = i;
for (int i1 = i + 1; i1 < localMat.m_Row; i1++)
{
if (abs(localMat.m_Matrix[i1][i]) > abs(localMat.m_Matrix[tempMaxRow][i]))
{
tempMaxRow = i1;
}
}
if (tempMaxRow != i)
{
//std::cout << i << " 行交换" << tempMaxRow << " 行" << std::endl;
//进行交换 将当前第i行与第tempMaxRow行进行互换 初等行变换
localMat = localMat.SwapRow(i, tempMaxRow);
*resMat = resMat->SwapRow(i, tempMaxRow);
//localMat.PrintMat();
}
//此对角线以下的元素通过初等变化为0
for (int i2 = i + 1; i2 < localMat.m_Row; i2++)
{
if (localMat.m_Matrix[i2][i] != 0)
{
//std::cout << "<" << localMat.m_Matrix[i][i] << "> *" << i2 << " 行 + <" << -1.0 * (localMat.m_Matrix[i2][i]) << "> *" << i << " 行" << std::endl;
*resMat = resMat->AddRow(i2, i, -1.0 * (localMat.m_Matrix[i2][i]) / localMat.m_Matrix[i][i]);
localMat = localMat.AddRow(i2, i, -1.0 * (localMat.m_Matrix[i2][i]) / localMat.m_Matrix[i][i]);
//localMat.PrintMat();
}
}
}
//错误判断
if (localMat.m_Matrix[localMat.m_Row - 1][localMat.m_Col - 1] == 0)
{
std::cout << "Error: <InverseRow> marix[" << localMat.m_Row - 1 << "][" << localMat.m_Col - 1 <<"] == 0" << std::endl;
return *this;
}
//*******************通过 num1*i + num2*j 实现上三角为0***************//
for (int i = localMat.m_Row - 1; i > 0; i--)
{
for (int i2 = i - 1; i2 >= 0; i2--)
{
if (localMat.m_Matrix[i2][i] != 0)
{
//std::cout << "<" << localMat.m_Matrix[i][i] << "> *" << i2 << " 行 + <" << -1.0 * (localMat.m_Matrix[i2][i]) << "> *" << i << " 行" << std::endl;
*resMat = resMat->AddRow(i2, i, -1.0 * (localMat.m_Matrix[i2][i]) / localMat.m_Matrix[i][i]);
localMat = localMat.AddRow(i2, i, -1.0 * (localMat.m_Matrix[i2][i]) / localMat.m_Matrix[i][i]);
//localMat.PrintMat();
}
}
}
//*******************通过 i*num 实现矩阵为单位矩阵***************//
for (int i = 0; i < localMat.m_Row; i++)
{
if (localMat.m_Matrix[i][i] == 0)
{
std::cout << "Error: <InverseRow> matrix[" << i << "]" << "[" << i << "] == 0" << std::endl;
return *this;
}
else
{
//std::cout << "<" << 1 / localMat.m_Matrix[i][i] << "> *" << i << " 行" << std::endl;
*resMat = resMat->MultRow(1 / localMat.m_Matrix[i][i], i);
localMat = localMat.MultRow(1 / localMat.m_Matrix[i][i], i);
//localMat.PrintMat();
}
}
return *resMat;
}
}
//矩阵求逆 下三角矩阵
Matrix& Matrix::InverseDownTriangle()
{
//错误判断 方阵检测
if (this->m_Row != this->m_Col)
{
std::cout << "Error: <InverseDownTriangle> Matrix Col != Row" << std::endl;
return *this;
}
//下三角求逆
Matrix* resMat = new Matrix(*this);
for (int i = 0; i < resMat->m_Row; i++)
{
for (int j = 0; j <= i; j++)
{
//分段求解 对角线为倒数
if (i == j)
{
resMat->m_Matrix[i][j] = 1 / resMat->m_Matrix[i][j];
}
else
{
//分段求解 非对角线元素
double tempSum = 0.0;
for (int k = j; k <= i - 1; k++)
{
tempSum += resMat->m_Matrix[i][k] * resMat->m_Matrix[k][j];
}
resMat->m_Matrix[i][j] = -1.0*tempSum / resMat->m_Matrix[i][i];
}
}
}
return *resMat;
}
//矩阵求逆 上三角矩阵
Matrix& Matrix::InverseUpTriangle()
{
//错误判断 方阵检测
if (this->m_Row != this->m_Col)
{
std::cout << "Error: <InverseUpTriangle> Matrix Col != Row" << std::endl;
return *this;
}
//上三角求逆
Matrix* resMat = new Matrix(*this);
for (int j = resMat->m_Col-1; j >=0; j--)
{
for (int i = j; i >=0; i--)
{
//分段求解 对角线为倒数
if (i == j)
{
resMat->m_Matrix[i][j] = 1 / resMat->m_Matrix[i][j];
}
else
{
//分段求解 非对角线元素
double tempSum = 0.0;
for (int k = j; k >= i+1; k--)
{
tempSum += resMat->m_Matrix[i][k] * resMat->m_Matrix[k][j];
}
resMat->m_Matrix[i][j] = -1.0 * tempSum / resMat->m_Matrix[i][i];
}
}
}
return *resMat;
}
//矩阵LU分解 顺序分解 对于病态矩阵可能存在精度问题
void Matrix::ResolveLU(Matrix& LMat, Matrix& UMat)
{
if (this->m_Col != this->m_Row)
{
std::cout << "Error: <ResolveLU> Is Not Square Matrix" << std::endl;
return;
}
//存在性判定 顺序主子式不为0
for (int i = 0; i < this->m_Row; i++)
{
if (this->Det(i) == 0)
{
std::cout << "Error: <ResolveLU> order Det = 0" << std::endl;
return;
}
}
//LU 分解
//L矩阵为单位矩阵
LMat = this->Uint();
//U矩阵初始化为空矩阵
Matrix temp(this->m_Row, this->m_Col);
UMat = temp;
for (int i = 0; i < this->m_Row; i++)
{
//计算U
for (int j1 = i; j1 < this->m_Col; j1++)
{
double tempSum1 = 0.0;
if (i != 0)
{
for (int j2 = 0; j2 <= i - 1; j2++)
{
tempSum1 += LMat.m_Matrix[i][j2] * UMat.m_Matrix[j2][j1];
}
}
UMat.m_Matrix[i][j1] = this->m_Matrix[i][j1] - tempSum1;
}
//计算L
for (int i1 = i; i1 < this->m_Row; i1++)
{
double tempSum2 = 0.0;
if (i != 0)
{
for (int j2 = 0; j2 <= i - 1; j2++)
{
tempSum2 += LMat.m_Matrix[i1][j2] * UMat.m_Matrix[j2][i];
}
}
LMat.m_Matrix[i1][i] = (this->m_Matrix[i1][i] - tempSum2)/UMat.m_Matrix[i][i];
}
}
}
//矩阵的LUP分解 P*A = L*U 添加了列主消元功能
//L为主对角线元素为1的下三角矩阵 U为上二角矩阵 P为行交换矩阵 P*A=L*U
void Matrix::ResolveLUP(Matrix& LMat, Matrix& UMat, Matrix& PMat)
{
//条件判断 矩阵行列式不为0
if (this->Det() == 0)
{
std::cout << "Error: <ResolveLUP> Can't Resolve Matrix To L U P" << std::endl;
return;
}
//初始化 L U P
LMat = this->Uint();
PMat = this->Uint();
UMat = *this;
//进行分解计算
for (int i = 0; i < UMat.m_Row - 1; i++)
{
//记录最大行所在行标
int tempMaxRow = i;
for (int i1 = i + 1; i1 < UMat.m_Row; i1++)
{
if (abs(UMat.m_Matrix[i1][i]) > abs(UMat.m_Matrix[tempMaxRow][i]))
{
tempMaxRow = i1;
}
}
//进行交换 将当前第i行与第tempMaxRow行进行互换 初等行变换
UMat = UMat.SwapRow(i, tempMaxRow);
//L矩阵做出对应交换 先交换<itempMaxRow>列再交换<itempMaxRow>行
LMat = LMat.SwapCol(i, tempMaxRow);
LMat = LMat.SwapRow(i, tempMaxRow);
//P矩阵做出对应变换 交换<itempMaxRow>行
PMat = PMat.SwapRow(i, tempMaxRow);
//高斯消元 V矩阵消除下三角区域，L矩阵添加下三角区域
for (int i1 = i + 1; i1 < UMat.m_Row; i1++)
{
//记录消元系数
double deleteVar = UMat.m_Matrix[i1][i] / UMat.m_Matrix[i][i];
//L矩阵列填充
LMat.m_Matrix[i1][i] = deleteVar;
//U矩阵列消除
UMat = UMat.MultRow(UMat.m_Matrix[i][i], i1).AddRow(i1, i, -1.0 * UMat.m_Matrix[i1][i]).MultRow(1 / UMat.m_Matrix[i][i], i1);
}
}
return;
}
- App下载
- 项目客服
- 培训客服
- 平台客服
TOP
13
4
7
