#include <CRMatrixRect.h>

Inheritance diagram for CRMatrixRect< T_Coeff, X, B >:

Collaboration diagram for CRMatrixRect< T_Coeff, X, B >:

Classes
class	PairWiseAssembler

Public Types
typedef T_Coeff	Coeff

typedef T_Coeff	Diag

typedef T_Coeff	OffDiag

typedef Array< Coeff >	CoeffArray

typedef Array< Coeff >	DiagArray

typedef Array< X >	XArray

typedef Array< B >	BArray

typedef map< const CRConnectivity , PairWiseAssembler >	PairWiseAssemblerMap

Public Member Functions
	CRMatrixRect (const CRConnectivity &conn)

	DEFINE_TYPENAME ("CRMatrixRect<"+NumTypeTraits< Coeff >::getTypeName()+","+NumTypeTraits< X >::getTypeName()+","+NumTypeTraits< B >::getTypeName()+">")

virtual void	initAssembly ()

virtual void	multiply (IContainer &yB, const IContainer &xB) const

virtual void	multiplyAndAdd (IContainer &yB, const IContainer &xB) const

const CRConnectivity &	getConnectivity () const

Array< Coeff > &	getDiag ()

Array< Coeff > &	getOffDiag ()

const Array< Coeff > &	getDiag () const

const Array< Coeff > &	getOffDiag () const

virtual	~CRMatrixRect ()

shared_ptr< CRConnectivity >	createCoarseConnectivity (const IContainer &gCoarseIndex, const CRConnectivity &coarseToFine, const StorageSite &coarseRowSite, const StorageSite &coarseColSite)

virtual shared_ptr< Matrix >	createCoarseMatrix (const IContainer &gCoarseIndex, const CRConnectivity &coarseToFine, const CRConnectivity &coarseConnectivity)

PairWiseAssembler &	getPairWiseAssembler (const CRConnectivity &pairs)

Public Member Functions inherited from Matrix
	Matrix ()

virtual	~Matrix ()

	DEFINE_TYPENAME ("Matrix")

virtual shared_ptr< ArrayBase >	quadProduct (const IContainer &xB) const

virtual void	forwardGS (IContainer &xB, IContainer &bB, IContainer &residual) const

virtual void	reverseGS (IContainer &xB, IContainer &bB, IContainer &residual) const

virtual void	Jacobi (IContainer &xnew, const IContainer &xold, const IContainer &b) const

virtual void	iluSolve (IContainer &xB, const IContainer &bB, const IContainer &residual) const

virtual void	spikeSolve (IContainer &xB, const IContainer &bB, const IContainer &residual, const SpikeStorage &spike_storage) const

virtual void	solveBoundary (IContainer &xB, IContainer &bB, IContainer &residual) const

virtual void	computeResidual (const IContainer &xB, const IContainer &bB, IContainer &residual) const

virtual int	createCoarsening (IContainer &coarseIndex, const int groupSize, const double weighRatioThreshold)

virtual void	eliminateBoundaryEquations (IContainer &xB)

virtual void	printRow (const int nr) const

virtual void *	getDiagData () const

virtual void *	getOffDiagData () const

virtual int	getDiagDataSize () const

virtual int	getOffDiagDataSize () const

virtual shared_ptr< Matrix >	createMergeMatrix (const LinearSystemMerger &mergeLS)

virtual void	setFlatMatrix (Matrix &fmg) const

virtual void	transpose ()

virtual bool	isInvertible ()

Private Attributes
const CRConnectivity &	_conn

const Array< int > &	_row

const Array< int > &	_col

Array< Coeff >	_diag

Array< Coeff >	_offDiag

PairWiseAssemblerMap	_pairWiseAssemblers

Detailed Description

template<class T_Coeff, class X, class B>
class CRMatrixRect< T_Coeff, X, B >

Sparse matrix stored using a compressed row format. The sparsity pattern is provided by a CRConnectivity object that is required at construction time. Note that we assume that diagonal is implicitly always present (i.e., col[row[i]] through col[row[i+1]] do not contain i) and store diagonal and off diagonal entries in separate arrays.

