DOLFIN
DOLFIN C++ interface
TAOLinearBoundSolver.h
1 // Copyright (C) 2012 Corrado Maurini
2 //
3 // This file is part of DOLFIN.
4 //
5 // DOLFIN is free software: you can redistribute it and/or modify
6 // it under the terms of the GNU Lesser General Public License as published by
7 // the Free Software Foundation, either version 3 of the License, or
8 // (at your option) any later version.
9 //
10 // DOLFIN is distributed in the hope that it will be useful,
11 // but WITHOUT ANY WARRANTY; without even the implied warranty of
12 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 // GNU Lesser General Public License for more details.
14 //
15 // You should have received a copy of the GNU Lesser General Public License
16 // along with DOLFIN. If not, see <http://www.gnu.org/licenses/>.
17 
18 #ifndef _TAOLinearBoundSolver_H
19 #define _TAOLinearBoundSolver_H
20 
21 #ifdef HAS_PETSC
22 
23 #include <map>
24 #include <memory>
25 #include <petscksp.h>
26 #include <petscpc.h>
27 
28 #include <petsctao.h>
29 #include <dolfin/common/NoDeleter.h>
30 #include <dolfin/common/types.h>
31 
32 #include <dolfin/la/PETScObject.h>
33 #include <dolfin/la/KrylovSolver.h>
34 
35 namespace dolfin
36 {
37 
38  // Forward declarations
39  class GenericMatrix;
40  class GenericVector;
41  class PETScMatrix;
42  class PETScVector;
43  class PETScPreconditioner;
44  class PETScKrylovSolver;
45  class PETScKSPDeleter;
46 
76  class TAOLinearBoundSolver : public Variable, public PETScObject
77  {
78  public:
79 
81  explicit TAOLinearBoundSolver(MPI_Comm comm);
82 
84  TAOLinearBoundSolver(const std::string method = "default",
85  const std::string ksp_type = "default",
86  const std::string pc_type = "default");
87 
90 
93  std::size_t solve(const GenericMatrix& A, GenericVector& x,
94  const GenericVector& b, const GenericVector& xl,
95  const GenericVector& xu);
96 
99  std::size_t solve(const PETScMatrix& A, PETScVector& x,
100  const PETScVector& b,
101  const PETScVector& xl, const PETScVector& xu);
102 
104  void set_solver(const std::string&);
105 
107  void set_ksp(const std::string ksp_type = "default");
108 
110  Tao tao() const;
111 
113  static std::map<std::string, std::string> methods();
114 
116  static std::map<std::string, std::string> krylov_solvers();
117 
119  static std::map<std::string, std::string> preconditioners();
120 
123  {
124  Parameters p("tao_solver");
125 
126  p.add("monitor_convergence" , false);
127  p.add("report" , false);
128  p.add("gradient_absolute_tol" , 1.0e-8);
129  p.add("gradient_relative_tol" , 1.0e-8);
130  p.add("gradient_t_tol" , 0.0);
131  p.add("error_on_nonconvergence", true);
132  p.add("maximum_iterations" , 100);
133  p.add("options_prefix" , "default");
134 
135  Parameters ksp("krylov_solver");
137  p.add(ksp);
138 
139  return p;
140  }
141 
143  std::shared_ptr<const PETScMatrix> get_matrix() const;
144 
146  std::shared_ptr<const PETScVector> get_vector() const;
147 
148  private:
149 
150  // Set operators with GenericMatrix and GenericVector
151  void set_operators(std::shared_ptr<const GenericMatrix> A,
152  std::shared_ptr<const GenericVector> b);
153 
154  // Set operators with shared pointer to PETSc objects
155  void set_operators(std::shared_ptr<const PETScMatrix> A,
156  std::shared_ptr<const PETScVector> b);
157 
158  // Callback for changes in parameter values
159  void read_parameters();
160 
161  // Available ksp solvers
162  static const std::map<std::string, const KSPType> _ksp_methods;
163 
164  // Available tao solvers descriptions
165  static const std::map<std::string, std::string> _methods_descr;
166 
167  // Set options
168  void set_ksp_options();
169 
170  // Initialize TAO solver
171  void init(const std::string& method);
172 
173  // Tao solver pointer
174  Tao _tao;
175 
176  // Petsc preconditioner
177  std::shared_ptr<PETScPreconditioner> _preconditioner;
178 
179  // Operator (the matrix) and the vector
180  std::shared_ptr<const PETScMatrix> _matA;
181  std::shared_ptr<const PETScVector> _b;
182 
183  bool _preconditioner_set;
184 
185  // Computes the value of the objective function and its gradient.
186  static PetscErrorCode
187  __TAOFormFunctionGradientQuadraticProblem(Tao tao, Vec X,
188  PetscReal *ener, Vec G,
189  void *ptr);
190 
191  // Computes the hessian of the quadratic objective function
192  static PetscErrorCode
193  __TAOFormHessianQuadraticProblem(Tao tao,Vec X, Mat H, Mat Hpre,
194  void *ptr);
195 
196  //-------------------------------------------------------------------------
197  // Monitor the state of the solution at each iteration. The
198  // output printed to the screen is:
199  //
200  // iterate - the current iterate number (>=0)
201  // f - the current function value
202  // gnorm - the square of the gradient norm, duality gap, or other
203  // measure
204  // indicating distance from optimality.
205  // cnorm - the infeasibility of the current solution with regard
206  // to the constraints.
207  // xdiff - the step length or trust region radius of the most
208  // recent iterate.
209  //-------------------------------------------------------------------------
210  static PetscErrorCode __TAOMonitor(Tao tao, void *ctx);
211 
212  };
213 
214 }
215 
216 #endif
217 #endif
void set_solver(const std::string &)
Set the TAO solver type.
Definition: TAOLinearBoundSolver.cpp:272
Tao tao() const
Return TAO solver pointer.
Definition: TAOLinearBoundSolver.cpp:327
Definition: PETScVector.h:60
Common base class for DOLFIN variables.
Definition: Variable.h:35
std::shared_ptr< const PETScVector > get_vector() const
Return load vector shared pointer.
Definition: TAOLinearBoundSolver.cpp:337
static std::map< std::string, std::string > preconditioners()
Return a list of available preconditioners.
Definition: TAOLinearBoundSolver.cpp:74
Definition: adapt.h:29
void add(std::string key)
Definition: Parameters.h:128
Definition: TAOLinearBoundSolver.h:76
~TAOLinearBoundSolver()
Destructor.
Definition: TAOLinearBoundSolver.cpp:119
Definition: PETScMatrix.h:58
TAOLinearBoundSolver(MPI_Comm comm)
Create TAO bound constrained solver.
Definition: TAOLinearBoundSolver.cpp:79
std::shared_ptr< const PETScMatrix > get_matrix() const
Return Matrix shared pointer.
Definition: TAOLinearBoundSolver.cpp:332
static std::map< std::string, std::string > methods()
Return a list of available Tao solver methods.
Definition: TAOLinearBoundSolver.cpp:64
Definition: Parameters.h:94
std::size_t solve(const GenericMatrix &A, GenericVector &x, const GenericVector &b, const GenericVector &xl, const GenericVector &xu)
Definition: TAOLinearBoundSolver.cpp:142
This class defines a common interface for matrices.
Definition: GenericMatrix.h:46
static std::map< std::string, std::string > krylov_solvers()
Return a list of available krylov solvers.
Definition: TAOLinearBoundSolver.cpp:69
void set_ksp(const std::string ksp_type="default")
Set PETSC Krylov Solver (ksp) used by TAO.
Definition: TAOLinearBoundSolver.cpp:302
Definition: PETScObject.h:33
static Parameters default_parameters()
Default parameter values.
Definition: KrylovSolver.cpp:32
This class defines a common interface for vectors.
Definition: GenericVector.h:47
static Parameters default_parameters()
Default parameter values.
Definition: TAOLinearBoundSolver.h:122