Reference documentation for deal.II version GIT relicensing-487-ge9eb5ab491 2024-04-25 07:20:02+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_tpetra_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) 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_tpetra_trilinos_block_sparse_matrix_h
16#define dealii_tpetra_trilinos_block_sparse_matrix_h
17
18
19#include <deal.II/base/config.h>
20
21#ifdef DEAL_II_TRILINOS_WITH_TPETRA
22
24
31
32# include <cmath>
33
35
36// forward declarations
37# ifndef DOXYGEN
39template <typename number>
41# endif
42
43namespace LinearAlgebra
44{
49 namespace TpetraWrappers
50 {
74 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
76 : public BlockMatrixBase<SparseMatrix<Number, MemorySpace>>
77 {
78 public:
83
88
93 using pointer = typename BaseClass::pointer;
98 using iterator = typename BaseClass::iterator;
100
112 BlockSparseMatrix() = default;
113
118
125
136 operator=(const Number d);
137
151 void
152 reinit(const size_type n_block_rows, const size_type n_block_columns);
153
159 template <typename BlockSparsityPatternType>
160 void
161 reinit(const std::vector<IndexSet> &input_maps,
162 const BlockSparsityPatternType &block_sparsity_pattern,
163 const MPI_Comm communicator = MPI_COMM_WORLD,
164 const bool exchange_data = false);
165
171 template <typename BlockSparsityPatternType>
172 void
173 reinit(const BlockSparsityPatternType &block_sparsity_pattern);
174
181 void
183 const std::vector<IndexSet> &parallel_partitioning,
184 const ::BlockSparseMatrix<double> &dealii_block_sparse_matrix,
185 const MPI_Comm communicator = MPI_COMM_WORLD,
186 const double drop_tolerance = 1e-13);
187
195 void
196 reinit(const ::BlockSparseMatrix<double> &deal_ii_sparse_matrix,
197 const double drop_tolerance = 1e-13);
198
206 bool
208
220 void
222
227 std::uint64_t
229
235
241 std::vector<IndexSet>
243
249 std::vector<IndexSet>
251
258 template <typename VectorType1, typename VectorType2>
259 void
260 vmult(VectorType1 &dst, const VectorType2 &src) const;
261
267 template <typename VectorType1, typename VectorType2>
268 void
269 Tvmult(VectorType1 &dst, const VectorType2 &src) const;
270
283 template <typename VectorType1,
284 typename VectorType2,
285 typename VectorType3>
286 Number
287 residual(VectorType1 &dst,
288 const VectorType2 &x,
289 const VectorType3 &b) const;
290
296
297 private:
301 template <typename VectorType1, typename VectorType2>
302 void
303 vmult(VectorType1 &dst,
304 const VectorType2 &src,
305 const bool transpose,
306 const std::bool_constant<true>,
307 const std::bool_constant<true>) const;
308
313 template <typename VectorType1, typename VectorType2>
314 void
315 vmult(VectorType1 &dst,
316 const VectorType2 &src,
317 const bool transpose,
318 const std::bool_constant<false>,
319 const std::bool_constant<true>) const;
320
325 template <typename VectorType1, typename VectorType2>
326 void
327 vmult(VectorType1 &dst,
328 const VectorType2 &src,
329 const bool transpose,
330 const std::bool_constant<true>,
331 const std::bool_constant<false>) const;
332
338 template <typename VectorType1, typename VectorType2>
339 void
340 vmult(VectorType1 &dst,
341 const VectorType2 &src,
342 const bool transpose,
343 const std::bool_constant<false>,
344 const std::bool_constant<false>) const;
345 };
346
347 } // namespace TpetraWrappers
348
351 // ------------- inline and template functions -----------------
352 namespace TpetraWrappers
353 {
354 template <typename Number, typename MemorySpace>
355 template <typename VectorType1, typename VectorType2>
356 inline void
358 const VectorType2 &src) const
359 {
360 vmult(dst,
361 src,
362 false,
363 std::bool_constant<IsBlockVector<VectorType1>::value>(),
364 std::bool_constant<IsBlockVector<VectorType2>::value>());
365 }
366
367
368
369 template <typename Number, typename MemorySpace>
370 template <typename VectorType1, typename VectorType2>
371 inline void
373 const VectorType2 &src) const
374 {
375 vmult(dst,
376 src,
377 true,
378 std::bool_constant<IsBlockVector<VectorType1>::value>(),
379 std::bool_constant<IsBlockVector<VectorType2>::value>());
380 }
381
382
383
384 template <typename Number, typename MemorySpace>
385 template <typename VectorType1, typename VectorType2>
386 inline void
388 VectorType1 &dst,
389 const VectorType2 &src,
390 const bool transpose,
391 std::bool_constant<true>,
392 std::bool_constant<true>) const
393 {
394 if (transpose == true)
395 BaseClass::Tvmult_block_block(dst, src);
396 else
397 BaseClass::vmult_block_block(dst, src);
398 }
399
400
401
402 template <typename Number, typename MemorySpace>
403 template <typename VectorType1, typename VectorType2>
404 inline void
406 VectorType1 &dst,
407 const VectorType2 &src,
408 const bool transpose,
409 std::bool_constant<false>,
410 std::bool_constant<true>) const
411 {
412 if (transpose == true)
413 BaseClass::Tvmult_nonblock_block(dst, src);
414 else
415 BaseClass::vmult_nonblock_block(dst, src);
416 }
417
418
419
420 template <typename Number, typename MemorySpace>
421 template <typename VectorType1, typename VectorType2>
422 inline void
424 VectorType1 &dst,
425 const VectorType2 &src,
426 const bool transpose,
427 std::bool_constant<true>,
428 std::bool_constant<false>) const
429 {
430 if (transpose == true)
431 BaseClass::Tvmult_block_nonblock(dst, src);
432 else
433 BaseClass::vmult_block_nonblock(dst, src);
434 }
435
436
437
438 template <typename Number, typename MemorySpace>
439 template <typename VectorType1, typename VectorType2>
440 inline void
442 VectorType1 &dst,
443 const VectorType2 &src,
444 const bool transpose,
445 std::bool_constant<false>,
446 std::bool_constant<false>) const
447 {
448 if (transpose == true)
449 BaseClass::Tvmult_nonblock_nonblock(dst, src);
450 else
451 BaseClass::vmult_nonblock_nonblock(dst, src);
452 }
453 } // namespace TpetraWrappers
454
455} // namespace LinearAlgebra
456
457
459
460#endif // DEAL_II_TRILINOS_WITH_TPETRA
461
462#endif // dealii_tpetra_trilinos_block_sparse_matrix_h
value_type & reference
const value_type & const_reference
unsigned int n_block_rows() const
MatrixIterator< BlockMatrixIterators::Accessor< BlockMatrixBase, true > > const_iterator
types::global_dof_index size_type
typename BlockType::value_type value_type
const value_type * const_pointer
MatrixIterator< BlockMatrixIterators::Accessor< BlockMatrixBase, false > > iterator
void reinit(const std::vector< IndexSet > &parallel_partitioning, const ::BlockSparseMatrix< double > &dealii_block_sparse_matrix, const MPI_Comm communicator=MPI_COMM_WORLD, const double drop_tolerance=1e-13)
void Tvmult(VectorType1 &dst, const VectorType2 &src) const
BlockSparseMatrix< Number, MemorySpace > & operator=(const BlockSparseMatrix< Number, MemorySpace > &)=default
void vmult(VectorType1 &dst, const VectorType2 &src) const
std::vector< IndexSet > locally_owned_domain_indices() const
std::vector< IndexSet > locally_owned_range_indices() const
void reinit(const ::BlockSparseMatrix< double > &deal_ii_sparse_matrix, const double drop_tolerance=1e-13)
Number residual(VectorType1 &dst, const VectorType2 &x, const VectorType3 &b) const
void reinit(const BlockSparsityPatternType &block_sparsity_pattern)
BlockSparseMatrix< Number, MemorySpace > & operator=(const Number d)
void reinit(const size_type n_block_rows, const size_type n_block_columns)
void reinit(const std::vector< IndexSet > &input_maps, const BlockSparsityPatternType &block_sparsity_pattern, const MPI_Comm communicator=MPI_COMM_WORLD, const bool exchange_data=false)
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
DerivativeForm< 1, spacedim, dim, Number > transpose(const DerivativeForm< 1, dim, spacedim, Number > &DF)