deal.II version GIT relicensing-2202-g5ed228918a 2024-12-01 17:10:00+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
trilinos_block_sparse_matrix.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2008 - 2024 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_trilinos_block_sparse_matrix_h
16#define dealii_trilinos_block_sparse_matrix_h
17
18
19#include <deal.II/base/config.h>
20
21#ifdef DEAL_II_WITH_TRILINOS
22
24
30
31# include <cmath>
32
34
35// forward declarations
36# ifndef DOXYGEN
38template <typename number>
40# endif
41
42namespace TrilinosWrappers
43{
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
216
221 std::uint64_t
222 n_nonzero_elements() const;
223
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
303 const MPI::BlockVector &x,
304 const MPI::Vector &b) const;
305
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.");
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::bool_constant<true>,
362 const std::bool_constant<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::bool_constant<false>,
374 const std::bool_constant<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::bool_constant<true>,
386 const std::bool_constant<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::bool_constant<false>,
399 const std::bool_constant<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::bool_constant<IsBlockVector<VectorType1>::value>(),
449 std::bool_constant<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::bool_constant<IsBlockVector<VectorType1>::value>(),
462 std::bool_constant<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::bool_constant<true>,
473 std::bool_constant<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::bool_constant<false>,
489 std::bool_constant<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::bool_constant<true>,
505 std::bool_constant<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::bool_constant<false>,
521 std::bool_constant<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_v<
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 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:498
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:499
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:580
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:178