Reference documentation for deal.II version 9.5.0
\(\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
trilinos_block_sparse_matrix.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2008 - 2023 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
31
32# include <cmath>
33
35
36// forward declarations
37# ifndef DOXYGEN
39template <typename number>
41# endif
42
43namespace TrilinosWrappers
44{
72 class BlockSparseMatrix : public BlockMatrixBase<SparseMatrix>
73 {
74 public:
79
84
96
108 BlockSparseMatrix() = default;
109
113 ~BlockSparseMatrix() override;
114
120 operator=(const BlockSparseMatrix &) = default;
121
132 operator=(const double d);
133
147 void
148 reinit(const size_type n_block_rows, const size_type n_block_columns);
149
155 template <typename BlockSparsityPatternType>
156 void
157 reinit(const std::vector<IndexSet> & input_maps,
158 const BlockSparsityPatternType &block_sparsity_pattern,
159 const MPI_Comm communicator = MPI_COMM_WORLD,
160 const bool exchange_data = false);
161
167 template <typename BlockSparsityPatternType>
168 void
169 reinit(const BlockSparsityPatternType &block_sparsity_pattern);
170
177 void
178 reinit(
179 const std::vector<IndexSet> & parallel_partitioning,
180 const ::BlockSparseMatrix<double> &dealii_block_sparse_matrix,
181 const MPI_Comm communicator = MPI_COMM_WORLD,
182 const double drop_tolerance = 1e-13);
183
191 void
192 reinit(const ::BlockSparseMatrix<double> &deal_ii_sparse_matrix,
193 const double drop_tolerance = 1e-13);
194
202 bool
203 is_compressed() const;
204
215 void
217
222 std::uint64_t
223 n_nonzero_elements() const;
224
229 get_mpi_communicator() const;
230
236 std::vector<IndexSet>
238
244 std::vector<IndexSet>
246
253 template <typename VectorType1, typename VectorType2>
254 void
255 vmult(VectorType1 &dst, const VectorType2 &src) const;
256
262 template <typename VectorType1, typename VectorType2>
263 void
264 Tvmult(VectorType1 &dst, const VectorType2 &src) const;
265
280 const MPI::BlockVector &x,
281 const MPI::BlockVector &b) const;
282
292 const MPI::Vector & x,
293 const MPI::BlockVector &b) const;
294
303 residual(MPI::Vector & dst,
304 const MPI::BlockVector &x,
305 const MPI::Vector & b) const;
306
315 residual(MPI::Vector & dst,
316 const MPI::Vector &x,
317 const MPI::Vector &b) const;
318
324
334 int,
335 int,
336 int,
337 int,
338 << "The blocks [" << arg1 << ',' << arg2 << "] and [" << arg3
339 << ',' << arg4 << "] have differing row numbers.");
340
345 int,
346 int,
347 int,
348 int,
349 << "The blocks [" << arg1 << ',' << arg2 << "] and [" << arg3
350 << ',' << arg4 << "] have differing column numbers.");
351 //** @} */
352
353 private:
357 template <typename VectorType1, typename VectorType2>
358 void
359 vmult(VectorType1 & dst,
360 const VectorType2 &src,
361 const bool transpose,
362 const std::integral_constant<bool, true>,
363 const std::integral_constant<bool, true>) const;
364
369 template <typename VectorType1, typename VectorType2>
370 void
371 vmult(VectorType1 & dst,
372 const VectorType2 &src,
373 const bool transpose,
374 const std::integral_constant<bool, false>,
375 const std::integral_constant<bool, true>) const;
376
381 template <typename VectorType1, typename VectorType2>
382 void
383 vmult(VectorType1 & dst,
384 const VectorType2 &src,
385 const bool transpose,
386 const std::integral_constant<bool, true>,
387 const std::integral_constant<bool, false>) const;
388
394 template <typename VectorType1, typename VectorType2>
395 void
396 vmult(VectorType1 & dst,
397 const VectorType2 &src,
398 const bool transpose,
399 const std::integral_constant<bool, false>,
400 const std::integral_constant<bool, false>) const;
401 };
402
403
404
407 // ------------- inline and template functions -----------------
408
409
410
411 inline BlockSparseMatrix &
413 {
415
416 for (size_type r = 0; r < this->n_block_rows(); ++r)
417 for (size_type c = 0; c < this->n_block_cols(); ++c)
418 this->block(r, c) = d;
419
420 return *this;
421 }
422
423
424
425 inline bool
427 {
428 bool compressed = true;
429 for (size_type row = 0; row < n_block_rows(); ++row)
430 for (size_type col = 0; col < n_block_cols(); ++col)
431 if (block(row, col).is_compressed() == false)
432 {
433 compressed = false;
434 break;
435 }
436
437 return compressed;
438 }
439
440
441
442 template <typename VectorType1, typename VectorType2>
443 inline void
444 BlockSparseMatrix::vmult(VectorType1 &dst, const VectorType2 &src) const
445 {
446 vmult(dst,
447 src,
448 false,
449 std::integral_constant<bool, IsBlockVector<VectorType1>::value>(),
450 std::integral_constant<bool, IsBlockVector<VectorType2>::value>());
451 }
452
453
454
455 template <typename VectorType1, typename VectorType2>
456 inline void
457 BlockSparseMatrix::Tvmult(VectorType1 &dst, const VectorType2 &src) const
458 {
459 vmult(dst,
460 src,
461 true,
462 std::integral_constant<bool, IsBlockVector<VectorType1>::value>(),
463 std::integral_constant<bool, IsBlockVector<VectorType2>::value>());
464 }
465
466
467
468 template <typename VectorType1, typename VectorType2>
469 inline void
470 BlockSparseMatrix::vmult(VectorType1 & dst,
471 const VectorType2 &src,
472 const bool transpose,
473 std::integral_constant<bool, true>,
474 std::integral_constant<bool, true>) const
475 {
476 if (transpose == true)
478 else
480 }
481
482
483
484 template <typename VectorType1, typename VectorType2>
485 inline void
486 BlockSparseMatrix::vmult(VectorType1 & dst,
487 const VectorType2 &src,
488 const bool transpose,
489 std::integral_constant<bool, false>,
490 std::integral_constant<bool, true>) const
491 {
492 if (transpose == true)
494 else
496 }
497
498
499
500 template <typename VectorType1, typename VectorType2>
501 inline void
502 BlockSparseMatrix::vmult(VectorType1 & dst,
503 const VectorType2 &src,
504 const bool transpose,
505 std::integral_constant<bool, true>,
506 std::integral_constant<bool, false>) const
507 {
508 if (transpose == true)
510 else
512 }
513
514
515
516 template <typename VectorType1, typename VectorType2>
517 inline void
518 BlockSparseMatrix::vmult(VectorType1 & dst,
519 const VectorType2 &src,
520 const bool transpose,
521 std::integral_constant<bool, false>,
522 std::integral_constant<bool, false>) const
523 {
524 if (transpose == true)
526 else
528 }
529
530
531
532 inline std::vector<IndexSet>
534 {
535 Assert(this->n_block_cols() != 0, ExcNotInitialized());
536 Assert(this->n_block_rows() != 0, ExcNotInitialized());
537
538 std::vector<IndexSet> domain_indices;
539 for (size_type c = 0; c < this->n_block_cols(); ++c)
540 domain_indices.push_back(
541 this->sub_objects[0][c]->locally_owned_domain_indices());
542
543 return domain_indices;
544 }
545
546
547
548 inline std::vector<IndexSet>
550 {
551 Assert(this->n_block_cols() != 0, ExcNotInitialized());
552 Assert(this->n_block_rows() != 0, ExcNotInitialized());
553
554 std::vector<IndexSet> range_indices;
555 for (size_type r = 0; r < this->n_block_rows(); ++r)
556 range_indices.push_back(
557 this->sub_objects[r][0]->locally_owned_range_indices());
558
559 return range_indices;
560 }
561
562
563
564 namespace internal
565 {
566 namespace BlockLinearOperatorImplementation
567 {
582 template <typename PayloadBlockType>
584 {
585 public:
589 using BlockType = PayloadBlockType;
590
600 template <typename... Args>
601 TrilinosBlockPayload(const Args &...)
602 {
603 static_assert(
604 std::is_same<
605 PayloadBlockType,
607 "TrilinosBlockPayload can only accept a payload of type TrilinosPayload.");
608 }
609 };
610
611 } // namespace BlockLinearOperatorImplementation
612 } /* namespace internal */
613
614
615} /* namespace TrilinosWrappers */
616
617
619
620#endif // DEAL_II_WITH_TRILINOS
621
622#endif // dealii_trilinos_block_sparse_matrix_h
void Tvmult_nonblock_nonblock(VectorType &dst, const VectorType &src) const
void vmult_nonblock_nonblock(VectorType &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
unsigned int n_block_rows() const
MatrixIterator< BlockMatrixIterators::Accessor< BlockMatrixBase, true > > const_iterator
void vmult_block_block(BlockVectorType &dst, const BlockVectorType &src) const
types::global_dof_index size_type
typename BlockType::value_type value_type
void vmult_nonblock_block(VectorType &dst, const BlockVectorType &src) const
unsigned int n_block_cols() const
MatrixIterator< BlockMatrixIterators::Accessor< BlockMatrixBase, false > > iterator
BlockType & block(const unsigned int row, const unsigned int column)
void vmult_block_nonblock(BlockVectorType &dst, const VectorType &src) const
void Tvmult_block_nonblock(BlockVectorType &dst, const VectorType &src) const
BlockSparseMatrix & operator=(const BlockSparseMatrix &)=default
TrilinosScalar residual(MPI::BlockVector &dst, const MPI::BlockVector &x, const MPI::BlockVector &b) const
void Tvmult(VectorType1 &dst, const VectorType2 &src) const
std::vector< IndexSet > locally_owned_domain_indices() const
void reinit(const size_type n_block_rows, const size_type n_block_columns)
std::vector< IndexSet > locally_owned_range_indices() const
void vmult(VectorType1 &dst, const VectorType2 &src) const
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)
static ::ExceptionBase & ExcIncompatibleRowNumbers(int arg1, int arg2, int arg3, int arg4)
#define DeclException4(Exception4, type1, type2, type3, type4, outsequence)
Definition exceptions.h:579
static ::ExceptionBase & ExcScalarAssignmentOnlyForZeroValue()
#define Assert(cond, exc)
static ::ExceptionBase & ExcIncompatibleColNumbers(int arg1, int arg2, int arg3, int arg4)
static ::ExceptionBase & ExcNotInitialized()
double TrilinosScalar
Definition types.h:175