Loading [MathJax]/extensions/TeX/newcommand.js
 Reference documentation for deal.II version 9.2.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\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
sparse_decomposition.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2002 - 2018 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_sparse_decomposition_h
17 #define dealii_sparse_decomposition_h
18 
19 #include <deal.II/base/config.h>
20 
22 
23 #include <cmath>
24 
26 
109 template <typename number>
110 class SparseLUDecomposition : protected SparseMatrix<number>,
111  public virtual Subscriptor
112 {
113 protected:
121 
122 public:
127 
132  virtual ~SparseLUDecomposition() override = 0;
133 
138  virtual void
139  clear() override;
140 
145  {
146  public:
150  AdditionalData(const double strengthen_diagonal = 0,
151  const unsigned int extra_off_diagonals = 0,
152  const bool use_previous_sparsity = false,
153  const SparsityPattern *use_this_sparsity = nullptr);
154 
162 
173  unsigned int extra_off_diagonals;
174 
184 
196  };
197 
213  template <typename somenumber>
214  void
216  const AdditionalData parameters);
217 
222  bool
223  empty() const;
224 
230  size_type
231  m() const;
232 
238  size_type
239  n() const;
240 
248  template <class OutVector, class InVector>
249  void
250  vmult_add(OutVector &dst, const InVector &src) const;
251 
259  template <class OutVector, class InVector>
260  void
261  Tvmult_add(OutVector &dst, const InVector &src) const;
262 
267  virtual std::size_t
268  memory_consumption() const;
269 
279  double,
280  << "The strengthening parameter " << arg1
281  << " is not greater or equal than zero!");
283 protected:
288  template <typename somenumber>
289  void
291 
298  virtual void
300 
310  virtual number
311  get_strengthen_diagonal(const number rowsum, const size_type row) const;
312 
317 
323  std::vector<const size_type *> prebuilt_lower_bound;
324 
328  void
330 
331 private:
342 };
343 
345 //---------------------------------------------------------------------------
346 
347 #ifndef DOXYGEN
348 
349 template <typename number>
350 inline number
352  const number /*rowsum*/,
353  const size_type /*row*/) const
354 {
355  return strengthen_diagonal;
356 }
357 
358 
359 
360 template <typename number>
361 inline bool
363 {
365 }
366 
367 
368 template <typename number>
371 {
372  return SparseMatrix<number>::m();
373 }
374 
375 
376 template <typename number>
379 {
380  return SparseMatrix<number>::n();
381 }
382 
383 // Note: This function is required for full compatibility with
384 // the LinearOperator class. ::MatrixInterfaceWithVmultAdd
385 // picks up the vmult_add function in the protected SparseMatrix
386 // base class.
387 template <typename number>
388 template <class OutVector, class InVector>
389 inline void
391  const InVector &src) const
392 {
393  OutVector tmp;
394  tmp.reinit(dst);
395  this->vmult(tmp, src);
396  dst += tmp;
397 }
398 
399 // Note: This function is required for full compatibility with
400 // the LinearOperator class. ::MatrixInterfaceWithVmultAdd
401 // picks up the vmult_add function in the protected SparseMatrix
402 // base class.
403 template <typename number>
404 template <class OutVector, class InVector>
405 inline void
407  const InVector &src) const
408 {
409  OutVector tmp;
410  tmp.reinit(dst);
411  this->Tvmult(tmp, src);
412  dst += tmp;
413 }
414 
415 //---------------------------------------------------------------------------
416 
417 
418 template <typename number>
420  const double strengthen_diag,
421  const unsigned int extra_off_diag,
422  const bool use_prev_sparsity,
423  const SparsityPattern *use_this_spars)
424  : strengthen_diagonal(strengthen_diag)
425  , extra_off_diagonals(extra_off_diag)
426  , use_previous_sparsity(use_prev_sparsity)
427  , use_this_sparsity(use_this_spars)
428 {}
429 
430 
431 #endif // DOXYGEN
432 
434 
435 #endif // dealii_sparse_decomposition_h
SparseLUDecomposition::~SparseLUDecomposition
virtual ~SparseLUDecomposition() override=0
SparseLUDecomposition::AdditionalData::AdditionalData
AdditionalData(const double strengthen_diagonal=0, const unsigned int extra_off_diagonals=0, const bool use_previous_sparsity=false, const SparsityPattern *use_this_sparsity=nullptr)
sparse_matrix.h
SparseLUDecomposition::AdditionalData::use_this_sparsity
const SparsityPattern * use_this_sparsity
Definition: sparse_decomposition.h:195
SparseLUDecomposition::memory_consumption
virtual std::size_t memory_consumption() const
SparseMatrix::n
size_type n() const
SparseLUDecomposition::Tvmult_add
void Tvmult_add(OutVector &dst, const InVector &src) const
SparseMatrix
Definition: sparse_matrix.h:497
SparseLUDecomposition::own_sparsity
SparsityPattern * own_sparsity
Definition: sparse_decomposition.h:341
SparseLUDecomposition::AdditionalData
Definition: sparse_decomposition.h:144
SparseLUDecomposition::m
size_type m() const
SparseLUDecomposition::prebuilt_lower_bound
std::vector< const size_type * > prebuilt_lower_bound
Definition: sparse_decomposition.h:323
SparseLUDecomposition::get_strengthen_diagonal
virtual number get_strengthen_diagonal(const number rowsum, const size_type row) const
Subscriptor
Definition: subscriptor.h:62
SparseLUDecomposition::AdditionalData::extra_off_diagonals
unsigned int extra_off_diagonals
Definition: sparse_decomposition.h:173
SparseLUDecomposition
Definition: sparse_decomposition.h:110
SparseLUDecomposition::ExcInvalidStrengthening
static ::ExceptionBase & ExcInvalidStrengthening(double arg1)
SparseMatrix::empty
bool empty() const
SparseLUDecomposition::AdditionalData::strengthen_diagonal
double strengthen_diagonal
Definition: sparse_decomposition.h:161
SparseLUDecomposition::AdditionalData::use_previous_sparsity
bool use_previous_sparsity
Definition: sparse_decomposition.h:183
SparseLUDecomposition::size_type
typename SparseMatrix< number >::size_type size_type
Definition: sparse_decomposition.h:126
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:358
SparseLUDecomposition::initialize
void initialize(const SparseMatrix< somenumber > &matrix, const AdditionalData parameters)
LAPACKSupport::matrix
@ matrix
Contents is actually a matrix.
Definition: lapack_support.h:60
SparseLUDecomposition::strengthen_diagonal_impl
virtual void strengthen_diagonal_impl()
SparsityPattern
Definition: sparsity_pattern.h:865
DeclException1
#define DeclException1(Exception1, type1, outsequence)
Definition: exceptions.h:518
unsigned int
LinearAlgebra::CUDAWrappers::kernel::size_type
types::global_dof_index size_type
Definition: cuda_kernels.h:45
SparseLUDecomposition::clear
virtual void clear() override
SparseLUDecomposition::strengthen_diagonal
double strengthen_diagonal
Definition: sparse_decomposition.h:316
SparseLUDecomposition::n
size_type n() const
SparseLUDecomposition::SparseLUDecomposition
SparseLUDecomposition()
SparseLUDecomposition::empty
bool empty() const
config.h
SparseLUDecomposition::prebuild_lower_bound
void prebuild_lower_bound()
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:359
SparseLUDecomposition::copy_from
void copy_from(const SparseMatrix< somenumber > &matrix)
SparseLUDecomposition::vmult_add
void vmult_add(OutVector &dst, const InVector &src) const
SparseMatrix::m
size_type m() const