Reference documentation for deal.II version Git f102a78320 2021-10-16 14:49:02 -0400
\(\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\}}\)
trilinos_block_sparse_matrix.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2008 - 2020 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 #ifndef dealii_trilinos_block_sparse_matrix_h
17 #define dealii_trilinos_block_sparse_matrix_h
18 
19 
20 #include <deal.II/base/config.h>
21 
22 #ifdef DEAL_II_WITH_TRILINOS
23 
25 
27 # include <deal.II/lac/exceptions.h>
28 # include <deal.II/lac/full_matrix.h>
31 
32 # include <cmath>
33 
35 
36 // forward declarations
37 # ifndef DOXYGEN
39 template <typename number>
40 class BlockSparseMatrix;
41 # endif
42 
43 namespace TrilinosWrappers
44 {
71  class BlockSparseMatrix : public BlockMatrixBase<SparseMatrix>
72  {
73  public:
78 
83 
95 
107  BlockSparseMatrix() = default;
108 
112  ~BlockSparseMatrix() override;
113 
119  operator=(const BlockSparseMatrix &) = default;
120 
131  operator=(const double d);
132 
146  void
147  reinit(const size_type n_block_rows, const size_type n_block_columns);
148 
154  template <typename BlockSparsityPatternType>
155  void
156  reinit(const std::vector<IndexSet> & input_maps,
157  const BlockSparsityPatternType &block_sparsity_pattern,
158  const MPI_Comm & communicator = MPI_COMM_WORLD,
159  const bool exchange_data = false);
160 
166  template <typename BlockSparsityPatternType>
167  void
168  reinit(const BlockSparsityPatternType &block_sparsity_pattern);
169 
176  void
177  reinit(
178  const std::vector<IndexSet> & parallel_partitioning,
179  const ::BlockSparseMatrix<double> &dealii_block_sparse_matrix,
180  const MPI_Comm & communicator = MPI_COMM_WORLD,
181  const double drop_tolerance = 1e-13);
182 
190  void
191  reinit(const ::BlockSparseMatrix<double> &deal_ii_sparse_matrix,
192  const double drop_tolerance = 1e-13);
193 
201  bool
202  is_compressed() const;
203 
214  void
215  collect_sizes();
216 
221  size_type
222  n_nonzero_elements() const;
223 
227  MPI_Comm
228  get_mpi_communicator() const;
229 
235  std::vector<IndexSet>
237 
243  std::vector<IndexSet>
245 
252  template <typename VectorType1, typename VectorType2>
253  void
254  vmult(VectorType1 &dst, const VectorType2 &src) const;
255 
261  template <typename VectorType1, typename VectorType2>
262  void
263  Tvmult(VectorType1 &dst, const VectorType2 &src) const;
264 
279  const MPI::BlockVector &x,
280  const MPI::BlockVector &b) const;
281 
291  const MPI::Vector & x,
292  const MPI::BlockVector &b) const;
293 
302  residual(MPI::Vector & dst,
303  const MPI::BlockVector &x,
304  const MPI::Vector & b) const;
305 
314  residual(MPI::Vector & dst,
315  const MPI::Vector &x,
316  const MPI::Vector &b) const;
317 
323 
333  int,
334  int,
335  int,
336  int,
337  << "The blocks [" << arg1 << ',' << arg2 << "] and [" << arg3
338  << ',' << arg4 << "] have differing row numbers.");
339 
344  int,
345  int,
346  int,
347  int,
348  << "The blocks [" << arg1 << ',' << arg2 << "] and [" << arg3
349  << ',' << arg4 << "] have differing column numbers.");
351 
352  private:
356  template <typename VectorType1, typename VectorType2>
357  void
358  vmult(VectorType1 & dst,
359  const VectorType2 &src,
360  const bool transpose,
361  const std::integral_constant<bool, true>,
362  const std::integral_constant<bool, true>) const;
363 
368  template <typename VectorType1, typename VectorType2>
369  void
370  vmult(VectorType1 & dst,
371  const VectorType2 &src,
372  const bool transpose,
373  const std::integral_constant<bool, false>,
374  const std::integral_constant<bool, true>) const;
375 
380  template <typename VectorType1, typename VectorType2>
381  void
382  vmult(VectorType1 & dst,
383  const VectorType2 &src,
384  const bool transpose,
385  const std::integral_constant<bool, true>,
386  const std::integral_constant<bool, false>) const;
387 
393  template <typename VectorType1, typename VectorType2>
394  void
395  vmult(VectorType1 & dst,
396  const VectorType2 &src,
397  const bool transpose,
398  const std::integral_constant<bool, false>,
399  const std::integral_constant<bool, false>) const;
400  };
401 
402 
403 
406  // ------------- inline and template functions -----------------
407 
408 
409 
410  inline BlockSparseMatrix &
412  {
414 
415  for (size_type r = 0; r < this->n_block_rows(); ++r)
416  for (size_type c = 0; c < this->n_block_cols(); ++c)
417  this->block(r, c) = d;
418 
419  return *this;
420  }
421 
422 
423 
424  inline bool
426  {
427  bool compressed = true;
428  for (size_type row = 0; row < n_block_rows(); ++row)
429  for (size_type col = 0; col < n_block_cols(); ++col)
430  if (block(row, col).is_compressed() == false)
431  {
432  compressed = false;
433  break;
434  }
435 
436  return compressed;
437  }
438 
439 
440 
441  template <typename VectorType1, typename VectorType2>
442  inline void
443  BlockSparseMatrix::vmult(VectorType1 &dst, const VectorType2 &src) const
444  {
445  vmult(dst,
446  src,
447  false,
448  std::integral_constant<bool, IsBlockVector<VectorType1>::value>(),
449  std::integral_constant<bool, IsBlockVector<VectorType2>::value>());
450  }
451 
452 
453 
454  template <typename VectorType1, typename VectorType2>
455  inline void
456  BlockSparseMatrix::Tvmult(VectorType1 &dst, const VectorType2 &src) const
457  {
458  vmult(dst,
459  src,
460  true,
461  std::integral_constant<bool, IsBlockVector<VectorType1>::value>(),
462  std::integral_constant<bool, IsBlockVector<VectorType2>::value>());
463  }
464 
465 
466 
467  template <typename VectorType1, typename VectorType2>
468  inline void
469  BlockSparseMatrix::vmult(VectorType1 & dst,
470  const VectorType2 &src,
471  const bool transpose,
472  std::integral_constant<bool, true>,
473  std::integral_constant<bool, true>) const
474  {
475  if (transpose == true)
477  else
479  }
480 
481 
482 
483  template <typename VectorType1, typename VectorType2>
484  inline void
485  BlockSparseMatrix::vmult(VectorType1 & dst,
486  const VectorType2 &src,
487  const bool transpose,
488  std::integral_constant<bool, false>,
489  std::integral_constant<bool, true>) const
490  {
491  if (transpose == true)
493  else
495  }
496 
497 
498 
499  template <typename VectorType1, typename VectorType2>
500  inline void
501  BlockSparseMatrix::vmult(VectorType1 & dst,
502  const VectorType2 &src,
503  const bool transpose,
504  std::integral_constant<bool, true>,
505  std::integral_constant<bool, false>) const
506  {
507  if (transpose == true)
509  else
511  }
512 
513 
514 
515  template <typename VectorType1, typename VectorType2>
516  inline void
517  BlockSparseMatrix::vmult(VectorType1 & dst,
518  const VectorType2 &src,
519  const bool transpose,
520  std::integral_constant<bool, false>,
521  std::integral_constant<bool, false>) const
522  {
523  if (transpose == true)
525  else
527  }
528 
529 
530 
531  inline std::vector<IndexSet>
533  {
534  Assert(this->n_block_cols() != 0, ExcNotInitialized());
535  Assert(this->n_block_rows() != 0, ExcNotInitialized());
536 
537  std::vector<IndexSet> domain_indices;
538  for (size_type c = 0; c < this->n_block_cols(); ++c)
539  domain_indices.push_back(
540  this->sub_objects[0][c]->locally_owned_domain_indices());
541 
542  return domain_indices;
543  }
544 
545 
546 
547  inline std::vector<IndexSet>
549  {
550  Assert(this->n_block_cols() != 0, ExcNotInitialized());
551  Assert(this->n_block_rows() != 0, ExcNotInitialized());
552 
553  std::vector<IndexSet> range_indices;
554  for (size_type r = 0; r < this->n_block_rows(); ++r)
555  range_indices.push_back(
556  this->sub_objects[r][0]->locally_owned_range_indices());
557 
558  return range_indices;
559  }
560 
561 
562 
563  namespace internal
564  {
565  namespace BlockLinearOperatorImplementation
566  {
581  template <typename PayloadBlockType>
583  {
584  public:
588  using BlockType = PayloadBlockType;
589 
599  template <typename... Args>
600  TrilinosBlockPayload(const Args &...)
601  {
602  static_assert(
603  std::is_same<
604  PayloadBlockType,
606  "TrilinosBlockPayload can only accept a payload of type TrilinosPayload.");
607  }
608  };
609 
610  } // namespace BlockLinearOperatorImplementation
611  } /* namespace internal */
612 
613 
614 } /* namespace TrilinosWrappers */
615 
616 
618 
619 #endif // DEAL_II_WITH_TRILINOS
620 
621 #endif // dealii_trilinos_block_sparse_matrix_h
void vmult_block_nonblock(BlockVectorType &dst, const VectorType &src) const
static ::ExceptionBase & ExcIncompatibleRowNumbers(int arg1, int arg2, int arg3, int arg4)
types::global_dof_index size_type
SymmetricTensor< 2, dim, Number > e(const Tensor< 2, dim, Number > &F)
static ::ExceptionBase & ExcScalarAssignmentOnlyForZeroValue()
void vmult_nonblock_block(VectorType &dst, const BlockVectorType &src) const
MatrixIterator< BlockMatrixIterators::Accessor< BlockMatrixBase, true > > const_iterator
static ::ExceptionBase & ExcNotInitialized()
void Tvmult_block_nonblock(BlockVectorType &dst, const VectorType &src) const
void Tvmult_block_block(BlockVectorType &dst, const BlockVectorType &src) const
void Tvmult_nonblock_block(VectorType &dst, const BlockVectorType &src) const
void vmult_block_block(BlockVectorType &dst, const BlockVectorType &src) const
typename BlockType::value_type value_type
void Tvmult_nonblock_nonblock(VectorType &dst, const VectorType &src) const
static ::ExceptionBase & ExcIncompatibleColNumbers(int arg1, int arg2, int arg3, int arg4)
unsigned int n_block_cols() const
void vmult(VectorType1 &dst, const VectorType2 &src) const
#define Assert(cond, exc)
Definition: exceptions.h:1461
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
void reinit(const size_type n_block_rows, const size_type n_block_columns)
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:401
SymmetricTensor< 2, dim, Number > d(const Tensor< 2, dim, Number > &F, const Tensor< 2, dim, Number > &dF_dt)
SymmetricTensor< 2, dim, Number > b(const Tensor< 2, dim, Number > &F)
void Tvmult(VectorType1 &dst, const VectorType2 &src) const
MatrixIterator< BlockMatrixIterators::Accessor< BlockMatrixBase, false > > iterator
void vmult_nonblock_nonblock(VectorType &dst, const VectorType &src) const
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:400
BlockType & block(const unsigned int row, const unsigned int column)
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
Definition: exceptions.h:578
unsigned int n_block_rows() const
std::vector< IndexSet > locally_owned_range_indices() const
TrilinosScalar residual(MPI::BlockVector &dst, const MPI::BlockVector &x, const MPI::BlockVector &b) const
BlockSparseMatrix & operator=(const BlockSparseMatrix &)=default
std::vector< IndexSet > locally_owned_domain_indices() const