libMesh
linear_solver.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 // C++ includes
21 
22 // Local Includes
23 #include "libmesh/libmesh_logging.h"
24 #include "libmesh/auto_ptr.h"
25 #include "libmesh/linear_solver.h"
26 #include "libmesh/laspack_linear_solver.h"
27 #include "libmesh/eigen_sparse_linear_solver.h"
28 #include "libmesh/petsc_linear_solver.h"
29 #include "libmesh/trilinos_aztec_linear_solver.h"
30 #include "libmesh/preconditioner.h"
31 #include "libmesh/sparse_matrix.h"
32 #include "libmesh/string_to_enum.h"
33 #include "libmesh/solver_configuration.h"
34 
35 namespace libMesh
36 {
37 
38 //------------------------------------------------------------------
39 // LinearSolver members
40 template <typename T>
41 UniquePtr<LinearSolver<T>>
43  const SolverPackage solver_package)
44 {
45  // Avoid unused parameter warnings when no solver packages are enabled.
46  libmesh_ignore(comm);
47 
48  // Build the appropriate solver
49  switch (solver_package)
50  {
51 #ifdef LIBMESH_HAVE_LASPACK
52  case LASPACK_SOLVERS:
54 #endif
55 
56 
57 #ifdef LIBMESH_HAVE_PETSC
58  case PETSC_SOLVERS:
60 #endif
61 
62 
63 #ifdef LIBMESH_TRILINOS_HAVE_AZTECOO
64  case TRILINOS_SOLVERS:
66 #endif
67 
68 
69 #ifdef LIBMESH_HAVE_EIGEN
70  case EIGEN_SOLVERS:
72 #endif
73 
74  default:
75  libmesh_error_msg("ERROR: Unrecognized solver package: " << solver_package);
76  }
77 
78  return UniquePtr<LinearSolver<T>>();
79 }
80 
81 template <typename T>
84 {
85  if (_preconditioner)
86  return _preconditioner->type();
87 
88  return _preconditioner_type;
89 }
90 
91 template <typename T>
92 void
94 {
95  if (_preconditioner)
96  _preconditioner->set_type(pct);
97  else
98  _preconditioner_type = pct;
99 }
100 
101 template <typename T>
102 void
104 {
105  if (this->_is_initialized)
106  libmesh_error_msg("Preconditioner must be attached before the solver is initialized!");
107 
108  _preconditioner_type = SHELL_PRECOND;
109  _preconditioner = preconditioner;
110 }
111 
112 template <typename T>
113 void
115 {
116  same_preconditioner = reuse_flag;
117 }
118 
119 template <typename T>
120 void
121 LinearSolver<T>::restrict_solve_to(const std::vector<unsigned int> * const dofs,
122  const SubsetSolveMode /*subset_solve_mode*/)
123 {
124  if (dofs != libmesh_nullptr)
125  libmesh_not_implemented();
126 }
127 
128 
129 template <typename T>
130 std::pair<unsigned int, Real> LinearSolver<T>::adjoint_solve (SparseMatrix<T> & mat,
131  NumericVector<T> & sol,
132  NumericVector<T> & rhs,
133  const double tol,
134  const unsigned int n_iter)
135 {
136  // Log how long the linear solve takes.
137  LOG_SCOPE("adjoint_solve()", "LinearSolver");
138 
139  // Take the discrete adjoint
140  mat.close();
141  mat.get_transpose(mat);
142 
143  // Call the solve function for the relevant linear algebra library and
144  // solve the transpose matrix
145  const std::pair<unsigned int, Real> totalrval = this->solve (mat, sol, rhs, tol, n_iter);
146 
147  // Now transpose back and restore the original matrix
148  // by taking the discrete adjoint
149  mat.get_transpose(mat);
150 
151  return totalrval;
152 }
153 
154 template <typename T>
156 {
157  LinearConvergenceReason reason = this->get_converged_reason();
158  libMesh::out << "Linear solver convergence/divergence reason: " << Utility::enum_to_string(reason) << std::endl;
159 }
160 
161 template <typename T>
163 {
164  _solver_configuration = &solver_configuration;
165 }
166 
167 //------------------------------------------------------------------
168 // Explicit instantiations
169 template class LinearSolver<Number>;
170 
171 
172 
173 } // namespace libMesh
virtual void get_transpose(SparseMatrix< T > &dest) const =0
Copies the transpose of the matrix into dest, which may be *this.
Encapsulates the MPI_Comm object.
Definition: parallel.h:657
static UniquePtr< LinearSolver< T > > build(const libMesh::Parallel::Communicator &comm_in, const SolverPackage solver_package=libMesh::default_solver_package())
Builds a LinearSolver using the linear solver package specified by solver_package.
Definition: linear_solver.C:42
SolverPackage
Defines an enum for various linear solver packages.
This class provides an interface to Eigen iterative solvers that is compatible with the libMesh Linea...
const class libmesh_nullptr_t libmesh_nullptr
Numeric vector.
Definition: dof_map.h:66
PreconditionerType preconditioner_type() const
Definition: linear_solver.C:83
The libMesh namespace provides an interface to certain functionality in the library.
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
Generic sparse matrix.
Definition: dof_map.h:65
This class provides an interface to Laspack iterative solvers that is compatible with the libMesh Lin...
PreconditionerType
Defines an enum for preconditioner types.
This class provides a uniform interface for preconditioners.
bool _is_initialized
Flag that tells if init() has been called.
Definition: libmesh.C:255
void set_preconditioner_type(const PreconditionerType pct)
Sets the type of preconditioner to use.
Definition: linear_solver.C:93
virtual std::pair< unsigned int, Real > adjoint_solve(SparseMatrix< T > &, NumericVector< T > &, NumericVector< T > &, const double, const unsigned int)
Function to solve the adjoint system.
This class stores solver configuration data, e.g.
void attach_preconditioner(Preconditioner< T > *preconditioner)
Attaches a Preconditioner object to be used.
void set_solver_configuration(SolverConfiguration &solver_configuration)
Set the solver configuration object.
std::string enum_to_string(const T e)
void libmesh_ignore(const T &)
This class provides an interface to AztecOO iterative solvers that is compatible with the libMesh Lin...
virtual void close()=0
Calls the SparseMatrix&#39;s internal assembly routines, ensuring that the values are consistent across p...
virtual void reuse_preconditioner(bool)
Set the same_preconditioner flag, which indicates if we reuse the same preconditioner for subsequent ...
This class provides an interface to PETSc iterative solvers that is compatible with the libMesh Linea...
OStreamProxy out
virtual void print_converged_reason() const
Prints a useful message about why the latest linear solve con(di)verged.
LinearConvergenceReason
Linear solver convergence flags (taken from the PETSc flags)
virtual void restrict_solve_to(const std::vector< unsigned int > *const dofs, const SubsetSolveMode subset_solve_mode=SUBSET_ZERO)
After calling this method, all successive solves will be restricted to the given set of dofs...
SubsetSolveMode
defines an enum for the question what happens to the dofs outside the given subset when a system is s...