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
matrix_tools.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 1998 - 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_matrix_tools_h
17#define dealii_matrix_tools_h
18
19
20#include <deal.II/base/config.h>
21
24
26
28
29#include <map>
30
31#ifdef DEAL_II_WITH_PETSC
32# include <petscsys.h>
33#endif
34
36
37
38// forward declarations
39#ifndef DOXYGEN
40template <int dim>
41class Quadrature;
42
43
44template <typename number>
45class Vector;
46template <typename number>
47class FullMatrix;
48template <typename number>
49class SparseMatrix;
50
51template <typename number>
53template <typename Number>
54class BlockVector;
55
56template <int dim, int spacedim>
57class Mapping;
58template <int dim, int spacedim>
60class DoFHandler;
61
62namespace hp
63{
64 template <int>
65 class QCollection;
66 template <int, int>
67 class MappingCollection;
68} // namespace hp
69
70
71# ifdef DEAL_II_WITH_PETSC
72namespace PETScWrappers
73{
74 class MatrixBase;
75 class VectorBase;
76 namespace MPI
77 {
79 class BlockVector;
80 } // namespace MPI
81} // namespace PETScWrappers
82# endif
83
84# ifdef DEAL_II_WITH_TRILINOS
85namespace TrilinosWrappers
86{
87 class SparseMatrix;
89 namespace MPI
90 {
91 class Vector;
92 class BlockVector;
93 } // namespace MPI
94} // namespace TrilinosWrappers
95# endif
96#endif
97
98
99
307namespace MatrixTools
308{
314 using namespace MatrixCreator;
315
320 template <typename number>
321 void
323 const std::map<types::global_dof_index, number> &boundary_values,
324 SparseMatrix<number> & matrix,
325 Vector<number> & solution,
326 Vector<number> & right_hand_side,
327 const bool eliminate_columns = true);
328
334 template <typename number>
335 void
337 const std::map<types::global_dof_index, number> &boundary_values,
339 BlockVector<number> & solution,
340 BlockVector<number> & right_hand_side,
341 const bool eliminate_columns = true);
342
343#ifdef DEAL_II_WITH_PETSC
351 void
353 const std::map<types::global_dof_index, PetscScalar> &boundary_values,
355 PETScWrappers::VectorBase & solution,
356 PETScWrappers::VectorBase & right_hand_side,
357 const bool eliminate_columns = true);
358
362 void
364 const std::map<types::global_dof_index, PetscScalar> &boundary_values,
367 PETScWrappers::MPI::BlockVector & right_hand_side,
368 const bool eliminate_columns = true);
369
370#endif
371
372#ifdef DEAL_II_WITH_TRILINOS
406 void
408 const std::map<types::global_dof_index, TrilinosScalar> &boundary_values,
411 TrilinosWrappers::MPI::Vector & right_hand_side,
412 const bool eliminate_columns = true);
413
418 void
420 const std::map<types::global_dof_index, TrilinosScalar> &boundary_values,
423 TrilinosWrappers::MPI::BlockVector & right_hand_side,
424 const bool eliminate_columns = true);
425#endif
426
445 template <typename number>
446 void
448 const std::map<types::global_dof_index, number> &boundary_values,
449 const std::vector<types::global_dof_index> & local_dof_indices,
450 FullMatrix<number> & local_matrix,
451 Vector<number> & local_rhs,
452 const bool eliminate_columns);
453
458 "You are providing a matrix whose subdivision into "
459 "blocks in either row or column direction does not use "
460 "the same blocks sizes as the solution vector or "
461 "right hand side vectors, respectively.");
462} // namespace MatrixTools
463
464
465
467
468#endif
Abstract base class for mapping classes.
Definition mapping.h:317
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:160
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcBlocksDontMatch()
#define DeclExceptionMsg(Exception, defaulttext)
Definition exceptions.h:488
void apply_boundary_values(const std::map< types::global_dof_index, number > &boundary_values, SparseMatrix< number > &matrix, Vector< number > &solution, Vector< number > &right_hand_side, const bool eliminate_columns=true)
void local_apply_boundary_values(const std::map< types::global_dof_index, number > &boundary_values, const std::vector< types::global_dof_index > &local_dof_indices, FullMatrix< number > &local_matrix, Vector< number > &local_rhs, const bool eliminate_columns)
Definition hp.h:118