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_solver_direct.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_trilinos_tpetra_solver_direct_h
16#define dealii_trilinos_tpetra_solver_direct_h
17
18#include <deal.II/base/config.h>
19
21
22#include <string>
23
24#ifdef DEAL_II_TRILINOS_WITH_TPETRA
25# ifdef DEAL_II_TRILINOS_WITH_AMESOS2
26
27# include "deal.II/base/types.h"
28
32
33# include <Amesos2.hpp>
34# include <Teuchos_ConfigDefs.hpp>
35# include <Teuchos_ParameterList.hpp>
36# include <Teuchos_RCPDecl.hpp>
37
38
40
41namespace LinearAlgebra
42{
43 namespace TpetraWrappers
44 {
45 // forward declarations
46# ifndef DOXYGEN
47 template <typename Number, typename MemorySpace>
48 class SparseMatrix;
49# endif
50
54 std::string,
55 << "You tried to select the solver type <" << arg1 << ">\n"
56 << "but this solver is not supported by Trilinos/Amesos2\n"
57 << "due to one of the following reasons:\n"
58 << "* This solver does not exist\n"
59 << "* This solver is not (yet) supported by Trilinos/Amesos2\n"
60 << "* Trilinos/Amesos2 was not configured for its use.");
61
62
69 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
71 {
72 public:
77
81# if DEAL_II_TRILINOS_VERSION_GTE(14, 2, 0)
82 using NodeType = Tpetra::KokkosCompat::KokkosDeviceWrapperNode<
83 typename MemorySpace::kokkos_space::execution_space,
84 typename MemorySpace::kokkos_space>;
85# else
86 using NodeType = Kokkos::Compat::KokkosDeviceWrapperNode<
87 typename MemorySpace::kokkos_space::execution_space,
88 typename MemorySpace::kokkos_space>;
89# endif
90
94 using LinearOperator = Tpetra::Operator<Number, int, size_type, NodeType>;
95
99 using MultiVector = Tpetra::MultiVector<Number, int, size_type, NodeType>;
100
104 virtual ~SolverDirectBase() = default;
105
112 void
114
120 void
123
130 void
134
135
140 control() const;
141
146 int,
147 << "An error with error number " << arg1
148 << " occurred while calling a Trilinos function");
149
150 protected:
160 void
162
170
174 Teuchos::RCP<
175 Amesos2::Solver<typename SparseMatrix<Number, MemorySpace>::MatrixType,
178
179 /*
180 * The set solver type to be handed to the solver factory of Amesos2.
181 */
182 std::string solver_type;
183
189 Teuchos::ParameterList parameter_list;
190 }; // Base
191
192
193
207 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
208 class SolverDirect : public SolverDirectBase<Number, MemorySpace>
209 {
210 public:
217 {
218 AdditionalData(const std::string &solver_name);
219
220
244 std::string solver_name;
245 };
246
251 const AdditionalData &additional_data = AdditionalData());
252
259 void
260 set_pararameter_list(Teuchos::ParameterList &parameter_list);
261 }; // SolverDirect
262
263
274 template <typename Number, typename MemorySpace = ::MemorySpace::Host>
275 class SolverDirectKLU2 : public SolverDirectBase<Number, MemorySpace>
276 {
277 public:
284 {
285 AdditionalData(const std::string &transpose_mode = "NOTRANS",
286 const bool symmetric_mode = false,
287 const bool equilibrate_matrix = true,
288 const std::string &column_permutation = "COLAMD",
289 const std::string &iterative_refinement = "NO");
296 std::string transpose_mode;
297
302
307
316
325 };
326
331 SolverControl &cn,
332 const AdditionalData &additional_data = AdditionalData());
333 }; // KLU2
334
335 } // namespace TpetraWrappers
336} // namespace LinearAlgebra
337
339
340# endif // DEAL_II_TRILINOS_WITH_AMESOS2
341#endif // DEAL_II_TRILINOS_WITH_TPETRA
342
343#endif
void solve(const SparseMatrix< Number, MemorySpace > &A, Vector< Number, MemorySpace > &x, const Vector< Number, MemorySpace > &b)
Teuchos::RCP< Amesos2::Solver< typename SparseMatrix< Number, MemorySpace >::MatrixType, MultiVector > > solver
void initialize(const SparseMatrix< Number, MemorySpace > &A)
Tpetra::KokkosCompat::KokkosDeviceWrapperNode< typename MemorySpace::kokkos_space::execution_space, typename MemorySpace::kokkos_space > NodeType
SolverDirectBase(SolverControl &cn, const std::string &solver_type)
Tpetra::MultiVector< Number, int, size_type, NodeType > MultiVector
void solve(Vector< Number, MemorySpace > &x, const Vector< Number, MemorySpace > &b)
SolverDirectKLU2(SolverControl &cn, const AdditionalData &additional_data=AdditionalData())
void set_pararameter_list(Teuchos::ParameterList &parameter_list)
SolverDirect(SolverControl &cn, const AdditionalData &additional_data=AdditionalData())
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
static ::ExceptionBase & ExcTrilinosError(int arg1)
static ::ExceptionBase & ExcTrilinosAmesos2SolverUnsupported(std::string arg1)
The chosen Solver is not supported or configured with Amesos2.
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:516
int signed_global_dof_index
Definition types.h:92
AdditionalData(const std::string &transpose_mode="NOTRANS", const bool symmetric_mode=false, const bool equilibrate_matrix=true, const std::string &column_permutation="COLAMD", const std::string &iterative_refinement="NO")