Reference documentation for deal.II version GIT relicensing-214-g6e74dec06b 2024-03-27 18:10:01+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
Loading...
Searching...
No Matches
block_sparse_matrix_ez.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2002 - 2023 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_block_sparse_matrix_ez_h
16#define dealii_block_sparse_matrix_ez_h
17
18
19// TODO: Derive BlockSparseMatrixEZ from BlockMatrixBase, like all the
20// other block matrices as well; this would allow to instantiate a few
21// functions with this template argument as well (in particular
22// AffineConstraints::distribute_local_to_global)
23
24#include <deal.II/base/config.h>
25
29#include <deal.II/base/table.h>
30
33
35
36// Forward declaration
37#ifndef DOXYGEN
38template <typename Number>
39class BlockVector;
40#endif
41
59template <typename Number>
61{
62public:
67
72
77 BlockSparseMatrixEZ(const unsigned int block_rows,
78 const unsigned int block_cols);
79
86
93
103 operator=(const double d);
104
105
109 void
111
119 void
120 reinit(const unsigned int n_block_rows, const unsigned int n_block_cols);
127 void
129
134 block(const unsigned int row, const unsigned int column);
135
136
142 block(const unsigned int row, const unsigned int column) const;
143
147 unsigned int
148 n_block_rows() const;
149
153 unsigned int
154 n_block_cols() const;
155
162 bool
163 empty() const;
164
172 m() const;
173
181 n() const;
182
188 void
189 set(const size_type i, const size_type j, const Number value);
190
196 void
197 add(const size_type i, const size_type j, const Number value);
198
199
204 template <typename somenumber>
205 void
207
213 template <typename somenumber>
214 void
216 const BlockVector<somenumber> &src) const;
217
222 template <typename somenumber>
223 void
225 const BlockVector<somenumber> &src) const;
226
232 template <typename somenumber>
233 void
235 const BlockVector<somenumber> &src) const;
236
237
243 template <typename StreamType>
244 void
245 print_statistics(StreamType &s, bool full = false);
246
247private:
253
259
264};
265
267/*----------------------------------------------------------------------*/
268
269
270template <typename Number>
271inline unsigned int
273{
274 return row_indices.size();
275}
276
277
278
279template <typename Number>
280inline unsigned int
282{
283 return column_indices.size();
284}
285
286
287
288template <typename Number>
291 const unsigned int column)
292{
293 AssertIndexRange(row, n_block_rows());
294 AssertIndexRange(column, n_block_cols());
295
296 return blocks[row][column];
297}
298
299
300
301template <typename Number>
302inline const SparseMatrixEZ<Number> &
304 const unsigned int column) const
305{
306 AssertIndexRange(row, n_block_rows());
307 AssertIndexRange(column, n_block_cols());
308
309 return blocks[row][column];
310}
311
312
313
314template <typename Number>
317{
318 return row_indices.total_size();
319}
320
321
322
323template <typename Number>
326{
327 return column_indices.total_size();
328}
329
330
331
332template <typename Number>
333inline void
335 const size_type j,
336 const Number value)
337{
338 AssertIsFinite(value);
339
340 const std::pair<size_type, size_type> row_index =
341 row_indices.global_to_local(i),
342 col_index =
343 column_indices.global_to_local(j);
344 block(row_index.first, col_index.first)
345 .set(row_index.second, col_index.second, value);
346}
347
348
349
350template <typename Number>
351inline void
353 const size_type j,
354 const Number value)
355{
356 AssertIsFinite(value);
357
358 const std::pair<unsigned int, size_type> row_index =
359 row_indices.global_to_local(i),
360 col_index =
361 column_indices.global_to_local(j);
362 block(row_index.first, col_index.first)
363 .add(row_index.second, col_index.second, value);
364}
365
366
367template <typename Number>
368template <typename somenumber>
369void
371 const BlockVector<somenumber> &src) const
372{
373 Assert(dst.n_blocks() == n_block_rows(),
374 ExcDimensionMismatch(dst.n_blocks(), n_block_rows()));
375 Assert(src.n_blocks() == n_block_cols(),
376 ExcDimensionMismatch(src.n_blocks(), n_block_cols()));
377
378 dst = 0.;
379
380 for (unsigned int row = 0; row < n_block_rows(); ++row)
381 for (unsigned int col = 0; col < n_block_cols(); ++col)
382 block(row, col).vmult_add(dst.block(row), src.block(col));
383}
384
385
386
387template <typename Number>
388template <typename somenumber>
389void
391 const BlockVector<somenumber> &src) const
392{
393 Assert(dst.n_blocks() == n_block_rows(),
394 ExcDimensionMismatch(dst.n_blocks(), n_block_rows()));
395 Assert(src.n_blocks() == n_block_cols(),
396 ExcDimensionMismatch(src.n_blocks(), n_block_cols()));
397
398 for (unsigned int row = 0; row < n_block_rows(); ++row)
399 for (unsigned int col = 0; col < n_block_cols(); ++col)
400 block(row, col).vmult_add(dst.block(row), src.block(col));
401}
402
403
404
405template <typename Number>
406template <typename somenumber>
407void
409 const BlockVector<somenumber> &src) const
410{
411 Assert(dst.n_blocks() == n_block_cols(),
412 ExcDimensionMismatch(dst.n_blocks(), n_block_cols()));
413 Assert(src.n_blocks() == n_block_rows(),
414 ExcDimensionMismatch(src.n_blocks(), n_block_rows()));
415
416 dst = 0.;
417
418 for (unsigned int row = 0; row < n_block_rows(); ++row)
419 for (unsigned int col = 0; col < n_block_cols(); ++col)
420 block(row, col).Tvmult_add(dst.block(col), src.block(row));
421}
422
423
424
425template <typename Number>
426template <typename somenumber>
427void
430 const BlockVector<somenumber> &src) const
431{
432 Assert(dst.n_blocks() == n_block_cols(),
433 ExcDimensionMismatch(dst.n_blocks(), n_block_cols()));
434 Assert(src.n_blocks() == n_block_rows(),
435 ExcDimensionMismatch(src.n_blocks(), n_block_rows()));
436
437 for (unsigned int row = 0; row < n_block_rows(); ++row)
438 for (unsigned int col = 0; col < n_block_cols(); ++col)
439 block(row, col).Tvmult_add(dst.block(col), src.block(row));
440}
441
442
443template <typename number>
444template <typename StreamType>
445inline void
447{
448 size_type used_total = 0;
449 size_type allocated_total = 0;
450 size_type reserved_total = 0;
451 std::vector<size_type> used_by_line_total;
452
453 size_type used;
454 size_type allocated;
455 size_type reserved;
456 std::vector<size_type> used_by_line;
457
458 for (size_type i = 0; i < n_block_rows(); ++i)
459 for (size_type j = 0; j < n_block_cols(); ++j)
460 {
461 used_by_line.clear();
462 out << "block:\t" << i << '\t' << j << std::endl;
463 block(i, j).compute_statistics(
464 used, allocated, reserved, used_by_line, full);
465
466 out << "used:" << used << std::endl
467 << "allocated:" << allocated << std::endl
468 << "reserved:" << reserved << std::endl;
469
470 used_total += used;
471 allocated_total += allocated;
472 reserved_total += reserved;
473
474 if (full)
475 {
476 used_by_line_total.resize(used_by_line.size());
477 for (size_type i = 0; i < used_by_line.size(); ++i)
478 if (used_by_line[i] != 0)
479 {
480 out << "row-entries\t" << i << "\trows\t" << used_by_line[i]
481 << std::endl;
482 used_by_line_total[i] += used_by_line[i];
483 }
484 }
485 }
486 out << "Total" << std::endl
487 << "used:" << used_total << std::endl
488 << "allocated:" << allocated_total << std::endl
489 << "reserved:" << reserved_total << std::endl;
490 for (size_type i = 0; i < used_by_line_total.size(); ++i)
491 if (used_by_line_total[i] != 0)
492 {
493 out << "row-entries\t" << i << "\trows\t" << used_by_line_total[i]
494 << std::endl;
495 }
496}
497
498
500
501#endif // dealii_block_sparse_matrix_ez_h
void print_statistics(StreamType &s, bool full=false)
unsigned int n_block_cols() const
unsigned int n_block_rows() const
BlockSparseMatrixEZ()=default
void reinit(const unsigned int n_block_rows, const unsigned int n_block_cols)
void Tvmult_add(BlockVector< somenumber > &dst, const BlockVector< somenumber > &src) const
BlockSparseMatrixEZ & operator=(const BlockSparseMatrixEZ< Number > &)
BlockSparseMatrixEZ & operator=(const double d)
void vmult(BlockVector< somenumber > &dst, const BlockVector< somenumber > &src) const
Table< 2, SparseMatrixEZ< Number > > blocks
void Tvmult(BlockVector< somenumber > &dst, const BlockVector< somenumber > &src) const
SparseMatrixEZ< Number > & block(const unsigned int row, const unsigned int column)
BlockSparseMatrixEZ(const unsigned int block_rows, const unsigned int block_cols)
void set(const size_type i, const size_type j, const Number value)
void add(const size_type i, const size_type j, const Number value)
bool empty() const
void vmult_add(BlockVector< somenumber > &dst, const BlockVector< somenumber > &src) const
BlockSparseMatrixEZ(const BlockSparseMatrixEZ< Number > &)
unsigned int n_blocks() const
BlockType & block(const unsigned int i)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
#define Assert(cond, exc)
#define AssertIsFinite(number)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcDimensionMismatch(std::size_t arg1, std::size_t arg2)
unsigned int global_dof_index
Definition types.h:81