Commit 503f606b authored by akbib's avatar akbib
Browse files

Merge branch 'feature/add_blockmatrixview' into 'master'

Add class to merge a vector of matrices into a big combined matrix.

This class is copied from dune-fufem to here

See merge request !1
parents 1374e3e2 389a9d44
......@@ -5,6 +5,7 @@ install(FILES
addtodiagonal.hh
axpy.hh
axy.hh
blockmatrixview.hh
componentwisematrixmap.hh
crossproduct.hh
genericvectortools.hh
......
// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set ts=4 sw=2 et sts=2:
#ifndef DUNE_MATRIX_VECTOR_BLOCK_MATRIX_VIEW_HH
#define DUNE_MATRIX_VECTOR_BLOCK_MATRIX_VIEW_HH
#include <array>
#include <vector>
#include <dune/istl/matrixindexset.hh>
namespace Dune {
namespace MatrixVector {
/** \brief Class that provides global indices for a block matrix given local indices of the blocks. Furthermore, it provides a method to construct a block matrix from given sub matrices.
*
* This class provides a method to combine matrices \f A_1,\ldots,A_n\f into one single block matrix
* \f B = \begin{pmatrix} A_1 & 0 & \ldots & 0 \\
* 0 & A_2 & \ldots & 0 \\
* \ldots \\
* 0 & \ldots & 0 & A_n \end{pmatrix} \f,
* by copying the corresponding entries.
* A BlockMatrixView object provides methods that return global indices for \f B\f given a local index of \f A_i\f and the block index \f i\f.
*
* \tparam MatrixType The matrix type of the matrix blocks.
*/
template <class MatrixType>
class BlockMatrixView
{
public:
//! Construct block view from vector of matrices corr. to the blocks
BlockMatrixView(const std::vector<const MatrixType*>& submat) {
offsets_.resize(submat.size());
offsets_[0] = {0, 0};
for (size_t i=0; i<submat.size()-1; i++)
offsets_[i+1] = {offsets_[i][0] + submat[i]->N(),
offsets_[i][1] + submat[i]->M()};
nRows_ = offsets_.back()[0] + submat.back()->N();
nCols_ = offsets_.back()[1] + submat.back()->M();
}
//! Get global row index for given block and a local row index
size_t row(size_t block, size_t index) const {
return offsets_[block][0] + index;
}
//! Get global column index given a block and a local column index
size_t col(size_t block, size_t index) const {
return offsets_[block][1] + index;
}
//! Return the total number of rows
size_t nRows() const {
return nRows_;
}
//! Return the total number of columns
size_t nCols() const {
return nCols_;
}
/** \brief Combine submatrices in a big matrix. */
static void setupBlockMatrix(const std::vector<MatrixType>& submat, MatrixType& mat) {
std::vector<const MatrixType*> dummy(submat.size());
for (size_t i=0; i<dummy.size(); i++)
dummy[i] = &submat[i];
setupBlockMatrix(dummy,mat);
}
/** \brief Combine submatrices in a big matrix. */
static void setupBlockMatrix(const std::vector<const MatrixType*>& submat, MatrixType& mat)
{
BlockMatrixView<MatrixType> blockView(submat);
MatrixIndexSet indexSet(blockView.nRows(), blockView.nCols());
for (size_t i=0; i<submat.size(); i++)
indexSet.import(*submat[i], blockView.row(i,0), blockView.col(i,0));
indexSet.exportIdx(mat);
for (size_t i=0; i<submat.size(); i++) {
auto row = submat[i]->begin();
auto rowEnd = submat[i]->end();
for (; row != rowEnd; row++) {
auto col = row->begin();
auto colEnd = row->end();
for (; col != colEnd; col++)
mat[blockView.row(i,row.index())][blockView.col(i,col.index())] = *col;
}
}
}
protected:
//! Offsets for the different blocks
std::vector<std::array<size_t, 2> > offsets_;
//! Total numbers of rows and columns
size_t nRows_, nCols_;
};
} // end namespace MatrixVector
} // end namespace Dune
#endif
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment