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
sparse_direct.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 2001 - 2022 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_direct_h
17#define dealii_sparse_direct_h
18
19
20#include <deal.II/base/config.h>
21
24
28#include <deal.II/lac/vector.h>
29
30#ifdef DEAL_II_WITH_UMFPACK
31# include <umfpack.h>
32#endif
33
35
36namespace types
37{
42#ifdef SuiteSparse_long
43 using suitesparse_index = SuiteSparse_long;
44#else
45 using suitesparse_index = long int;
46#endif
47} // namespace types
48
88{
89public:
94
100 {};
101
102
108
112 ~SparseDirectUMFPACK() override;
113
125 void
126 initialize(const SparsityPattern &sparsity_pattern);
127
145 template <class Matrix>
146 void
147 factorize(const Matrix &matrix);
148
152 template <class Matrix>
153 void
154 initialize(const Matrix & matrix,
155 const AdditionalData additional_data = AdditionalData());
156
177 void
178 vmult(Vector<double> &dst, const Vector<double> &src) const;
179
183 void
184 vmult(BlockVector<double> &dst, const BlockVector<double> &src) const;
185
190 void
191 Tvmult(Vector<double> &dst, const Vector<double> &src) const;
192
196 void
197 Tvmult(BlockVector<double> &dst, const BlockVector<double> &src) const;
198
204 m() const;
205
211 n() const;
212
241 void
242 solve(Vector<double> &rhs_and_solution, const bool transpose = false) const;
243
255 void
256 solve(Vector<std::complex<double>> &rhs_and_solution,
257 const bool transpose = false) const;
258
262 void
263 solve(BlockVector<double> &rhs_and_solution,
264 const bool transpose = false) const;
265
269 void
270 solve(BlockVector<std::complex<double>> &rhs_and_solution,
271 const bool transpose = false) const;
272
279 template <class Matrix>
280 void
281 solve(const Matrix & matrix,
282 Vector<double> &rhs_and_solution,
283 const bool transpose = false);
284
288 template <class Matrix>
289 void
290 solve(const Matrix & matrix,
291 Vector<std::complex<double>> &rhs_and_solution,
292 const bool transpose = false);
293
297 template <class Matrix>
298 void
299 solve(const Matrix & matrix,
300 BlockVector<double> &rhs_and_solution,
301 const bool transpose = false);
302
306 template <class Matrix>
307 void
308 solve(const Matrix & matrix,
309 BlockVector<std::complex<double>> &rhs_and_solution,
310 const bool transpose = false);
311
323 std::string,
324 int,
325 << "UMFPACK routine " << arg1 << " returned error status " << arg2 << '.'
326 << "\n\n"
327 << ("A complete list of error codes can be found in the file "
328 "<bundled/umfpack/UMFPACK/Include/umfpack.h>."
329 "\n\n"
330 "That said, the two most common errors that can happen are "
331 "that your matrix cannot be factorized because it is "
332 "rank deficient, and that UMFPACK runs out of memory "
333 "because your problem is too large."
334 "\n\n"
335 "The first of these cases most often happens if you "
336 "forget terms in your bilinear form necessary to ensure "
337 "that the matrix has full rank, or if your equation has a "
338 "spatially variable coefficient (or nonlinearity) that is "
339 "supposed to be strictly positive but, for whatever "
340 "reasons, is negative or zero. In either case, you probably "
341 "want to check your assembly procedure. Similarly, a "
342 "matrix can be rank deficient if you forgot to apply the "
343 "appropriate boundary conditions. For example, the "
344 "Laplace equation for a problem where only Neumann boundary "
345 "conditions are posed (or where you forget to apply Dirichlet "
346 "boundary conditions) has exactly one eigenvalue equal to zero "
347 "and its rank is therefore deficient by one. Finally, the matrix "
348 "may be rank deficient because you are using a quadrature "
349 "formula with too few quadrature points."
350 "\n\n"
351 "The other common situation is that you run out of memory. "
352 "On a typical laptop or desktop, it should easily be possible "
353 "to solve problems with 100,000 unknowns in 2d. If you are "
354 "solving problems with many more unknowns than that, in "
355 "particular if you are in 3d, then you may be running out "
356 "of memory and you will need to consider iterative "
357 "solvers instead of the direct solver employed by "
358 "UMFPACK."));
359
360private:
365
371
379
383 void
384 clear();
385
392 template <typename number>
393 void
395
396 template <typename number>
397 void
399
400 template <typename number>
401 void
403
415 std::vector<types::suitesparse_index> Ap;
416 std::vector<types::suitesparse_index> Ai;
417 std::vector<double> Ax;
418 std::vector<double> Az;
419
423 std::vector<double> control;
424};
425
427
428#endif // dealii_sparse_direct_h
std::vector< double > Az
~SparseDirectUMFPACK() override
void initialize(const SparsityPattern &sparsity_pattern)
void Tvmult(Vector< double > &dst, const Vector< double > &src) const
size_type m() const
void solve(Vector< double > &rhs_and_solution, const bool transpose=false) const
std::vector< double > Ax
void sort_arrays(const SparseMatrixEZ< number > &)
size_type n() const
void factorize(const Matrix &matrix)
std::vector< double > control
void vmult(Vector< double > &dst, const Vector< double > &src) const
std::vector< types::suitesparse_index > Ap
std::vector< types::suitesparse_index > Ai
#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 & ExcUMFPACKError(std::string arg1, int arg2)
#define DeclException2(Exception2, type1, type2, outsequence)
Definition exceptions.h:533
Definition types.h:33
unsigned int global_dof_index
Definition types.h:82
long int suitesparse_index