matrix.h

偏微分方程数值解- 这是矩阵的紧致存储及在紧致存储的矩阵上用共轭梯度(CG)方法解Poisson方程的实现, gcc mingw 3.4.4下通过

#include<iostream.h>

class Matrix {
public:
  Matrix(unsigned rows, unsigned cols);
  double& operator() (unsigned row, unsigned col);
  double  operator() (unsigned row, unsigned col) const;
  // ...
 ~Matrix();
                             // 析构函数
//  Matrix(Matrix& m){
//  };               // 拷贝构造函数
  

  Matrix& operator= (Matrix& m){
        unsigned a = m.getrow(), b = m.getcol();
        Matrix::Matrix(a, b);
    	if(data_==NULL){
			cout<<"Not enough memory for "<<a<<"*"<<b<<" matrix.";
			exit(-1);
		}
        for(unsigned i=1;i<=a;i++)
            for(unsigned j=1;j<=b;j++)
                data_[b*(i-1) + (j-1)] = m(i,j);

  };   // 赋值运算符
  const unsigned getrow() {return rows_; };
  const unsigned getcol() {return cols_; };
  void clear(){	unsigned d = rows_ * cols_;	for(unsigned i=0;i<d;i++) data_[i]=0.0;}
  bool rescale(unsigned rows, unsigned cols);

  // ...
protected:
  unsigned rows_, cols_;
  double* data_;
};

inline
Matrix::Matrix(unsigned rows, unsigned cols)
  : rows_ (rows),
    cols_ (cols),
    data_ (new double[rows * cols])
{
  if (rows == 0 || cols == 0)
  void();
	if(data_==NULL){
		cout<<"Not enough memory for "<<rows<<"*"<<cols<<" matrix.";
		exit(-1);
	}

//     throw BadIndex("Matrix constructor has 0 size");
	clear();
}

inline
bool Matrix::rescale(unsigned rows, unsigned cols)
{
	if((rows_==rows) && (cols_==cols)) return false;
	rows_ = rows;
    cols_ = cols;
	delete[] data_;
	data_ = new double[rows * cols];
	if(data_==NULL){
		cout<<"Not enough memory for "<<rows<<"*"<<cols<<" matrix.";
		exit(-1);
	}
	clear();
	return true;
}
inline
Matrix::~Matrix()
{
  delete[] data_;
}

inline
double& Matrix::operator() (unsigned row, unsigned col)
{
  if (row > rows_ || col > cols_)    void();
//     throw BadIndex("Matrix subscript out of bounds");
  return data_[cols_*(row-1) + (col-1)];
}

inline
double Matrix::operator() (unsigned row, unsigned col) const
{
  if (row > rows_ || col > cols_)    void();
//     throw BadIndex("const Matrix subscript out of bounds");
  return data_[cols_*(row-1) + (col-1)];
} 

/*
int main()
{
  Matrix m(10,10);
  m(5,8) = 106.15;
  cout << m(5,8);
   getchar();
} 
*/