/***************************************************************
 * PrMatrix.h
 *
 * This class provides a rowxcol matrix.
 *
 ****************************************************************/

/*
 * modification history
 *----------------------
 *
 * 11/12/97: K.C. Chang: created.
 */

#ifndef _PrMatrix_h
#define _PrMatrix_h

#include "PrGlobalDefn.h"
#include "PrVector.h"

class PrMatrix
{
  public:
    PrMatrix(int row,int col); // initialized to [0]
    PrMatrix(const PrMatrix &m);
    ~PrMatrix();

    PrMatrix &operator=(const PrMatrix &m);
    Float *operator[](int row);
    const Float *operator[](int row) const;

    Float &elementAt(int i,int j);
    const Float &elementAt(int i,int j) const;

    int row() const { return row_; }
    int column() const { return col_; }

    void zero();

    PrMatrix operator-() const;
    void negate(PrMatrix &dest) const;

    PrMatrix operator~() const;
    void inverse(PrMatrix &dest) const;
    //(symmetric and positive definite)
    void inverseSPD(PrMatrix &dest) const;

    // y = (*this) x
    // solve for x given y
    void solve(const PrVector &y,PrVector &x) const;
    //(symmetric and positive definite)
    void solveSPD(const PrVector &y,PrVector &x) const;

    PrMatrix operator+(const PrMatrix &m) const;
    void add(const PrMatrix &m,PrMatrix &dest) const;
    PrMatrix operator-(const PrMatrix &m) const;
    void subtract(const PrMatrix &m,PrMatrix &dest) const;
    PrMatrix operator*(const PrMatrix &m) const;
    void multiply(const PrMatrix &m,PrMatrix &dest) const;
    PrVector operator*(const PrVector &v) const;
    void multiply(const PrVector &v,PrVector &dest) const;
    PrMatrix operator*(Float s) const;
    void multiply(Float s,PrMatrix &dest) const;

    PrMatrix &operator+=(const PrMatrix &m);
    PrMatrix &operator-=(const PrMatrix &m);
    PrMatrix &operator*=(const PrMatrix &m);
    PrMatrix &operator*=(Float s);

    PrMatrix transpose() const;
    void transpose(PrMatrix &dest) const;

    void diagonal(Float *data); // set the diagonal
    PrVector diagonal() const;
    void diagonal(PrVector &dest) const;

    void display(char *name=NULL);

  private:
    int row_;
    int col_;
    int n_;   // =row_*col_
    Float *data_;

    void LUdecomp(PrMatrix &lu) const;
    //(symmetric and positive definite)
    void LUdecompSPD(PrMatrix &lu) const;

    // y = LU x
    // find x given LU and y 
    void BackSub(const PrMatrix &lu,const PrVector &y,PrVector &x) const;
    //(symmetric and positive definite)
    void BackSubSPD(const PrMatrix &lu,const PrVector &y,PrVector &x) const;

    void AssertDimension(const PrMatrix &m) const;
    void AssertDimension(const PrVector &v) const;
};

#endif // _PrMatrix_h
