#ifndef _MATRIX_H
#define _MATRIX_H
class Matrix
{
private :
int row; // 矩阵的行数
int col; // 矩阵的列数
int n; // 矩阵元素个数
double * mtx; // 动态分配用来存放数组的空间
public :
Matrix( int row = 1 , int col = 1 ); // 带默认参数的构造函数
Matrix( int row, int col, double mtx[]); // 用数组创建一个矩阵
Matrix( const Matrix & obj); // copy构造函数
~ Matrix() { delete[] this -> mtx; }
void print() const ; // 格式化输出矩阵
int getRow() const { return this -> row; } // 访问矩阵行数
int getCol() const { return this -> col; } // 访问矩阵列数
int getN() const { return this -> n; } // 访问矩阵元素个数
double * getMtx() const { return this -> mtx; } // 获取该矩阵的数组
// 用下标访问矩阵元素
double get ( const int i, const int j) const ;
// 用下标修改矩阵元素值
void set ( const int i, const int j, const double e);
// 重载了一些常用操作符,包括 +,-,x,=,负号,正号,
// A = B
Matrix & operator = ( const Matrix & obj);
// +A
Matrix operator + () const { return * this ; }
// -A
Matrix operator - () const ;
// A + B
friend Matrix operator + ( const Matrix & A, const Matrix & B);
// A - B
friend Matrix operator - ( const Matrix & A, const Matrix & B);
// A * B 两矩阵相乘
friend Matrix operator * ( const Matrix & A, const Matrix & B);
// a * B 实数与矩阵相乘
friend Matrix operator * ( const double & a, const Matrix & B);
// A 的转置
friend Matrix trv( const Matrix & A);
// A 的行列式值,采用列主元消去法
// 求行列式须将矩阵化为三角阵,此处为了防止修改原矩阵,采用传值调用
friend double det(Matrix A);
// A 的逆矩阵,采用高斯-若当列主元消去法
friend Matrix inv(Matrix A);
};
#endif
#define _MATRIX_H
class Matrix
{
private :
int row; // 矩阵的行数
int col; // 矩阵的列数
int n; // 矩阵元素个数
double * mtx; // 动态分配用来存放数组的空间
public :
Matrix( int row = 1 , int col = 1 ); // 带默认参数的构造函数
Matrix( int row, int col, double mtx[]); // 用数组创建一个矩阵
Matrix( const Matrix & obj); // copy构造函数
~ Matrix() { delete[] this -> mtx; }
void print() const ; // 格式化输出矩阵
int getRow() const { return this -> row; } // 访问矩阵行数
int getCol() const { return this -> col; } // 访问矩阵列数
int getN() const { return this -> n; } // 访问矩阵元素个数
double * getMtx() const { return this -> mtx; } // 获取该矩阵的数组
// 用下标访问矩阵元素
double get ( const int i, const int j) const ;
// 用下标修改矩阵元素值
void set ( const int i, const int j, const double e);
// 重载了一些常用操作符,包括 +,-,x,=,负号,正号,
// A = B
Matrix & operator = ( const Matrix & obj);
// +A
Matrix operator + () const { return * this ; }
// -A
Matrix operator - () const ;
// A + B
friend Matrix operator + ( const Matrix & A, const Matrix & B);
// A - B
friend Matrix operator - ( const Matrix & A, const Matrix & B);
// A * B 两矩阵相乘
friend Matrix operator * ( const Matrix & A, const Matrix & B);
// a * B 实数与矩阵相乘
friend Matrix operator * ( const double & a, const Matrix & B);
// A 的转置
friend Matrix trv( const Matrix & A);
// A 的行列式值,采用列主元消去法
// 求行列式须将矩阵化为三角阵,此处为了防止修改原矩阵,采用传值调用
friend double det(Matrix A);
// A 的逆矩阵,采用高斯-若当列主元消去法
friend Matrix inv(Matrix A);
};
#endif
实现
#include < iostream.h >
#include < math.h >
#include < stdlib.h >
#include < iomanip.h >
#include " matrix.h "
// 带默认参数值的构造函数
// 构造一个row行col列的零矩阵
Matrix::Matrix( int row, int col) {
this->row = row;
this->col = col;
this->n = row * col;
this->mtx = new double[n];
for(int i=0; i<n; i++)
this->mtx[i] = 0.0;
}
// 用一个数组初始化矩阵
Matrix::Matrix( int row, int col, double mtx[])
{
this->row = row;
this->col = col;
this->n = row * col;
this->mtx = new double[n];
for(int i=0; i<n; i++)
this->mtx[i] = mtx[i];
}
// 拷贝构造函数,因为成员变量含有动态空间,防止传递参数
// 等操作发生错误
Matrix::Matrix( const Matrix & obj) {
this->row = obj.getRow();
this->col = obj.getCol();
this->n = obj.getN();
this->mtx = new double[n];
for(int i=0; i<n; i++)
this->mtx[i] = obj.getMtx()[i];
}
// 格式化输出矩阵所有元素
void Matrix::print() const {
for(int i=0; i<this->row; i++){
for(int j=0; j<this->col; j++)
if(fabs(this->get(i,j)) <= 1.0e-10)
cout << setiosflags(ios::left) << setw(12) << 0.0 << ' ';
else
cout << setiosflags(ios::left)
<< setw(12) << this->get(i,j) << ' ';
cout <<endl;
}
}
// 获取矩阵元素
// 注意这里矩阵下标从(0,0)开始
double Matrix:: get ( const int i, const int j) const {
return this->mtx[i*this->col + j];
}
// 修改矩阵元素
void Matrix:: set ( const int i, const int j, const double e) {
this->mtx[i*this->col + j] = e;
}
// 重载赋值操作符,由于成员变量中含有动态分配
Matrix & Matrix:: operator = ( const Matrix & obj) {
if(this == &obj) // 将一个矩阵赋给它自己时简单做返回即可
return *this;
delete[] this->mtx; // 首先删除目的操作数的动态空间
this->row = obj.getRow();
this->col = obj.getCol();
this->n = obj.getN();
this->mtx = new double[n]; // 重新分配空间保存obj的数组
for(int i=0; i<n; i++)
this->mtx[i] = obj.getMtx()[i];
return *this;
}
// 负号操作符,返回值为该矩阵的负矩阵,原矩阵不变
Matrix Matrix:: operator - () const {
// 为了不改变原来的矩阵,此处从新构造一个矩阵
Matrix _A(this->row, this->col);
for(int i=0; i<_A.n; i++)
_A.mtx[i] = -(this->mtx[i]);
return _A;
}
// 矩阵求和,对应位置元素相加
Matrix operator + ( const Matrix & A, const Matrix & B) {
Matrix AB(A.row, A.col);
if(A.row!=B.row || A.col!=B.col){
cout << "Can't do A+B\n"; // 如果矩阵A和B行列数不一致则不可相加
exit(0);
}
for(int i=0; i<AB.n; i++)
AB.mtx[i] = A.mtx[i] + B.mtx[i];
return AB;
}
// 矩阵减法,用加上一个负矩阵来实现
Matrix operator - ( const Matrix & A, const Matrix & B) {
return (A + (-B));
}
// 矩阵乘法
Matrix operator * ( const Matrix & A, const Matrix & B) {
if(A.col != B.row){ // A的列数必须和B的行数一致
cout << "Can't multiply\n";
exit(0);
}
Matrix AB(A.row, B.col); // AB用于保存乘积
for(int i=0; i<AB.row; i++)
for(int j=0; j<AB.col; j++)
for(int k=0; k<A.col; k++)
AB.set(i, j, AB.get(i,j) + A.get(i,k)*B.get(k,j));
return AB;
}
// 矩阵与实数相乘
Matrix operator * ( const double & a, const Matrix & B) {
Matrix aB(B);
for(int i=0; i<aB.row; i++)
for(int j=0; j<aB.col; j++)
aB.set(i,j, a*B.get(i,j));
return aB;
}
// 矩阵的转置 将(i,j)与(j,i)互换
// 此函数返回一个矩阵的转置矩阵,并不改变原来的矩阵
Matrix trv( const Matrix & A) {
Matrix AT(A.col, A.row);
for(int i=0; i<AT.row; i++)
for(int j=0; j<AT.col; j++)
AT.set(i, j, A.get(j,i));
return AT;
}
// 矩阵行列式值,采用列主元消去法
double det(Matrix A) {
if(A.row != A.col) { // 矩阵必须为n*n的才可进行行列式求值
cout << "error" << endl;
return 0.0; //如果不满足行列数相等返回0.0 }
double detValue = 1.0; // 用于保存行列式值
for(int i=0; i<A.getRow()-1; i++){ // 需要n-1步列化零操作
//------------------ 选主元 ---------------------------------
double max = fabs(A.get(i,i)); // 主元初始默认为右下方矩阵首个元素
int ind = i;
// 主元行号默认为右下方矩阵首行
for(int j=i+1; j<A.getRow(); j++){ // 选择列主元
if(fabs(A.get(j,i)) > max){ // 遇到绝对值更大的元素
max = fabs(A.get(j,i)); // 更新主元值
ind = j; // 更新主元行号
}
}//loop j
//------------------- 移动主元行 -----------------------------
if(max <= 1.0e-10) return 0.0; // 右下方矩阵首行为零,显然行列式值为零
if(ind != i){// 主元行非右下方矩阵首行
for(int k=i; k<A.getRow(); k++){ // 将主元行与右下方矩阵首行互换
double temp = A.get(i,k);
A.set(i,k,A.get(ind,k));
A.set(ind,k,temp);
}
detValue = -detValue; // 互换行列式两行,行列式值反号
}
//------------------- 消元 ----------------------------------
for(j=i+1; j<A.getRow(); j++){ // 遍历行
double temp = A.get(j,i)/A.get(i,i);
for(int k=i; k<A.getRow(); k++) // 遍历行中每个元素,行首置0
A.set(j, k, A.get(j,k)-A.get(i,k)*temp); }
detValue *= A.get(i,i); // 每步消元都会产生一个对角线上元素,将其累乘
}// loop i
// 注意矩阵最后一个元素在消元的过程中没有被累乘到
return detValue * A.get(A.getRow()-1, A.getRow()-1);}//det()
// A的逆矩阵 高斯-若当消去法,按列选主元
Matrix inv(Matrix A){
if(A.row != A.col){ // 只可求狭义逆矩阵,即行列数相同
cout << "Matrix should be N x N\n";
exit(0);
}
// 构造一个与A行列相同的单位阵B
Matrix B(A.row,A.col);
for(int r=0; r<A.row; r++)
for(int c=0; c<A.col; c++)
if(r == c) B.set(r,c,1.0);
// 对矩阵A进行A.row次消元运算,每次保证第K列只有对角线上非零
// 同时以同样的操作施与矩阵B,结果A变为单位阵B为所求逆阵
for(int k=0; k<A.row; k++){
//------------------ 选主元 --------------------------------------
double max = fabs(A.get(k,k)); // 主元初始默认为右下方矩阵首个元素
int ind = k; // 主元行号默认为右下方矩阵首行
// 结果第ind行为列主元行
for(int n=k+1; n<A.getRow(); n++){
if(fabs(A.get(n,k)) > max){// 遇到绝对值更大的元素
max = fabs(A.get(n,k)); // 更新主元值
ind = n;// 更新主元行号
}
}
//------------------- 移动主元行 --------------------------------
if(ind != k){// 主元行不是右下方矩阵首行
for(int m=k; m<A.row; m++){// 将主元行与右下方矩阵首行互换
double tempa = A.get(k,m);
A.set(k, m, A.get(ind,m));
A.set(ind, m, tempa);
}
for(m=0; m<B.row; m++){
double tempb = B.get(k,m); // 对矩阵B施以相同操作
B.set(k, m, B.get(ind,m)); // B与A阶数相同,可在一个循环中
B.set(ind, m, tempb);
}
}
//--------------------- 消元 -----------------------------------
// 第k次消元操作,以第k行作为主元行,将其上下各行的第k列元素化为零
// 同时以同样的参数对B施以同样的操作,此时可以将B看作A矩阵的一部分
for(int i=0; i<A.col; i++){
if(i != k){
double Mik = -A.get(i,k)/A.get(k,k);
for(int j=k+1; j<A.row; j++)
A.set(i, j, A.get(i,j) + Mik*A.get(k,j));
for(j=0; j<B.row; j++)
B.set(i, j, B.get(i,j) + Mik*B.get(k,j));
}//end if
}//loop i
double Mkk = 1.0/A.get(k,k);
for(int j=0; j<A.row; j++)
A.set(k, j, A.get(k,j) * Mkk);
for(j=0; j<B.row; j++)
B.set(k, j, B.get(k,j) * Mkk);
}//loop k
return B;
}//inv()
#include < iostream.h >
#include < math.h >
#include < stdlib.h >
#include < iomanip.h >
#include " matrix.h "
// 带默认参数值的构造函数
// 构造一个row行col列的零矩阵
Matrix::Matrix( int row, int col) {
this->row = row;
this->col = col;
this->n = row * col;
this->mtx = new double[n];
for(int i=0; i<n; i++)
this->mtx[i] = 0.0;
}
// 用一个数组初始化矩阵
Matrix::Matrix( int row, int col, double mtx[])
{
this->row = row;
this->col = col;
this->n = row * col;
this->mtx = new double[n];
for(int i=0; i<n; i++)
this->mtx[i] = mtx[i];
}
// 拷贝构造函数,因为成员变量含有动态空间,防止传递参数
// 等操作发生错误
Matrix::Matrix( const Matrix & obj) {
this->row = obj.getRow();
this->col = obj.getCol();
this->n = obj.getN();
this->mtx = new double[n];
for(int i=0; i<n; i++)
this->mtx[i] = obj.getMtx()[i];
}
// 格式化输出矩阵所有元素
void Matrix::print() const {
for(int i=0; i<this->row; i++){
for(int j=0; j<this->col; j++)
if(fabs(this->get(i,j)) <= 1.0e-10)
cout << setiosflags(ios::left) << setw(12) << 0.0 << ' ';
else
cout << setiosflags(ios::left)
<< setw(12) << this->get(i,j) << ' ';
cout <<endl;
}
}
// 获取矩阵元素
// 注意这里矩阵下标从(0,0)开始
double Matrix:: get ( const int i, const int j) const {
return this->mtx[i*this->col + j];
}
// 修改矩阵元素
void Matrix:: set ( const int i, const int j, const double e) {
this->mtx[i*this->col + j] = e;
}
// 重载赋值操作符,由于成员变量中含有动态分配
Matrix & Matrix:: operator = ( const Matrix & obj) {
if(this == &obj) // 将一个矩阵赋给它自己时简单做返回即可
return *this;
delete[] this->mtx; // 首先删除目的操作数的动态空间
this->row = obj.getRow();
this->col = obj.getCol();
this->n = obj.getN();
this->mtx = new double[n]; // 重新分配空间保存obj的数组
for(int i=0; i<n; i++)
this->mtx[i] = obj.getMtx()[i];
return *this;
}
// 负号操作符,返回值为该矩阵的负矩阵,原矩阵不变
Matrix Matrix:: operator - () const {
// 为了不改变原来的矩阵,此处从新构造一个矩阵
Matrix _A(this->row, this->col);
for(int i=0; i<_A.n; i++)
_A.mtx[i] = -(this->mtx[i]);
return _A;
}
// 矩阵求和,对应位置元素相加
Matrix operator + ( const Matrix & A, const Matrix & B) {
Matrix AB(A.row, A.col);
if(A.row!=B.row || A.col!=B.col){
cout << "Can't do A+B\n"; // 如果矩阵A和B行列数不一致则不可相加
exit(0);
}
for(int i=0; i<AB.n; i++)
AB.mtx[i] = A.mtx[i] + B.mtx[i];
return AB;
}
// 矩阵减法,用加上一个负矩阵来实现
Matrix operator - ( const Matrix & A, const Matrix & B) {
return (A + (-B));
}
// 矩阵乘法
Matrix operator * ( const Matrix & A, const Matrix & B) {
if(A.col != B.row){ // A的列数必须和B的行数一致
cout << "Can't multiply\n";
exit(0);
}
Matrix AB(A.row, B.col); // AB用于保存乘积
for(int i=0; i<AB.row; i++)
for(int j=0; j<AB.col; j++)
for(int k=0; k<A.col; k++)
AB.set(i, j, AB.get(i,j) + A.get(i,k)*B.get(k,j));
return AB;
}
// 矩阵与实数相乘
Matrix operator * ( const double & a, const Matrix & B) {
Matrix aB(B);
for(int i=0; i<aB.row; i++)
for(int j=0; j<aB.col; j++)
aB.set(i,j, a*B.get(i,j));
return aB;
}
// 矩阵的转置 将(i,j)与(j,i)互换
// 此函数返回一个矩阵的转置矩阵,并不改变原来的矩阵
Matrix trv( const Matrix & A) {
Matrix AT(A.col, A.row);
for(int i=0; i<AT.row; i++)
for(int j=0; j<AT.col; j++)
AT.set(i, j, A.get(j,i));
return AT;
}
// 矩阵行列式值,采用列主元消去法
double det(Matrix A) {
if(A.row != A.col) { // 矩阵必须为n*n的才可进行行列式求值
cout << "error" << endl;
return 0.0; //如果不满足行列数相等返回0.0 }
double detValue = 1.0; // 用于保存行列式值
for(int i=0; i<A.getRow()-1; i++){ // 需要n-1步列化零操作
//------------------ 选主元 ---------------------------------
double max = fabs(A.get(i,i)); // 主元初始默认为右下方矩阵首个元素
int ind = i;
// 主元行号默认为右下方矩阵首行
for(int j=i+1; j<A.getRow(); j++){ // 选择列主元
if(fabs(A.get(j,i)) > max){ // 遇到绝对值更大的元素
max = fabs(A.get(j,i)); // 更新主元值
ind = j; // 更新主元行号
}
}//loop j
//------------------- 移动主元行 -----------------------------
if(max <= 1.0e-10) return 0.0; // 右下方矩阵首行为零,显然行列式值为零
if(ind != i){// 主元行非右下方矩阵首行
for(int k=i; k<A.getRow(); k++){ // 将主元行与右下方矩阵首行互换
double temp = A.get(i,k);
A.set(i,k,A.get(ind,k));
A.set(ind,k,temp);
}
detValue = -detValue; // 互换行列式两行,行列式值反号
}
//------------------- 消元 ----------------------------------
for(j=i+1; j<A.getRow(); j++){ // 遍历行
double temp = A.get(j,i)/A.get(i,i);
for(int k=i; k<A.getRow(); k++) // 遍历行中每个元素,行首置0
A.set(j, k, A.get(j,k)-A.get(i,k)*temp); }
detValue *= A.get(i,i); // 每步消元都会产生一个对角线上元素,将其累乘
}// loop i
// 注意矩阵最后一个元素在消元的过程中没有被累乘到
return detValue * A.get(A.getRow()-1, A.getRow()-1);}//det()
// A的逆矩阵 高斯-若当消去法,按列选主元
Matrix inv(Matrix A){
if(A.row != A.col){ // 只可求狭义逆矩阵,即行列数相同
cout << "Matrix should be N x N\n";
exit(0);
}
// 构造一个与A行列相同的单位阵B
Matrix B(A.row,A.col);
for(int r=0; r<A.row; r++)
for(int c=0; c<A.col; c++)
if(r == c) B.set(r,c,1.0);
// 对矩阵A进行A.row次消元运算,每次保证第K列只有对角线上非零
// 同时以同样的操作施与矩阵B,结果A变为单位阵B为所求逆阵
for(int k=0; k<A.row; k++){
//------------------ 选主元 --------------------------------------
double max = fabs(A.get(k,k)); // 主元初始默认为右下方矩阵首个元素
int ind = k; // 主元行号默认为右下方矩阵首行
// 结果第ind行为列主元行
for(int n=k+1; n<A.getRow(); n++){
if(fabs(A.get(n,k)) > max){// 遇到绝对值更大的元素
max = fabs(A.get(n,k)); // 更新主元值
ind = n;// 更新主元行号
}
}
//------------------- 移动主元行 --------------------------------
if(ind != k){// 主元行不是右下方矩阵首行
for(int m=k; m<A.row; m++){// 将主元行与右下方矩阵首行互换
double tempa = A.get(k,m);
A.set(k, m, A.get(ind,m));
A.set(ind, m, tempa);
}
for(m=0; m<B.row; m++){
double tempb = B.get(k,m); // 对矩阵B施以相同操作
B.set(k, m, B.get(ind,m)); // B与A阶数相同,可在一个循环中
B.set(ind, m, tempb);
}
}
//--------------------- 消元 -----------------------------------
// 第k次消元操作,以第k行作为主元行,将其上下各行的第k列元素化为零
// 同时以同样的参数对B施以同样的操作,此时可以将B看作A矩阵的一部分
for(int i=0; i<A.col; i++){
if(i != k){
double Mik = -A.get(i,k)/A.get(k,k);
for(int j=k+1; j<A.row; j++)
A.set(i, j, A.get(i,j) + Mik*A.get(k,j));
for(j=0; j<B.row; j++)
B.set(i, j, B.get(i,j) + Mik*B.get(k,j));
}//end if
}//loop i
double Mkk = 1.0/A.get(k,k);
for(int j=0; j<A.row; j++)
A.set(k, j, A.get(k,j) * Mkk);
for(j=0; j<B.row; j++)
B.set(k, j, B.get(k,j) * Mkk);
}//loop k
return B;
}//inv()