LAP_DoubleMatrix.h

Go to the documentation of this file.

#ifndef LAP_DoubleMatrix_H
#define LAP_DoubleMatrix_H

#include "LAP_DoubleMatrixStorage.h"
#include "LAP_DoubleVector.h"
#include "LAP_ConstDoubleVector.h"


DEFINE_SPTR(LAP_DoubleMatrix);

class LAP_DoubleMatrix: public  LAP_Matrix {
    
    SP_OBJECT(LAP_DoubleMatrix);
    
private:
    // ATTRIBUTES
 
        
   
  
    // ASSOCIATIONS

    SP::LAP_DoubleMatrixStorage mStorage;
    
    // METHOD
  
  
public:
    // CONSTRUCTORS
    LAP_DoubleMatrix();
    LAP_DoubleMatrix(const tLVectorIndex& n,const tLVectorIndex& p);
    
   
  
    // DESTRUCTORS
    virtual ~LAP_DoubleMatrix();
    
 

public:
    inline const double& operator[](const tLVectorIndex& index)  const {
        return getStorage()[index];
    }
    inline double& operator[](const tLVectorIndex& index)  {
        return getStorage()[index];
    }
    
  
    virtual const double& operator()(const tLVectorIndex& i,
                                     const tLVectorIndex& j)  const=0;
    
    virtual double& operator()(const tLVectorIndex&  i,
                               const tLVectorIndex& j)=0;

    // NEW
    // ---

    virtual SP::LAP_DoubleMatrix NewInstance() const=0;
    
    // ----------------
    // MATRIX SETTINGS
    // ----------------
public:
    virtual tBoolean copy(const LAP_DoubleMatrix& x);
    inline tBoolean copy(SPC::LAP_DoubleMatrix x) {
        if (x.get()!=null) return copy(*x.get());
        return true;
    };
    inline tBoolean copy(LAP_DoubleMatrix* x) {
        if (x!=null) return copy(*x);
        return true;
    };

   
public:

    inline void init(const double v) {
        getStorage().init(v);
    }
    
    inline void set(const tLVectorIndex& i,
                     const tLVectorIndex& j,
                     const double& v) {
        (*this)(i,j)=v;
    }
    
    inline  void add(const tLVectorIndex& i,
                     const tLVectorIndex& j,
                     const double& v) {
        (*this)(i,j)+=v;
    }
    
    inline void setStorage(SP::LAP_DoubleMatrixStorage s) {
        mStorage=s;
    }

    virtual void setValues(SP::LAP_DoubleVector v) {
        mStorage->setValuesPointer(v);
        
    }
    virtual void setValues(const tLVectorIndex& n,const double* values) {
        mStorage->setValues(n,values);
        
    }
    
    // GET Methods
    // ===========
public:
   
   

    inline double get(const tLVectorIndex& i,
                      const tLVectorIndex& j) const {
        return (*this)(i,j);
    }
    
    inline LAP_DoubleMatrixStorage& getStorage() {
        if (mStorage.get()==null) {
            throw LAP_Exception("math/linalg/core","LAP_DoubleMatrix::getStorage()",
                                "a matrix must have a storage");
        }
        return *mStorage.get();
    } 
    inline const LAP_DoubleMatrixStorage& getStorage() const {
        if (mStorage.get()==null) {
            throw LAP_Exception("math/linalg/core","LAP_DoubleMatrix::getStorage()",
                                "a matrix must have a storage");
        }
        return *mStorage.get();
    }
    
    virtual tBoolean isSymmetric() const=0;
    virtual tBoolean isUpper() const=0;
    
public:
    virtual void getColumn(const tLVectorIndex& j,LAP_DoubleVector& v) const=0;
    virtual void getColumn(const tLVectorIndex& j,SP::LAP_DoubleVector v) const {
        if (v.get()!=null) {
            getColumn(j,*v.get());
        }
    };
    virtual void getRow(const tLVectorIndex& i,LAP_DoubleVector& v) const=0;
    virtual void getRow(const tLVectorIndex& i,SP::LAP_DoubleVector v) const {
        if (v.get()!=null) {
            getRow(i,*v.get());
        }
    };

