"Fossies" - the Fresh Open Source Software Archive  

Source code changes of the file "dune/istl/colcompmatrix.hh" between
dune-istl-2.6.0.tar.gz and dune-istl-2.7.0.tar.gz

About: dune-istl - DUNE (Distributed and Unified Numerics Environment) is a modular toolbox for solving partial differential equations (PDEs) with grid-based methods: DUNE Iterative Solver Template Library.

colcompmatrix.hh  (dune-istl-2.6.0):colcompmatrix.hh  (dune-istl-2.7.0)
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*- // -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2: // vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_ISTL_COLCOMPMATRIX_HH #ifndef DUNE_ISTL_COLCOMPMATRIX_HH
#define DUNE_ISTL_COLCOMPMATRIX_HH #define DUNE_ISTL_COLCOMPMATRIX_HH
#include "bcrsmatrix.hh" #include "bcrsmatrix.hh"
#include "bvector.hh" #include "bvector.hh"
#include <dune/common/fmatrix.hh> #include <dune/common/fmatrix.hh>
#include <dune/common/fvector.hh> #include <dune/common/fvector.hh>
#include <dune/common/typetraits.hh> #include <dune/common/typetraits.hh>
#include <dune/common/unused.hh> #include <dune/common/unused.hh>
#include <dune/common/scalarmatrixview.hh>
#include <limits> #include <limits>
namespace Dune namespace Dune
{ {
/** /**
* @brief Provides access to an iterator over all matrix rows. * @brief Provides access to an iterator over all matrix rows.
* *
* @tparam M The type of the matrix. * @tparam M The type of the matrix.
*/ */
template<class M> template<class M>
skipping to change at line 139 skipping to change at line 140
} }
private: private:
//! @brief The matrix for which we manage the row subset. //! @brief The matrix for which we manage the row subset.
const Matrix& m_; const Matrix& m_;
//! @brief The set of row indices we manage. //! @brief The set of row indices we manage.
const RowIndexSet& s_; const RowIndexSet& s_;
}; };
/** /**
* @brief Utility class for converting an ISTL Matrix into a column-compressed
matrix
* @tparam M the matrix type
*/
template<class M>
struct ColCompMatrix
{};
/**
* @brief Inititializer for the ColCompMatrix * @brief Inititializer for the ColCompMatrix
* as needed by OverlappingSchwarz * as needed by OverlappingSchwarz
* @tparam M the matrix type * @tparam M the matrix type
* @tparam I the internal index type
*/ */
template<class M> template<class M, class I = int>
struct ColCompMatrixInitializer class ColCompMatrixInitializer;
{};
template<class M, class X, class TM, class TD, class T1> template<class M, class X, class TM, class TD, class T1>
class SeqOverlappingSchwarz; class SeqOverlappingSchwarz;
template<class T, bool flag> template<class T, bool flag>
struct SeqOverlappingSchwarzAssemblerHelper; struct SeqOverlappingSchwarzAssemblerHelper;
/** /**
* @brief Converter for BCRSMatrix to column-compressed Matrix. * @brief Utility class for converting an ISTL Matrix into a column-compressed
* specialization for BCRSMatrix matrix
* @tparam M The matrix type
* @tparam I the internal index type
*/ */
template<class B, class TA, int n, int m> template<class Mat, class I = int>
class ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> > class ColCompMatrix
{ {
friend struct ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<B,n,m>,TA> >; friend class ColCompMatrixInitializer<Mat, I>;
using B = typename Mat::field_type;
public: public:
/** @brief The type of the matrix to convert. */ /** @brief The type of the matrix to convert. */
typedef BCRSMatrix<FieldMatrix<B,n,m>,TA> Matrix; using Matrix = Mat;
friend struct SeqOverlappingSchwarzAssemblerHelper<ColCompMatrix<Matrix>, tr ue>; friend struct SeqOverlappingSchwarzAssemblerHelper<ColCompMatrix<Matrix>, tr ue>;
typedef typename Matrix::size_type size_type; typedef typename Matrix::size_type size_type;
using Index = I;
/** /**
* @brief Constructor that initializes the data. * @brief Constructor that initializes the data.
* @param mat The matrix to convert. * @param mat The matrix to convert.
*/ */
explicit ColCompMatrix(const Matrix& mat); explicit ColCompMatrix(const Matrix& mat);
ColCompMatrix(); ColCompMatrix();
/** @brief Destructor */ /** @brief Destructor */
virtual ~ColCompMatrix(); virtual ~ColCompMatrix();
skipping to change at line 199 skipping to change at line 198
* @brief Get the number of rows. * @brief Get the number of rows.
* @return The number of rows. * @return The number of rows.
*/ */
size_type N() const size_type N() const
{ {
return N_; return N_;
} }
size_type nnz() const size_type nnz() const
{ {
// TODO: The following code assumes that the blocks are dense
// and that they all have the same dimensions.
typename Matrix::block_type dummy;
const auto n = MatrixDimension<typename Matrix::block_type>::rowdim(dummy)
;
const auto m = MatrixDimension<typename Matrix::block_type>::coldim(dummy)
;
return Nnz_/n/m; return Nnz_/n/m;
} }
/** /**
* @brief Get the number of columns. * @brief Get the number of columns.
* @return The number of columns. * @return The number of columns.
*/ */
size_type M() const size_type M() const
{ {
return M_; return M_;
} }
B* getValues() const B* getValues() const
{ {
return values; return values;
} }
int* getRowIndex() const Index* getRowIndex() const
{ {
return rowindex; return rowindex;
} }
int* getColStart() const Index* getColStart() const
{ {
return colstart; return colstart;
} }
ColCompMatrix& operator=(const Matrix& mat); ColCompMatrix& operator=(const Matrix& mat);
ColCompMatrix& operator=(const ColCompMatrix& mat); ColCompMatrix& operator=(const ColCompMatrix& mat);
/** /**
* @brief Initialize data from a given set of matrix rows and columns * @brief Initialize data from a given set of matrix rows and columns
* @param mat the matrix with the values * @param mat the matrix with the values
* @param mrs The set of row (and column) indices to remove * @param mrs The set of row (and column) indices to remove
*/ */
virtual void setMatrix(const Matrix& mat, const std::set<std::size_t>& mrs); virtual void setMatrix(const Matrix& mat, const std::set<std::size_t>& mrs);
/** @brief free allocated space. */ /** @brief free allocated space. */
virtual void free(); virtual void free();
/** @brief Initialize data from given matrix. */ /** @brief Initialize data from given matrix. */
virtual void setMatrix(const Matrix& mat); virtual void setMatrix(const Matrix& mat);
public: public:
int N_, M_, Nnz_; size_type N_, M_, Nnz_;
B* values; B* values;
int* rowindex; Index* rowindex;
int* colstart; Index* colstart;
}; };
template<class T, class A, int n, int m> template<class M, class I>
class ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> > class ColCompMatrixInitializer
{ {
template<class I, class S, class D> template<class IList, class S, class D>
friend class OverlappingSchwarzInitializer; friend class OverlappingSchwarzInitializer;
public: public:
typedef Dune::BCRSMatrix<FieldMatrix<T,n,m>,A> Matrix; using Matrix = M;
typedef Dune::ColCompMatrix<Matrix> ColCompMatrix; using Index = I;
typedef Dune::ColCompMatrix<Matrix, Index> ColCompMatrix;
typedef typename Matrix::row_type::const_iterator CIter; typedef typename Matrix::row_type::const_iterator CIter;
typedef typename Matrix::size_type size_type; typedef typename Matrix::size_type size_type;
ColCompMatrixInitializer(ColCompMatrix& lum); /** \brief Constructor for scalar-valued matrices
*
* \tparam Block Dummy parameter to make SFINAE work
*/
template <class Block=typename M::block_type>
ColCompMatrixInitializer(ColCompMatrix& lum,
typename std::enable_if_t<Dune::IsNumber<Block>::va
lue>* sfinae = nullptr);
/** \brief Constructor for dense matrix-valued matrices
*
* \tparam Block Dummy parameter to make SFINAE work
*/
template <class Block=typename M::block_type>
ColCompMatrixInitializer(ColCompMatrix& lum,
typename std::enable_if_t<!Dune::IsNumber<Block>::v
alue>* sfinae = nullptr);
ColCompMatrixInitializer(); ColCompMatrixInitializer();
virtual ~ColCompMatrixInitializer(); virtual ~ColCompMatrixInitializer()
{}
template<typename Iter> template<typename Iter>
void addRowNnz(const Iter& row) const; void addRowNnz(const Iter& row) const;
template<typename Iter, typename Set> template<typename Iter, typename FullMatrixIndex>
void addRowNnz(const Iter& row, const Set& s) const; void addRowNnz(const Iter& row, const std::set<FullMatrixIndex>& indices) co
nst;
template<typename Iter, typename SubMatrixIndex>
void addRowNnz(const Iter& row, const std::vector<SubMatrixIndex>& indices)
const;
void allocate(); void allocate();
template<typename Iter> template<typename Iter>
void countEntries(const Iter& row, const CIter& col) const; void countEntries(const Iter& row, const CIter& col) const;
void countEntries(size_type colidx) const; void countEntries(size_type colidx) const;
void calcColstart() const; void calcColstart() const;
skipping to change at line 295 skipping to change at line 318
virtual void createMatrix() const; virtual void createMatrix() const;
protected: protected:
void allocateMatrixStorage() const; void allocateMatrixStorage() const;
void allocateMarker(); void allocateMarker();
ColCompMatrix* mat; ColCompMatrix* mat;
size_type cols; size_type cols;
mutable size_type *marker;
// Number of rows/columns of the matrix entries
// (assumed to be scalars or dense matrices)
size_type n, m;
mutable std::vector<size_type> marker;
}; };
template<class T, class A, int n, int m> template<class M, class I>
ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::ColCompMatrixInit template <class Block>
ializer(ColCompMatrix& mat_) ColCompMatrixInitializer<M, I>::ColCompMatrixInitializer(ColCompMatrix& mat_,t
: mat(&mat_), cols(mat_.M()), marker(0) ypename std::enable_if_t<Dune::IsNumber<Block>::value>* sfinae)
: mat(&mat_), cols(mat_.M())
{ {
n = 1;
m = 1;
mat->Nnz_=0; mat->Nnz_=0;
} }
template<class T, class A, int n, int m> template<class M, class I>
ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::ColCompMatrixInit template <class Block>
ializer() ColCompMatrixInitializer<M, I>::ColCompMatrixInitializer(ColCompMatrix& mat_,t
: mat(0), cols(0), marker(0) ypename std::enable_if_t<!Dune::IsNumber<Block>::value>* sfinae)
{} : mat(&mat_), cols(mat_.M())
{
// WARNING: This assumes that all blocks are dense and identical
n = M::block_type::rows;
m = M::block_type::cols;
template<class T, class A, int n, int m> mat->Nnz_=0;
ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::~ColCompMatrixIni
tializer()
{
if(marker)
delete[] marker;
} }
template<class T, class A, int n, int m> template<class M, class I>
ColCompMatrixInitializer<M, I>::ColCompMatrixInitializer()
: mat(0), cols(0), n(0), m(0)
{}
template<class M, class I>
template<typename Iter> template<typename Iter>
void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::addRowNnz(co nst Iter& row) const void ColCompMatrixInitializer<M, I>::addRowNnz(const Iter& row) const
{ {
mat->Nnz_+=row->getsize(); mat->Nnz_+=row->getsize();
} }
template<class T, class A, int n, int m> template<class M, class I>
template<typename Iter, typename Map> template<typename Iter, typename FullMatrixIndex>
void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::addRowNnz(co void ColCompMatrixInitializer<M, I>::addRowNnz(const Iter& row,
nst Iter& row, cons
cons t std::set<FullMatrixIndex>& indices) const
t Map& indices) const
{ {
typedef typename Iter::value_type::const_iterator RIter; typedef typename Iter::value_type::const_iterator RIter;
typedef typename Map::const_iterator MIter; typedef typename std::set<FullMatrixIndex>::const_iterator MIter;
MIter siter =indices.begin(); MIter siter =indices.begin();
for(RIter entry=row->begin(); entry!=row->end(); ++entry) for(RIter entry=row->begin(); entry!=row->end(); ++entry)
{ {
for(; siter!=indices.end() && *siter<entry.index(); ++siter) ; for(; siter!=indices.end() && *siter<entry.index(); ++siter) ;
if(siter==indices.end()) if(siter==indices.end())
break; break;
if(*siter==entry.index()) if(*siter==entry.index())
// index is in subdomain // index is in subdomain
++mat->Nnz_; ++mat->Nnz_;
} }
} }
template<class T, class A, int n, int m> template<class M, class I>
void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::allocate() template<typename Iter, typename SubMatrixIndex>
void ColCompMatrixInitializer<M, I>::addRowNnz(const Iter& row,
cons
t std::vector<SubMatrixIndex>& indices) const
{
using RIter = typename Iter::value_type::const_iterator;
for(RIter entry=row->begin(); entry!=row->end(); ++entry)
if (indices[entry.index()]!=std::numeric_limits<SubMatrixIndex>::max())
++mat->Nnz_;
}
template<class M, class I>
void ColCompMatrixInitializer<M, I>::allocate()
{ {
allocateMatrixStorage(); allocateMatrixStorage();
allocateMarker(); allocateMarker();
} }
template<class T, class A, int n, int m> template<class M, class I>
void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::allocateMatr void ColCompMatrixInitializer<M, I>::allocateMatrixStorage() const
ixStorage() const
{ {
mat->Nnz_*=n*m; mat->Nnz_*=n*m;
// initialize data // initialize data
mat->values=new T[mat->Nnz_]; mat->values=new typename M::field_type[mat->Nnz_];
mat->rowindex=new int[mat->Nnz_]; mat->rowindex=new I[mat->Nnz_];
mat->colstart=new int[cols+1]; mat->colstart=new I[cols+1];
} }
template<class T, class A, int n, int m> template<class M, class I>
void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::allocateMark void ColCompMatrixInitializer<M, I>::allocateMarker()
er()
{ {
marker = new typename Matrix::size_type[cols]; marker.resize(cols);
std::fill(marker.begin(), marker.end(), 0);
for(size_type i=0; i < cols; ++i)
marker[i]=0;
} }
template<class T, class A, int n, int m> template<class M, class I>
template<typename Iter> template<typename Iter>
void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::countEntries (const Iter& row, const CIter& col) const void ColCompMatrixInitializer<M, I>::countEntries(const Iter& row, const CIter & col) const
{ {
DUNE_UNUSED_PARAMETER(row); DUNE_UNUSED_PARAMETER(row);
countEntries(col.index()); countEntries(col.index());
} }
template<class T, class A, int n, int m> template<class M, class I>
void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::countEntries void ColCompMatrixInitializer<M, I>::countEntries(size_type colindex) const
(size_type colindex) const
{ {
for(size_type i=0; i < m; ++i) for(size_type i=0; i < m; ++i)
{ {
assert(colindex*m+i<cols); assert(colindex*m+i<cols);
marker[colindex*m+i]+=n; marker[colindex*m+i]+=n;
} }
} }
template<class T, class A, int n, int m> template<class M, class I>
void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::calcColstart void ColCompMatrixInitializer<M, I>::calcColstart() const
() const
{ {
mat->colstart[0]=0; mat->colstart[0]=0;
for(size_type i=0; i < cols; ++i) { for(size_type i=0; i < cols; ++i) {
assert(i<cols); assert(i<cols);
mat->colstart[i+1]=mat->colstart[i]+marker[i]; mat->colstart[i+1]=mat->colstart[i]+marker[i];
marker[i]=mat->colstart[i]; marker[i]=mat->colstart[i];
} }
} }
template<class T, class A, int n, int m> template<class M, class I>
template<typename Iter> template<typename Iter>
void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::copyValue(co nst Iter& row, const CIter& col) const void ColCompMatrixInitializer<M, I>::copyValue(const Iter& row, const CIter& c ol) const
{ {
copyValue(col, row.index(), col.index()); copyValue(col, row.index(), col.index());
} }
template<class T, class A, int n, int m> template<class M, class I>
void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::copyValue(co void ColCompMatrixInitializer<M, I>::copyValue(const CIter& col, size_type row
nst CIter& col, size_type rowindex, size_type colindex) const index, size_type colindex) const
{ {
for(size_type i=0; i<n; i++) { for(size_type i=0; i<n; i++) {
for(size_type j=0; j<m; j++) { for(size_type j=0; j<m; j++) {
assert(colindex*m+j<cols-1 || (int)marker[colindex*m+j]<mat->colstart[co assert(colindex*m+j<cols-1 || (size_type)marker[colindex*m+j]<(size_type
lindex*m+j+1]); )mat->colstart[colindex*m+j+1]);
assert((int)marker[colindex*m+j]<mat->Nnz_); assert((size_type)marker[colindex*m+j]<mat->Nnz_);
mat->rowindex[marker[colindex*m+j]]=rowindex*n+i; mat->rowindex[marker[colindex*m+j]]=rowindex*n+i;
mat->values[marker[colindex*m+j]]=(*col)[i][j]; mat->values[marker[colindex*m+j]]=Impl::asMatrix(*col)[i][j];
++marker[colindex*m+j]; // index for next entry in column ++marker[colindex*m+j]; // index for next entry in column
} }
} }
} }
template<class T, class A, int n, int m> template<class M, class I>
void ColCompMatrixInitializer<BCRSMatrix<FieldMatrix<T,n,m>,A> >::createMatrix void ColCompMatrixInitializer<M, I>::createMatrix() const
() const
{ {
delete[] marker; marker.clear();
marker=0;
} }
template<class F, class MRS> template<class F, class MRS>
void copyToColCompMatrix(F& initializer, const MRS& mrs) void copyToColCompMatrix(F& initializer, const MRS& mrs)
{ {
typedef typename MRS::const_iterator Iter; typedef typename MRS::const_iterator Iter;
typedef typename std::iterator_traits<Iter>::value_type::const_iterator CIt er; typedef typename std::iterator_traits<Iter>::value_type::const_iterator CIt er;
for(Iter row=mrs.begin(); row!= mrs.end(); ++row) for(Iter row=mrs.begin(); row!= mrs.end(); ++row)
initializer.addRowNnz(row); initializer.addRowNnz(row);
skipping to change at line 463 skipping to change at line 508
template<class F, class M,class S> template<class F, class M,class S>
void copyToColCompMatrix(F& initializer, const MatrixRowSubset<M,S>& mrs) void copyToColCompMatrix(F& initializer, const MatrixRowSubset<M,S>& mrs)
{ {
typedef MatrixRowSubset<M,S> MRS; typedef MatrixRowSubset<M,S> MRS;
typedef typename MRS::RowIndexSet SIS; typedef typename MRS::RowIndexSet SIS;
typedef typename SIS::const_iterator SIter; typedef typename SIS::const_iterator SIter;
typedef typename MRS::const_iterator Iter; typedef typename MRS::const_iterator Iter;
typedef typename std::iterator_traits<Iter>::value_type row_type; typedef typename std::iterator_traits<Iter>::value_type row_type;
typedef typename row_type::const_iterator CIter; typedef typename row_type::const_iterator CIter;
// Calculate upper Bound for nonzeros
for(Iter row=mrs.begin(); row!= mrs.end(); ++row)
initializer.addRowNnz(row, mrs.rowIndexSet());
initializer.allocate();
typedef typename MRS::Matrix::size_type size_type; typedef typename MRS::Matrix::size_type size_type;
// A vector containing the corresponding indices in // A vector containing the corresponding indices in
// the to create submatrix. // the to create submatrix.
// If an entry is the maximum of size_type then this index will not appear i n // If an entry is the maximum of size_type then this index will not appear i n
// the submatrix. // the submatrix.
std::vector<size_type> subMatrixIndex(mrs.matrix().N(), std::vector<size_type> subMatrixIndex(mrs.matrix().N(),
std::numeric_limits<size_type>::max()) ; std::numeric_limits<size_type>::max()) ;
size_type s=0; size_type s=0;
for(SIter index = mrs.rowIndexSet().begin(); index!=mrs.rowIndexSet().end(); ++index) for(SIter index = mrs.rowIndexSet().begin(); index!=mrs.rowIndexSet().end(); ++index)
subMatrixIndex[*index]=s++; subMatrixIndex[*index]=s++;
// Calculate upper Bound for nonzeros
for(Iter row=mrs.begin(); row!= mrs.end(); ++row)
initializer.addRowNnz(row, subMatrixIndex);
initializer.allocate();
for(Iter row=mrs.begin(); row!= mrs.end(); ++row) for(Iter row=mrs.begin(); row!= mrs.end(); ++row)
for(CIter col=row->begin(); col != row->end(); ++col) { for(CIter col=row->begin(); col != row->end(); ++col) {
if(subMatrixIndex[col.index()]!=std::numeric_limits<size_type>::max()) if(subMatrixIndex[col.index()]!=std::numeric_limits<size_type>::max())
// This column is in our subset (use submatrix column index) // This column is in our subset (use submatrix column index)
initializer.countEntries(subMatrixIndex[col.index()]); initializer.countEntries(subMatrixIndex[col.index()]);
} }
initializer.calcColstart(); initializer.calcColstart();
for(Iter row=mrs.begin(); row!= mrs.end(); ++row) for(Iter row=mrs.begin(); row!= mrs.end(); ++row)
for(CIter col=row->begin(); col != row->end(); ++col) { for(CIter col=row->begin(); col != row->end(); ++col) {
if(subMatrixIndex[col.index()]!=std::numeric_limits<size_type>::max()) if(subMatrixIndex[col.index()]!=std::numeric_limits<size_type>::max())
// This value is in our submatrix -> copy (use submatrix indices // This value is in our submatrix -> copy (use submatrix indices
initializer.copyValue(col, subMatrixIndex[row.index()], subMatrixIndex [col.index()]); initializer.copyValue(col, subMatrixIndex[row.index()], subMatrixIndex [col.index()]);
} }
initializer.createMatrix(); initializer.createMatrix();
} }
#ifndef DOXYGEN #ifndef DOXYGEN
template<class B, class TA, int n, int m> template<class Mat, class I>
ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::ColCompMatrix() ColCompMatrix<Mat, I>::ColCompMatrix()
: N_(0), M_(0), Nnz_(0), values(0), rowindex(0), colstart(0) : N_(0), M_(0), Nnz_(0), values(0), rowindex(0), colstart(0)
{} {}
template<class B, class TA, int n, int m> template<class Mat, class I>
ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> > ColCompMatrix<Mat, I>
::ColCompMatrix(const Matrix& mat) ::ColCompMatrix(const Matrix& mat)
: N_(n*mat.N()), M_(m*mat.M()), Nnz_(n*m*mat.nonzeroes()) {
{} // WARNING: This assumes that all blocks are dense and identical
typename Matrix::block_type dummy;
const auto n = MatrixDimension<typename Matrix::block_type>::rowdim(dummy);
const auto m = MatrixDimension<typename Matrix::block_type>::coldim(dummy);
N_ = n*mat.N();
M_ = m*mat.N();
Nnz_ = n*m*mat.nonzeroes();
}
template<class B, class TA, int n, int m> template<class Mat, class I>
ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >& ColCompMatrix<Mat, I>&
ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::operator=(const Matrix& mat ColCompMatrix<Mat, I>::operator=(const Matrix& mat)
)
{ {
if(N_+M_+Nnz_!=0) if(N_+M_+Nnz_!=0)
free(); free();
setMatrix(mat); setMatrix(mat);
return *this; return *this;
} }
template<class B, class TA, int n, int m> template<class Mat, class I>
ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >& ColCompMatrix<Mat, I>&
ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::operator=(const ColCompMatr ColCompMatrix<Mat, I>::operator=(const ColCompMatrix& mat)
ix& mat)
{ {
if(N_+M_+Nnz_!=0) if(N_+M_+Nnz_!=0)
free(); free();
N_=mat.N_; N_=mat.N_;
M_=mat.M_; M_=mat.M_;
Nnz_= mat.Nnz_; Nnz_= mat.Nnz_;
if(M_>0) { if(M_>0) {
colstart=new int[M_+1]; colstart=new size_type[M_+1];
for(int i=0; i<=M_; ++i) for(size_type i=0; i<=M_; ++i)
colstart[i]=mat.colstart[i]; colstart[i]=mat.colstart[i];
} }
if(Nnz_>0) { if(Nnz_>0) {
values = new B[Nnz_]; values = new B[Nnz_];
rowindex = new int[Nnz_]; rowindex = new size_type[Nnz_];
for(int i=0; i<Nnz_; ++i) for(size_type i=0; i<Nnz_; ++i)
values[i]=mat.values[i]; values[i]=mat.values[i];
for(int i=0; i<Nnz_; ++i) for(size_type i=0; i<Nnz_; ++i)
rowindex[i]=mat.rowindex[i]; rowindex[i]=mat.rowindex[i];
} }
return *this; return *this;
} }
template<class B, class TA, int n, int m> template<class Mat, class I>
void ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> > void ColCompMatrix<Mat, I>
::setMatrix(const Matrix& mat) ::setMatrix(const Matrix& mat)
{ {
N_=n*mat.N(); N_=MatrixDimension<Mat>::rowdim(mat);
M_=m*mat.M(); M_=MatrixDimension<Mat>::coldim(mat);
ColCompMatrixInitializer<Matrix> initializer(*this); ColCompMatrixInitializer<Mat, I> initializer(*this);
copyToColCompMatrix(initializer, MatrixRowSet<Matrix>(mat)); copyToColCompMatrix(initializer, MatrixRowSet<Matrix>(mat));
} }
template<class B, class TA, int n, int m> template<class Mat, class I>
void ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> > void ColCompMatrix<Mat, I>
::setMatrix(const Matrix& mat, const std::set<std::size_t>& mrs) ::setMatrix(const Matrix& mat, const std::set<std::size_t>& mrs)
{ {
if(N_+M_+Nnz_!=0) if(N_+M_+Nnz_!=0)
free(); free();
N_=mrs.size()*n;
M_=mrs.size()*m;
ColCompMatrixInitializer<Matrix> initializer(*this);
copyToColCompMatrix(initializer, MatrixRowSubset<Matrix,std::set<std::size_t N_=mrs.size()*MatrixDimension<Mat>::rowdim(mat) / mat.N();
> >(mat,mrs)); M_=mrs.size()*MatrixDimension<Mat>::coldim(mat) / mat.M();
ColCompMatrixInitializer<Mat, I> initializer(*this);
copyToColCompMatrix(initializer, MatrixRowSubset<Mat,std::set<std::size_t> >
(mat,mrs));
} }
template<class B, class TA, int n, int m> template<class Mat, class I>
ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::~ColCompMatrix() ColCompMatrix<Mat, I>::~ColCompMatrix()
{ {
if(N_+M_+Nnz_!=0) if(N_+M_+Nnz_!=0)
free(); free();
} }
template<class B, class TA, int n, int m> template<class Mat, class I>
void ColCompMatrix<BCRSMatrix<FieldMatrix<B,n,m>,TA> >::free() void ColCompMatrix<Mat, I>::free()
{ {
delete[] values; delete[] values;
delete[] rowindex; delete[] rowindex;
delete[] colstart; delete[] colstart;
N_ = 0; N_ = 0;
M_ = 0; M_ = 0;
Nnz_ = 0; Nnz_ = 0;
} }
#endif // DOXYGEN #endif // DOXYGEN
 End of changes. 63 change blocks. 
135 lines changed or deleted 186 lines changed or added

Home  |  About  |  Features  |  All  |  Newest  |  Dox  |  Diffs  |  RSS Feeds  |  Screenshots  |  Comments  |  Imprint  |  Privacy  |  HTTP(S)