34 , chunk_size(s.chunk_size)
35 , sparsity_pattern(s.sparsity_pattern)
39 "This constructor can only be called if the provided argument "
40 "is the sparsity pattern for an empty matrix. This constructor can "
41 "not be used to copy-construct a non-empty sparsity pattern."));
63 const std::vector<size_type> &row_lengths,
84 const std::vector<size_type> &row_lengths,
99 "This operator can only be called if the provided argument "
100 "is the sparsity pattern for an empty matrix. This operator can "
101 "not be used to copy a non-empty sparsity pattern."));
104 ExcMessage(
"This operator can only be called if the current object is "
124 const std::vector<size_type> row_lengths(m, max_per_row);
157 std::vector<unsigned int> chunk_row_lengths(m_chunks, 0);
159 chunk_row_lengths[i /
chunk_size] += row_lengths[i];
164 if (m != n && m_chunks == n_chunks)
165 for (
unsigned int i = 0; i < m_chunks; ++i)
166 ++chunk_row_lengths[i];
181template <
typename SparsityPatternType>
205 for (
size_type row = 0; row < dsp.n_rows(); ++row)
212 for (
typename SparsityPatternType::iterator col_num = dsp.begin(row);
213 col_num != dsp.end(row);
215 temporary_sp.
add(reduced_row, col_num->column() /
chunk_size);
223template <
typename number>
234 std::vector<size_type> entries_per_row(matrix.m(), 0);
235 for (
size_type row = 0; row < matrix.m(); ++row)
237 for (
size_type col = 0; col < matrix.n(); ++col)
238 if (matrix(row, col) != 0)
239 ++entries_per_row[row];
241 if ((matrix.m() == matrix.n()) && (matrix(row, row) == 0))
242 ++entries_per_row[row];
248 for (
size_type row = 0; row < matrix.m(); ++row)
249 for (
size_type col = 0; col < matrix.n(); ++col)
250 if (matrix(row, col) != 0)
262 const std::vector<size_type> &row_lengths,
276 template <
typename SparsityPatternType>
293template <
typename Sparsity>
297 const Sparsity &sparsity_pattern_for_chunks,
301 Assert(m > (sparsity_pattern_for_chunks.n_rows() - 1) * chunk_size_in &&
302 m <= sparsity_pattern_for_chunks.n_rows() * chunk_size_in,
303 ExcMessage(
"Number of rows m is not compatible with chunk size "
304 "and number of rows in sparsity pattern for the chunks."));
305 Assert(n > (sparsity_pattern_for_chunks.n_cols() - 1) * chunk_size_in &&
306 n <= sparsity_pattern_for_chunks.n_cols() * chunk_size_in,
308 "Number of columns m is not compatible with chunk size "
309 "and number of columns in sparsity pattern for the chunks."));
386 for (; p !=
end; ++p)
477 out <<
']' << std::endl;
512 << -
static_cast<signed int>(i *
chunk_size + e) << std::endl;
598ChunkSparsityPattern::copy_from<DynamicSparsityPattern>(
602ChunkSparsityPattern::create_from<SparsityPattern>(
const size_type,
608ChunkSparsityPattern::create_from<DynamicSparsityPattern>(
ArrayView< typename std::remove_reference< typename std::iterator_traits< Iterator >::reference >::type, MemorySpaceType > make_array_view(const Iterator begin, const Iterator end)
void create_from(const size_type m, const size_type n, const Sparsity &sparsity_pattern_for_chunks, const size_type chunk_size, const bool optimize_diagonal=true)
void add(const size_type i, const size_type j)
void block_write(std::ostream &out) const
bool stores_only_added_elements() const
std::size_t memory_consumption() const
bool exists(const size_type i, const size_type j) const
void reinit(const size_type m, const size_type n, const size_type max_per_row, const size_type chunk_size)
void print_gnuplot(std::ostream &out) const
void block_read(std::istream &in)
ChunkSparsityPattern & operator=(const ChunkSparsityPattern &)
size_type max_entries_per_row() const
size_type bandwidth() const
size_type n_nonzero_elements() const
SparsityPattern sparsity_pattern
size_type row_length(const size_type row) const
void print(std::ostream &out) const
void copy_from(const size_type n_rows, const size_type n_cols, const ForwardIterator begin, const ForwardIterator end, const size_type chunk_size)
void add(const size_type i, const size_type j)
void block_write(std::ostream &out) const
void reinit(const size_type m, const size_type n, const ArrayView< const unsigned int > &row_lengths)
bool stores_only_added_elements() const
size_type bandwidth() const
std::size_t n_nonzero_elements() const
bool exists(const size_type i, const size_type j) const
std::unique_ptr< size_type[]> colnums
void copy_from(const size_type n_rows, const size_type n_cols, const ForwardIterator begin, const ForwardIterator end)
static constexpr size_type invalid_entry
std::unique_ptr< std::size_t[]> rowstart
void add(const size_type i, const size_type j)
size_type max_entries_per_row() const
unsigned int row_length(const size_type row) const
std::size_t memory_consumption() const
void block_read(std::istream &in)
#define DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_CLOSE
static ::ExceptionBase & ExcIO()
#define Assert(cond, exc)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcEmptyObject()
static ::ExceptionBase & ExcNotQuadratic()
static ::ExceptionBase & ExcInvalidIndex(size_type arg1, size_type arg2)
static ::ExceptionBase & ExcMessage(std::string arg1)
#define AssertThrow(cond, exc)
static ::ExceptionBase & ExcInvalidNumber(size_type arg1)
::VectorizedArray< Number, width > min(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)