optimization_solver.h
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 #ifndef LIBMESH_OPTIMIZATION_SOLVER_H
21 #define LIBMESH_OPTIMIZATION_SOLVER_H
22 
23 // Local includes
24 #include "libmesh/libmesh_common.h"
27 #include "libmesh/libmesh.h"
29 #include "libmesh/auto_ptr.h"
31 
32 // C++ includes
33 #include <cstddef>
34 
35 namespace libMesh
36 {
37 
38 // forward declarations
39 template <typename T> class SparseMatrix;
40 template <typename T> class NumericVector;
41 template <typename T> class Preconditioner;
42 
50 template <typename T>
51 class OptimizationSolver : public ReferenceCountedObject<OptimizationSolver<T> >,
52  public ParallelObject
53 {
54 public:
59 
63  explicit
64  OptimizationSolver (sys_type & s);
65 
69  virtual ~OptimizationSolver ();
70 
75  static UniquePtr<OptimizationSolver<T> > build(sys_type & s,
76  const SolverPackage solver_package = libMesh::default_solver_package());
77 
82  bool initialized () const { return _is_initialized; }
83 
87  virtual void clear () {}
88 
92  virtual void init () = 0;
93 
97  virtual void solve () = 0;
98 
104  virtual void get_dual_variables()
105  { libmesh_not_implemented(); }
106 
111  virtual void print_converged_reason() { libmesh_not_implemented(); }
112 
120  virtual int get_converged_reason() { return 0; }
121 
127 
133 
139 
145 
150 
156 
161 
166 
171  const sys_type & system () const { return _system; }
172 
177  sys_type & system () { return _system; }
178 
183 
188 
192  bool verbose;
193 
194 protected:
195 
199  sys_type & _system;
200 
205 
206 };
207 
208 } // namespace libMesh
209 
210 
211 #endif // LIBMESH_OPTIMIZATION_SOLVER_H
OptimizationSystem::ComputeObjective * objective_object
static UniquePtr< OptimizationSolver< T > > build(sys_type &s, const SolverPackage solver_package=libMesh::default_solver_package())
OptimizationSystem::ComputeGradient * gradient_object
const sys_type & system() const
OptimizationSystem::ComputeHessian * hessian_object
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
OptimizationSystem::ComputeInequalityConstraints * inequality_constraints_object
SolverPackage default_solver_package()
Definition: libmesh.C:988
OptimizationSystem::ComputeEqualityConstraintsJacobian * equality_constraints_jacobian_object
unsigned int max_objective_function_evaluations
OptimizationSystem::ComputeLowerAndUpperBounds * lower_and_upper_bounds_object
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
OptimizationSystem::ComputeInequalityConstraintsJacobian * inequality_constraints_jacobian_object
OptimizationSystem::ComputeEqualityConstraints * equality_constraints_object