The class is templated with the types of diagonal, off-diagonal and the corresponsing x arrays as template parameters.

Definition at line 28 of file CRMatrixRect.h.

Member Typedef Documentation

template<class T_Coeff , class X , class B >

typedef Array<B> CRMatrixRect< T_Coeff, X, B >::BArray

Definition at line 37 of file CRMatrixRect.h.

template<class T_Coeff , class X , class B >

typedef T_Coeff CRMatrixRect< T_Coeff, X, B >::Coeff

Definition at line 31 of file CRMatrixRect.h.

template<class T_Coeff , class X , class B >

typedef Array<Coeff> CRMatrixRect< T_Coeff, X, B >::CoeffArray

Definition at line 34 of file CRMatrixRect.h.

template<class T_Coeff , class X , class B >

typedef T_Coeff CRMatrixRect< T_Coeff, X, B >::Diag

Definition at line 32 of file CRMatrixRect.h.

template<class T_Coeff , class X , class B >

typedef Array<Coeff> CRMatrixRect< T_Coeff, X, B >::DiagArray

Definition at line 35 of file CRMatrixRect.h.

template<class T_Coeff , class X , class B >

typedef T_Coeff CRMatrixRect< T_Coeff, X, B >::OffDiag

Definition at line 33 of file CRMatrixRect.h.

template<class T_Coeff , class X , class B >

typedef map<const CRConnectivity*,PairWiseAssembler*> CRMatrixRect< T_Coeff, X, B >::PairWiseAssemblerMap

Definition at line 92 of file CRMatrixRect.h.

template<class T_Coeff , class X , class B >

typedef Array<X> CRMatrixRect< T_Coeff, X, B >::XArray

Definition at line 36 of file CRMatrixRect.h.

Constructor & Destructor Documentation

template<class T_Coeff , class X , class B >

CRMatrixRect< T_Coeff, X, B >::CRMatrixRect ( const CRConnectivity & conn )

inline

Definition at line 95 of file CRMatrixRect.h.

References logCtor.

Referenced by CRMatrixRect< T_Coeff, X, B >::createCoarseMatrix().

                                            :
     Matrix(),
     _conn(conn),
     _row(_conn.getRow()),
     _col(_conn.getCol()),
     _diag(_conn.getRowDim()),
     _offDiag(_col.getLength()),
     _pairWiseAssemblers()
   {
     logCtor();
   }

template<class T_Coeff , class X , class B >

virtual CRMatrixRect< T_Coeff, X, B >::~CRMatrixRect ( )

inlinevirtual

Definition at line 173 of file CRMatrixRect.h.

References logDtor.

   {
     logDtor();
   }

Member Function Documentation

template<class T_Coeff , class X , class B >

shared_ptr<CRConnectivity> CRMatrixRect< T_Coeff, X, B >::createCoarseConnectivity	(	const IContainer &	gCoarseIndex,
		const CRConnectivity &	coarseToFine,
		const StorageSite &	coarseRowSite,
		const StorageSite &	coarseColSite
	)

inlinevirtual

Create the connectivity for a coarse level matrix given the coarsening (ie. the fine to coarse index mapping) and its transpose (ie coarse to fine mapping, provided as a CRConnectivity object)

Reimplemented from Matrix.

Definition at line 187 of file CRMatrixRect.h.