    // Computation methods
    // ===================

    
    
    virtual tReal sum(const tFlag& d,LAP_DoubleVector& s) const=0;
    virtual tReal sum(const tFlag& d,const tLVectorIndex& index) const=0;
   
    virtual tReal norm2() const=0;
    virtual tReal norm2(LAP_DoubleVector& v) const=0;
  
    inline double trace() const {
        return getStorage().trace();
    }
    
    inline void addDiagonal(const double& alpha) {
        mStorage->addDiagonal(alpha);
    }
    
    /* \brief compute the bands number
     * @param eps:absolute module to detect null value
     * @return a pair whose firest element is the number of upper bands and the second the lower bands number
     */
    inline pair<tLVectorIndex,tLVectorIndex> computeBandsNumber() const {
        return mStorage->computeBandsNumber();
    }

    inline double getFillRate() const {
        return mStorage->getFillRate();
    }
    inline tLVectorIndex getNullValuesNumber() const {
        return mStorage->getNullValuesNumber();
    }
    inline tLVectorIndex getValuesNumber() const {
        return mStorage->getValuesNumber();
    }
  
    
public:
    //convenience methods
    //=====================
    
    inline static double zpow(const double& f,const int& n) {
        
        //x^0=1
        if (n==0) return 1;
        //0^n=0 n>0
        if (f==0) return 0;
        //n<0
        if (n<0) return zpow(1./f,-n);
        //n=1
        if (n==1) return f;
        //n>1
        double ret=f;
        for (int i=0;i<n;i++) ret*=ret;
        return ret;
    }
    
    // LINEAR Method
    //===============
    
    

  
    
    // ----------------
    // vector product
    // ----------------

public:
    virtual SP::LAP_DoubleVector vectorProduct(const LAP_DoubleVector& X) const {
        SP::LAP_DoubleVector Y=LAP_DoubleVector::New();
        vectorProduct(false,X,*Y.get());
        return Y;
    }
    virtual void vectorProduct(const LAP_DoubleVector& X,
                               LAP_DoubleVector& Y) const{
        vectorProduct(false,X,Y);
    }
    
    virtual void vectorProduct(SPC::LAP_DoubleVector& X,
                               SP::LAP_DoubleVector& Y) const {
        if ((X.get()!=null) &&
            (Y.get()!=null))
            vectorProduct(false,*X.get(),*Y.get());
        
    };
     
    
    virtual void vectorProduct(const tBoolean& isTrans,
                               const LAP_DoubleVector& X,
                               LAP_DoubleVector& Y) const{
        tLVectorIndex nRows=getRowsNumber();
        tLVectorIndex nCols=getColumnsNumber();
        tLVectorIndex nX=X.getSize();
        if (isTrans) std::swap(nRows,nCols);
        //no size
        if (nRows*nCols==0) return;

        //verify X vector dimension
        if (nCols!=nX) {
            throw LAP_Exception("math/linal/core","LAP_DoubleMatrix::product","incompatible size of vector");
        }
        //set the size of Y
        Y.setSize(nRows);
        //compute the product
        vectorProduct(false,
                      nX,X.getIncrement(),&X(0),
                      1.,0.,
                      Y.getSize(),Y.getIncrement(),&Y(0));
    };
   
     
   
  

    
    virtual void vectorProduct(const tBoolean& isTrans,
                               const tLVectorIndex& nX,const tLVectorIncrement& incX,const double *X,
                               const double& alpha,const double& beta,
                               const tLVectorIndex& nY,const tLVectorIncrement& incY,double *Y) const=0;
    
    // ----------------
    // Matrix product
    // ----------------
    virtual SP::LAP_DoubleMatrix matrixProduct(const LAP_DoubleMatrix& B) const=0;
    
    // ----------------------
    // eigen values & vectors
    // ------------------------
public:
  
    virtual tBoolean computeEigenValues(LAP_DoubleVector& U) const=0;

   
    
    
    // -----------------------
    // linear system solvers
    // -----------------------
  


    
    // -----------------------
    // STRING Representation
    // ----------------------
public:
    virtual tString toString() const;
    
  
  
};
#endif