qmatvec.h

来自「算断裂的」· C头文件 代码 · 共 178 行

// ------------------------------------------------------------------
// qmatvec.h
//
// This file contains the definitions for vectors and matrices, including
// several linear algebra ops. .
// ------------------------------------------------------------------
// Author: Stephen A. Vavasis
// Copyright (c) 1999 by Cornell University.  All rights reserved.
// 
// See the accompanying file 'Copyright' for authorship information,
// the terms of the license governing this software, and disclaimers
// concerning this software.
// ------------------------------------------------------------------
// This file is part of the QMG software.  
// Version 2.0 of QMG, release date September 3, 1999.
// ------------------------------------------------------------------



#include "qnamesp.h"

#ifndef QMGMATVEC_H
#define QMGMATVEC_H



namespace QMG {
  template<class T, int sz>
    class StatVector {
    private:
      T data_[sz];
    public:
      const T& operator[](int i) const {
#ifdef RANGECHECK
	if (i < 0 || i >= sz) {
	  throw_error("Range error in StatVector 1");
	}
#endif
	return data_[i];
      }
      
      T& operator[](int i) {
#ifdef RANGECHECK
        if (i < 0 || i >= sz) {
          throw_error("Range error in StatVector 2");
        }
#endif
        return data_[i];
      }
    };
}


class QMG::Matrix {
private:
  Real* entries_;
  int numrows_;
  int numcols_;
  bool heap_allocated_;
  // No assignment.
  void operator=(const Matrix&) { }
public:
  Matrix(int m, int n) : entries_(new Real[m*n]), numrows_(m), numcols_(n), 
    heap_allocated_(true) { }
  Matrix(int m, int n, Real* data) : entries_(data), numrows_(m), numcols_(n),
    heap_allocated_(false) { }
  ~Matrix() {if (heap_allocated_) delete entries_;}
  const Real& operator()(int i, int j) const {
#ifdef RANGECHECK
    if (i < 0 || i >= numrows_ || j < 0 || j >= numcols_) {
      throw_error("Range error in Matrix 1");
    }
#endif
    return entries_[i*numcols_+j];
  }
  Real& operator()(int i, int j) {
#ifdef RANGECHECK
    if (i < 0 || i >= numrows_ || j < 0 || j >= numcols_) {
      throw_error("Range error in Matrix 2");
    }
#endif
    return entries_[i*numcols_+j];
  }
  Real linear_solve(Real rhs[], Real soln[]);
  Real determinant() const;
  static Real compute_hh_transform(Real* hhtrans, 
    const Real* vec, 
    int vec_stride, 
    int n);

  // Matrix copying steals a pointer (like autoptr);
  // Destructive copy semantics

private:
  void init_from_mtx_(Matrix& m);
  class Matrix_Returnval {
  private:
    friend class Matrix;
    Matrix& sref_;
    Matrix_Returnval(Matrix& s, int) :
      sref_(s) { }
  };
public:
  Matrix(Matrix& m) {init_from_mtx_(m);}
  Matrix(const Matrix_Returnval& r) {init_from_mtx_(r.sref_);}
  operator Matrix_Returnval() {
    return Matrix_Returnval(*this, 0);
  }

  int numrows() const { return numrows_;}
  int numcols() const { return numcols_;}
};


extern QMG::ostream& operator<<(QMG::ostream& os, const QMG::Matrix& mt);



class QMG::IntMatrix {
private:
  int* entries_;
  int numrows_;
  int numcols_;
  bool heap_allocated_;
  // No assignment.
  void operator=(const IntMatrix&) { }
public:
  IntMatrix(int m, int n) : entries_(new int[m*n]), numrows_(m), numcols_(n), 
    heap_allocated_(true) { }
  ~IntMatrix() {if (heap_allocated_) delete entries_;}
  const int& operator()(int i, int j) const {
#ifdef RANGECHECK
    if (i < 0 || i >= numrows_ || j < 0 || j >= numcols_) {
      throw_error("Range error in IntMatrix 1");
    }
#endif
    return entries_[i*numcols_+j];
  }
  int& operator()(int i, int j) {
#ifdef RANGECHECK
    if (i < 0 || i >= numrows_ || j < 0 || j >= numcols_) {
      throw_error("Range error in IntMatrix 2");
    }
#endif
    return entries_[i*numcols_+j];
  }
  
  // Matrix copying steals a pointer (like autoptr);
  
  // Destructive copy semantics

private:
  void init_from_mtx_(IntMatrix& m);
  class IntMatrix_Returnval {
  private:
    friend class IntMatrix;
    IntMatrix& sref_;
    IntMatrix_Returnval(IntMatrix& s, int) :
      sref_(s) { }
  };
public:
  IntMatrix(IntMatrix& m) {init_from_mtx_(m);}
  IntMatrix(const IntMatrix_Returnval& r) {init_from_mtx_(r.sref_);}
  operator IntMatrix_Returnval() {
    return IntMatrix_Returnval(*this, 0);
  }
  
  int numrows() const { return numrows_;}
  int numcols() const { return numcols_;}
};





#endif