References CRMatrixRect< T_Coeff, X, B >::_col, CRMatrixRect< T_Coeff, X, B >::_row, StorageSite::getCount(), and CRConnectivity::getCount().

   {
     const Array<int>&  coarseIndex =
       dynamic_cast<const Array<int>& >(gCoarseIndex);
 
     const int nCoarseRows = coarseRowSite.getCount();
     
     shared_ptr<CRConnectivity> coarseCR(new CRConnectivity(coarseRowSite,coarseColSite));
 
 
     coarseCR->initCount();
 
     Array<bool> coarseCounted(nCoarseRows);
 
     coarseCounted = false;
     
     for(int nrCoarse=0; nrCoarse<nCoarseRows; nrCoarse++)
     {
         const int nFine = coarseToFine.getCount(nrCoarse);
         for(int nfr=0; nfr<nFine; nfr++)
         {
             const int nrFine = coarseToFine(nrCoarse,nfr);
 
             for (int nb = _row[nrFine]; nb<_row[nrFine+1]; nb++)
             {
                 const int nc = _col[nb];
                 const int ncCoarse = coarseIndex[nc];
                 if (ncCoarse>=0 && nrCoarse!=ncCoarse && !coarseCounted[ncCoarse])
                 {
                     coarseCounted[ncCoarse] = true;
                     coarseCR->addCount(nrCoarse,1);
                 }
             }
         }
 
         // reset counted 
         for(int nfr=0; nfr<nFine; nfr++)
         {
             const int nrFine = coarseToFine(nrCoarse,nfr);
             
             for (int nb = _row[nrFine]; nb<_row[nrFine+1]; nb++)
             {
                 const int nc = _col[nb];
                 const int ncCoarse = coarseIndex[nc];
                 if (ncCoarse>=0)
                   coarseCounted[ncCoarse] = false;
             }
         }
     }
 
     coarseCR->finishCount();
 
     for(int nrCoarse=0; nrCoarse<nCoarseRows; nrCoarse++)
     {
         const int nFine = coarseToFine.getCount(nrCoarse);
         for(int nfr=0; nfr<nFine; nfr++)
         {
             const int nrFine = coarseToFine(nrCoarse,nfr);
 
             for (int nb = _row[nrFine]; nb<_row[nrFine+1]; nb++)
             {
                 const int nc = _col[nb];
                 const int ncCoarse = coarseIndex[nc];
                 if (ncCoarse>=0 && nrCoarse!=ncCoarse && !coarseCounted[ncCoarse])
                 {
                     coarseCounted[ncCoarse] = true;
                     coarseCR->add(nrCoarse,ncCoarse);
                 }
             }
         }
 
         // reset counted 
         for(int nfr=0; nfr<nFine; nfr++)
         {
             const int nrFine = coarseToFine(nrCoarse,nfr);
             
             for (int nb = _row[nrFine]; nb<_row[nrFine+1]; nb++)
             {
                 const int nc = _col[nb];
                 const int ncCoarse = coarseIndex[nc];
                 if (ncCoarse>=0)
                 coarseCounted[ncCoarse] = false;
             }
         }
     }
 
     coarseCR->finishAdd();
     return coarseCR;
   }

template<class T_Coeff , class X , class B >

virtual shared_ptr<Matrix> CRMatrixRect< T_Coeff, X, B >::createCoarseMatrix	(	const IContainer &	gCoarseIndex,
		const CRConnectivity &	coarseToFine,
		const CRConnectivity &	coarseConnectivity
	)

inlinevirtual

create the coarse matrix given the coarsening and the coarse level connectivity.

Reimplemented from Matrix.

Definition at line 288 of file CRMatrixRect.h.

