3#ifndef DUNE_ISTL_BCCSMATRIX_INITIALIZER_HH
4#define DUNE_ISTL_BCCSMATRIX_INITIALIZER_HH
9#include <dune/common/typetraits.hh>
10#include <dune/common/scalarmatrixview.hh>
16 template<
class I,
class S,
class D>
17 class OverlappingSchwarzInitializer;
20namespace Dune::ISTL::Impl
29 template<
class M,
class S>
36 typedef S RowIndexSet;
43 MatrixRowSubset(
const Matrix& m,
const RowIndexSet& s)
47 const Matrix& matrix()
const
52 const RowIndexSet& rowIndexSet()
const
59 :
public ForwardIteratorFacade<const_iterator, const typename Matrix::row_type>
62 const_iterator(
typename Matrix::const_iterator firstRow,
63 typename RowIndexSet::const_iterator pos)
64 : firstRow_(firstRow), pos_(pos)
70 return *(firstRow_+ *pos_);
72 bool equals(
const const_iterator& o)
const
80 typename RowIndexSet::value_type index()
const
87 typename Matrix::const_iterator firstRow_;
89 typename RowIndexSet::const_iterator pos_;
93 const_iterator begin()
const
95 return const_iterator(m_.begin(), s_.begin());
98 const_iterator end()
const
100 return const_iterator(m_.begin(), s_.end());
107 const RowIndexSet& s_;
116 template<
class M,
class I =
typename M::
size_type>
117 class BCCSMatrixInitializer
119 template<
class IList,
class S,
class D>
124 typedef Dune::ISTL::Impl::BCCSMatrix<typename Matrix::field_type, I> OutputMatrix;
129 BCCSMatrixInitializer(OutputMatrix& mat_)
130 :
mat(&mat_), cols(mat_.M())
132 if constexpr (Dune::IsNumber<typename M::block_type>::value)
139 n = M::block_type::rows;
140 m = M::block_type::cols;
146 BCCSMatrixInitializer()
147 :
mat(0), cols(0), n(0), m(0)
150 virtual ~BCCSMatrixInitializer()
153 template<
typename Iter>
156 mat->Nnz_+=row->getsize();
159 template<
typename Iter,
typename FullMatrixIndex>
160 void addRowNnz(
const Iter& row,
const std::set<FullMatrixIndex>& indices)
const
162 auto siter =indices.begin();
163 for (
auto entry=row->begin(); entry!=row->end(); ++entry)
165 for(; siter!=indices.end() && *siter<entry.index(); ++siter) ;
166 if(siter==indices.end())
168 if(*siter==entry.index())
174 template<
typename Iter,
typename SubMatrixIndex>
175 void addRowNnz(
const Iter& row,
const std::vector<SubMatrixIndex>& indices)
const
177 for (
auto entry=row->begin(); entry!=row->end(); ++entry)
178 if (indices[entry.index()]!=std::numeric_limits<SubMatrixIndex>::max())
184 allocateMatrixStorage();
188 template<
typename Iter,
typename CIter>
198 assert(colindex*m+i<cols);
199 marker[colindex*m+i]+=n;
208 mat->colstart[i+1]=
mat->colstart[i]+marker[i];
209 marker[i]=
mat->colstart[i];
213 template<
typename Iter,
typename CIter>
219 template<
typename CIter>
224 assert(colindex*m+j<cols-1 || (
size_type)marker[colindex*m+j]<(
size_type)
mat->colstart[colindex*m+j+1]);
226 mat->rowindex[marker[colindex*m+j]]=rowindex*n+i;
227 mat->values[marker[colindex*m+j]] = Dune::Impl::asMatrix(*
col)[i][j];
228 ++marker[colindex*m+j];
240 void allocateMatrixStorage()
const
244 mat->values=
new typename M::field_type[
mat->Nnz_];
245 mat->rowindex=
new I[
mat->Nnz_];
246 mat->colstart=
new I[cols+1];
249 void allocateMarker()
252 std::fill(marker.begin(), marker.end(), 0);
262 mutable std::vector<size_type> marker;
265 template<
class F,
class Matrix>
266 void copyToBCCSMatrix(F& initializer,
const Matrix& matrix)
268 for (
auto row=matrix.begin(); row!= matrix.end(); ++row)
269 initializer.addRowNnz(row);
271 initializer.allocate();
273 for (
auto row=matrix.begin(); row!= matrix.end(); ++row) {
275 for (
auto col=row->begin();
col != row->end(); ++
col)
276 initializer.countEntries(row,
col);
279 initializer.calcColstart();
281 for (
auto row=matrix.begin(); row!= matrix.end(); ++row) {
282 for (
auto col=row->begin();
col != row->end(); ++
col) {
283 initializer.copyValue(row,
col);
287 initializer.createMatrix();
290 template<
class F,
class M,
class S>
291 void copyToBCCSMatrix(F& initializer,
const MatrixRowSubset<M,S>& mrs)
293 typedef MatrixRowSubset<M,S> MRS;
294 typedef typename MRS::RowIndexSet SIS;
295 typedef typename SIS::const_iterator SIter;
296 typedef typename MRS::const_iterator
Iter;
297 typedef typename std::iterator_traits<Iter>::value_type row_type;
298 typedef typename row_type::const_iterator
CIter;
300 typedef typename MRS::Matrix::size_type
size_type;
306 std::vector<size_type> subMatrixIndex(mrs.matrix().N(),
307 std::numeric_limits<size_type>::max());
309 for(SIter index = mrs.rowIndexSet().begin(); index!=mrs.rowIndexSet().end(); ++index)
310 subMatrixIndex[*index]=s++;
313 for(
Iter row=mrs.begin(); row!= mrs.end(); ++row)
314 initializer.addRowNnz(row, subMatrixIndex);
316 initializer.allocate();
318 for(
Iter row=mrs.begin(); row!= mrs.end(); ++row)
320 if(subMatrixIndex[
col.index()]!=std::numeric_limits<size_type>::max())
322 initializer.countEntries(subMatrixIndex[
col.index()]);
325 initializer.calcColstart();
327 for(
Iter row=mrs.begin(); row!= mrs.end(); ++row)
329 if(subMatrixIndex[
col.index()]!=std::numeric_limits<size_type>::max())
331 initializer.copyValue(
col, subMatrixIndex[row.index()], subMatrixIndex[
col.index()]);
333 initializer.createMatrix();
Col col
Definition: matrixmatrix.hh:349
Matrix & mat
Definition: matrixmatrix.hh:345
void addRowNnz(const Iter &row)
Definition: overlappingschwarz.hh:893
void calcColstart() const
Definition: overlappingschwarz.hh:924
void copyValue(const Iter &row, const CIter &col) const
Definition: overlappingschwarz.hh:931
void createMatrix() const
Definition: overlappingschwarz.hh:945
void allocate()
Definition: overlappingschwarz.hh:903
std::size_t countEntries(const BlockVector< T, A > &vector)
Definition: matrixmarket.hh:1060
Definition: allocator.hh:9
Initializer for SuperLU Matrices representing the subdomains.
Definition: overlappingschwarz.hh:45
Matrix::row_type::const_iterator CIter
Definition: overlappingschwarz.hh:54
Matrix::const_iterator Iter
Definition: overlappingschwarz.hh:53
IndexSet::size_type size_type
Definition: overlappingschwarz.hh:57
AtomInitializer::Matrix Matrix
Definition: overlappingschwarz.hh:52
A::size_type size_type
Type for indices and sizes.
Definition: matrix.hh:575
MatrixImp::DenseMatrixBase< T, A >::window_type row_type
The type implementing a matrix row.
Definition: matrix.hh:572