References CRMatrixRect< T_Coeff, X, B >::_col, CRMatrixRect< T_Coeff, X, B >::_diag, CRMatrixRect< T_Coeff, X, B >::_offDiag, CRMatrixRect< T_Coeff, X, B >::_row, CRMatrixRect< T_Coeff, X, B >::CRMatrixRect(), CRConnectivity::getCol(), CRConnectivity::getCount(), CRConnectivity::getRow(), CRConnectivity::getRowDim(), and Array< T >::zero().

   {
     const Array<int>&  coarseIndex =
       dynamic_cast<const Array<int>& >(gCoarseIndex);
     
     const int nCoarseRows = coarseConnectivity.getRowDim();
 
     shared_ptr<CRMatrixRect> coarseMatrix(new CRMatrixRect(coarseConnectivity));
 
     Array<Diag>& coarseDiag = coarseMatrix->getDiag();
     Array<OffDiag>& coarseOffDiag = coarseMatrix->getOffDiag();
 
     const Array<int>& coarseConnRow = coarseConnectivity.getRow();
     const Array<int>& coarseConnCol = coarseConnectivity.getCol();
 
     coarseDiag.zero();
     coarseOffDiag.zero();
 
     //used to avoid searches when inserting coeffs
     Array<int> coarseCoeffPos(nCoarseRows);
 
 
     for(int nrCoarse=0; nrCoarse<nCoarseRows; nrCoarse++)
     {
         // for easy indexing when inserting coefficients set col
         // positions into the coarse connectivity
         for(int nb=coarseConnRow[nrCoarse]; nb<coarseConnRow[nrCoarse+1]; nb++)
           coarseCoeffPos[coarseConnCol[nb]] = nb;
 
         // loop over the fine rows that make up this coarse row
         const int nFine = coarseToFine.getCount(nrCoarse);
         for(int nfr=0; nfr<nFine; nfr++)
         {
             const int nrFine = coarseToFine(nrCoarse,nfr);
 
             coarseDiag[nrCoarse] += _diag[nrFine];
             
             for (int nb = _row[nrFine]; nb<_row[nrFine+1]; nb++)
             {
                 const int nc = _col[nb];
                 const int ncCoarse = coarseIndex[nc];
 
                 if (ncCoarse<0) continue;
                 
                 if (nrCoarse!=ncCoarse)
                 {
                     const int pos = coarseCoeffPos[ncCoarse];
                     coarseOffDiag[pos] += _offDiag[nb];
                 }
                 else 
                 {
                     coarseDiag[nrCoarse] += _offDiag[nb];
                 }
             }
         }
     }
 
     return coarseMatrix;
   }

template<class T_Coeff , class X , class B >

CRMatrixRect< T_Coeff, X, B >::DEFINE_TYPENAME ( "CRMatrixRect<"+NumTypeTraits< Coeff >::getTypeName()+","+NumTypeTraits< X >::getTypeName()+","+NumTypeTraits< B >::getTypeName()+">" )

template<class T_Coeff , class X , class B >

const CRConnectivity& CRMatrixRect< T_Coeff, X, B >::getConnectivity ( ) const

inlinevirtual

Reimplemented from Matrix.

Definition at line 165 of file CRMatrixRect.h.

References CRMatrixRect< T_Coeff, X, B >::_conn.

165 {return _conn;}

CRMatrixRect::_conn

const CRConnectivity & _conn

Definition: CRMatrixRect.h:363

template<class T_Coeff , class X , class B >

Array<Coeff>& CRMatrixRect< T_Coeff, X, B >::getDiag ( )

inline

Definition at line 167 of file CRMatrixRect.h.

References CRMatrixRect< T_Coeff, X, B >::_diag.

167 {return _diag;}

CRMatrixRect::_diag

Array< Coeff > _diag

Definition: CRMatrixRect.h:366

template<class T_Coeff , class X , class B >

const Array<Coeff>& CRMatrixRect< T_Coeff, X, B >::getDiag ( ) const

inline

Definition at line 170 of file CRMatrixRect.h.

References CRMatrixRect< T_Coeff, X, B >::_diag.

170 {return _diag;}

CRMatrixRect::_diag

Array< Coeff > _diag

Definition: CRMatrixRect.h:366

template<class T_Coeff , class X , class B >

Array<Coeff>& CRMatrixRect< T_Coeff, X, B >::getOffDiag ( )

inline

Definition at line 168 of file CRMatrixRect.h.

References CRMatrixRect< T_Coeff, X, B >::_offDiag.

168 {return _offDiag;}

CRMatrixRect::_offDiag

Array< Coeff > _offDiag

Definition: CRMatrixRect.h:367

template<class T_Coeff , class X , class B >

const Array<Coeff>& CRMatrixRect< T_Coeff, X, B >::getOffDiag ( ) const

inline

Definition at line 171 of file CRMatrixRect.h.

References CRMatrixRect< T_Coeff, X, B >::_offDiag.

171 {return _offDiag;}

CRMatrixRect::_offDiag

Array< Coeff > _offDiag

Definition: CRMatrixRect.h:367

template<class T_Coeff , class X , class B >

PairWiseAssembler& CRMatrixRect< T_Coeff, X, B >::getPairWiseAssembler ( const CRConnectivity & pairs )

inline

Definition at line 350 of file CRMatrixRect.h.

References CRMatrixRect< T_Coeff, X, B >::_conn, CRMatrixRect< T_Coeff, X, B >::_offDiag, CRMatrixRect< T_Coeff, X, B >::_pairWiseAssemblers, and CRConnectivity::getPairToColMapping().

   {
     if (_pairWiseAssemblers.find(&pairs) == _pairWiseAssemblers.end())
     {
         _pairWiseAssemblers[&pairs] =
           new PairWiseAssembler(_offDiag,
                                 _conn.getPairToColMapping(pairs));
     }
     return *_pairWiseAssemblers[&pairs];
   }

template<class T_Coeff , class X , class B >

virtual void CRMatrixRect< T_Coeff, X, B >::initAssembly ( )

inlinevirtual

Implements Matrix.

Definition at line 114 of file CRMatrixRect.h.

References CRMatrixRect< T_Coeff, X, B >::_diag, CRMatrixRect< T_Coeff, X, B >::_offDiag, and Array< T >::zero().

   {
     _diag.zero();
     _offDiag.zero();
   }

template<class T_Coeff , class X , class B >

virtual void CRMatrixRect< T_Coeff, X, B >::multiply	(	IContainer &	yB,
		const IContainer &	xB
	)		const

inlinevirtual

y = this * x

Reimplemented from Matrix.

Definition at line 125 of file CRMatrixRect.h.

References CRMatrixRect< T_Coeff, X, B >::_col, CRMatrixRect< T_Coeff, X, B >::_conn, CRMatrixRect< T_Coeff, X, B >::_diag, CRMatrixRect< T_Coeff, X, B >::_offDiag, CRMatrixRect< T_Coeff, X, B >::_row, CRConnectivity::getRowSite(), and StorageSite::getSelfCount().

   {
     BArray& y = dynamic_cast<BArray&>(yB);
     const XArray& x = dynamic_cast<const XArray&>(xB);
     
     const int nRows = _conn.getRowSite().getSelfCount();
     for(int nr=0; nr<nRows; nr++)
     {
         y[nr] = _diag[nr]*x[nr];
         for (int nb = _row[nr]; nb<_row[nr+1]; nb++)
         {
             const int j = _col[nb];
             y[nr] += _offDiag[nb]*x[j];
         }
     }
   }

template<class T_Coeff , class X , class B >

virtual void CRMatrixRect< T_Coeff, X, B >::multiplyAndAdd	(	IContainer &	yB,
		const IContainer &	xB
	)		const

inlinevirtual

y += this * x

Reimplemented from Matrix.

Definition at line 147 of file CRMatrixRect.h.

References CRMatrixRect< T_Coeff, X, B >::_col, CRMatrixRect< T_Coeff, X, B >::_conn, CRMatrixRect< T_Coeff, X, B >::_diag, CRMatrixRect< T_Coeff, X, B >::_offDiag, CRMatrixRect< T_Coeff, X, B >::_row, CRConnectivity::getRowSite(), and StorageSite::getSelfCount().

   {
     BArray& y = dynamic_cast<BArray&>(yB);
     const XArray& x = dynamic_cast<const XArray&>(xB);
     
     const int nRows = _conn.getRowSite().getSelfCount();
     for(int nr=0; nr<nRows; nr++)
     {
         y[nr] += _diag[nr]*x[nr];
         for (int nb = _row[nr]; nb<_row[nr+1]; nb++)
         {
             const int j = _col[nb];
             y[nr] += _offDiag[nb]*x[j];
         }
     }
   }