libMesh
tao_optimization_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 #include "libmesh/libmesh_common.h"
21 
22 #if defined(LIBMESH_HAVE_PETSC_TAO) && !defined(LIBMESH_USE_COMPLEX_NUMBERS)
23 
24 
25 // C++ includes
26 
27 // Local Includes
28 #include "libmesh/petsc_vector.h"
29 #include "libmesh/petsc_matrix.h"
30 #include "libmesh/dof_map.h"
31 #include "libmesh/tao_optimization_solver.h"
32 #include "libmesh/equation_systems.h"
33 
34 namespace libMesh
35 {
36 
37 //--------------------------------------------------------------------
38 // Functions with C linkage to pass to Tao. Tao will call these
39 // methods as needed.
40 //
41 // Since they must have C linkage they have no knowledge of a namespace.
42 // Give them an obscure name to avoid namespace pollution.
43 extern "C"
44 {
45 
46  //---------------------------------------------------------------
47  // This function is called by Tao to evaluate objective function at x
48  PetscErrorCode
49  __libmesh_tao_objective (Tao /*tao*/, Vec x, PetscReal * objective, void * ctx)
50  {
51  LOG_SCOPE("objective()", "TaoOptimizationSolver");
52 
53  PetscErrorCode ierr = 0;
54 
55  libmesh_assert(x);
56  libmesh_assert(objective);
57  libmesh_assert(ctx);
58 
59  // ctx should be a pointer to the solver (it was passed in as void *)
61  static_cast<TaoOptimizationSolver<Number> *> (ctx);
62 
63  OptimizationSystem & sys = solver->system();
64 
65  // We'll use current_local_solution below, so let's ensure that it's consistent
66  // with the vector x that was passed in.
67  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
68  PetscVector<Number> X(x, sys.comm());
69 
70  // Perform a swap so that sys.solution points to X
71  X.swap(X_sys);
72  // Impose constraints on X
73  sys.get_dof_map().enforce_constraints_exactly(sys);
74  // Update sys.current_local_solution based on X
75  sys.update();
76  // Swap back
77  X.swap(X_sys);
78 
79  if (solver->objective_object != libmesh_nullptr)
80  (*objective) = solver->objective_object->objective(*(sys.current_local_solution), sys);
81  else
82  libmesh_error_msg("Objective function not defined in __libmesh_tao_objective");
83 
84  return ierr;
85  }
86 
87 
88 
89  //---------------------------------------------------------------
90  // This function is called by Tao to evaluate the gradient at x
91  PetscErrorCode
92  __libmesh_tao_gradient(Tao /*tao*/, Vec x, Vec g, void * ctx)
93  {
94  LOG_SCOPE("gradient()", "TaoOptimizationSolver");
95 
96  PetscErrorCode ierr = 0;
97 
98  libmesh_assert(x);
99  libmesh_assert(g);
100  libmesh_assert(ctx);
101 
102  // ctx should be a pointer to the solver (it was passed in as void *)
104  static_cast<TaoOptimizationSolver<Number> *> (ctx);
105 
106  OptimizationSystem & sys = solver->system();
107 
108  // We'll use current_local_solution below, so let's ensure that it's consistent
109  // with the vector x that was passed in.
110  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
111  PetscVector<Number> X(x, sys.comm());
112 
113  // Perform a swap so that sys.solution points to X
114  X.swap(X_sys);
115  // Impose constraints on X
116  sys.get_dof_map().enforce_constraints_exactly(sys);
117  // Update sys.current_local_solution based on X
118  sys.update();
119  // Swap back
120  X.swap(X_sys);
121 
122  // We'll also pass the gradient in to the assembly routine
123  // so let's make a PETSc vector for that too.
124  PetscVector<Number> gradient(g, sys.comm());
125 
126  // Clear the gradient prior to assembly
127  gradient.zero();
128 
129  if (solver->gradient_object != libmesh_nullptr)
130  solver->gradient_object->gradient(*(sys.current_local_solution), gradient, sys);
131  else
132  libmesh_error_msg("Gradient function not defined in __libmesh_tao_gradient");
133 
134  gradient.close();
135 
136  return ierr;
137  }
138 
139  //---------------------------------------------------------------
140  // This function is called by Tao to evaluate the Hessian at x
141  PetscErrorCode
142  __libmesh_tao_hessian(Tao /*tao*/, Vec x, Mat h, Mat pc, void * ctx)
143  {
144  LOG_SCOPE("hessian()", "TaoOptimizationSolver");
145 
146  PetscErrorCode ierr = 0;
147 
148  libmesh_assert(x);
149  libmesh_assert(h);
150  libmesh_assert(pc);
151  libmesh_assert(ctx);
152 
153  // ctx should be a pointer to the solver (it was passed in as void *)
155  static_cast<TaoOptimizationSolver<Number> *> (ctx);
156 
157  OptimizationSystem & sys = solver->system();
158 
159  // We'll use current_local_solution below, so let's ensure that it's consistent
160  // with the vector x that was passed in.
161  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
162  PetscVector<Number> X(x, sys.comm());
163 
164  // Perform a swap so that sys.solution points to X
165  X.swap(X_sys);
166  // Impose constraints on X
167  sys.get_dof_map().enforce_constraints_exactly(sys);
168  // Update sys.current_local_solution based on X
169  sys.update();
170  // Swap back
171  X.swap(X_sys);
172 
173  // Let's also wrap pc and h in PetscMatrix objects for convenience
174  PetscMatrix<Number> PC(pc, sys.comm());
175  PetscMatrix<Number> hessian(h, sys.comm());
176  PC.attach_dof_map(sys.get_dof_map());
177  hessian.attach_dof_map(sys.get_dof_map());
178 
179  if (solver->hessian_object != libmesh_nullptr)
180  {
181  // Following PetscNonlinearSolver by passing in PC. It's not clear
182  // why we pass in PC and not hessian though?
183  solver->hessian_object->hessian(*(sys.current_local_solution), PC, sys);
184  }
185  else
186  libmesh_error_msg("Hessian function not defined in __libmesh_tao_hessian");
187 
188  PC.close();
189  hessian.close();
190 
191  return ierr;
192  }
193 
194 
195  //---------------------------------------------------------------
196  // This function is called by Tao to evaluate the equality constraints at x
197  PetscErrorCode
198  __libmesh_tao_equality_constraints(Tao /*tao*/, Vec x, Vec ce, void * ctx)
199  {
200  LOG_SCOPE("equality_constraints()", "TaoOptimizationSolver");
201 
202  PetscErrorCode ierr = 0;
203 
204  libmesh_assert(x);
205  libmesh_assert(ce);
206  libmesh_assert(ctx);
207 
208  // ctx should be a pointer to the solver (it was passed in as void *)
210  static_cast<TaoOptimizationSolver<Number> *> (ctx);
211 
212  OptimizationSystem & sys = solver->system();
213 
214  // We'll use current_local_solution below, so let's ensure that it's consistent
215  // with the vector x that was passed in.
216  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
217  PetscVector<Number> X(x, sys.comm());
218 
219  // Perform a swap so that sys.solution points to X
220  X.swap(X_sys);
221  // Impose constraints on X
222  sys.get_dof_map().enforce_constraints_exactly(sys);
223  // Update sys.current_local_solution based on X
224  sys.update();
225  // Swap back
226  X.swap(X_sys);
227 
228  // We'll also pass the constraints vector ce into the assembly routine
229  // so let's make a PETSc vector for that too.
230  PetscVector<Number> eq_constraints(ce, sys.comm());
231 
232  // Clear the gradient prior to assembly
233  eq_constraints.zero();
234 
236  solver->equality_constraints_object->equality_constraints(*(sys.current_local_solution), eq_constraints, sys);
237  else
238  libmesh_error_msg("Constraints function not defined in __libmesh_tao_equality_constraints");
239 
240  eq_constraints.close();
241 
242  return ierr;
243  }
244 
245  //---------------------------------------------------------------
246  // This function is called by Tao to evaluate the Jacobian of the
247  // equality constraints at x
248  PetscErrorCode
249  __libmesh_tao_equality_constraints_jacobian(Tao /*tao*/, Vec x, Mat J, Mat Jpre, void * ctx)
250  {
251  LOG_SCOPE("equality_constraints_jacobian()", "TaoOptimizationSolver");
252 
253  PetscErrorCode ierr = 0;
254 
255  libmesh_assert(x);
256  libmesh_assert(J);
257  libmesh_assert(Jpre);
258 
259  // ctx should be a pointer to the solver (it was passed in as void *)
261  static_cast<TaoOptimizationSolver<Number> *> (ctx);
262 
263  OptimizationSystem & sys = solver->system();
264 
265  // We'll use current_local_solution below, so let's ensure that it's consistent
266  // with the vector x that was passed in.
267  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
268  PetscVector<Number> X(x, sys.comm());
269 
270  // Perform a swap so that sys.solution points to X
271  X.swap(X_sys);
272  // Impose constraints on X
273  sys.get_dof_map().enforce_constraints_exactly(sys);
274  // Update sys.current_local_solution based on X
275  sys.update();
276  // Swap back
277  X.swap(X_sys);
278 
279  // Let's also wrap J and Jpre in PetscMatrix objects for convenience
280  PetscMatrix<Number> J_petsc(J, sys.comm());
281  PetscMatrix<Number> Jpre_petsc(Jpre, sys.comm());
282 
284  solver->equality_constraints_jacobian_object->equality_constraints_jacobian(*(sys.current_local_solution), J_petsc, sys);
285  else
286  libmesh_error_msg("Constraints function not defined in __libmesh_tao_equality_constraints_jacobian");
287 
288  J_petsc.close();
289  Jpre_petsc.close();
290 
291  return ierr;
292  }
293 
294  //---------------------------------------------------------------
295  // This function is called by Tao to evaluate the inequality constraints at x
296  PetscErrorCode
297  __libmesh_tao_inequality_constraints(Tao /*tao*/, Vec x, Vec cineq, void * ctx)
298  {
299  LOG_SCOPE("inequality_constraints()", "TaoOptimizationSolver");
300 
301  PetscErrorCode ierr = 0;
302 
303  libmesh_assert(x);
304  libmesh_assert(cineq);
305  libmesh_assert(ctx);
306 
307  // ctx should be a pointer to the solver (it was passed in as void *)
309  static_cast<TaoOptimizationSolver<Number> *> (ctx);
310 
311  OptimizationSystem & sys = solver->system();
312 
313  // We'll use current_local_solution below, so let's ensure that it's consistent
314  // with the vector x that was passed in.
315  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
316  PetscVector<Number> X(x, sys.comm());
317 
318  // Perform a swap so that sys.solution points to X
319  X.swap(X_sys);
320  // Impose constraints on X
321  sys.get_dof_map().enforce_constraints_exactly(sys);
322  // Update sys.current_local_solution based on X
323  sys.update();
324  // Swap back
325  X.swap(X_sys);
326 
327  // We'll also pass the constraints vector ce into the assembly routine
328  // so let's make a PETSc vector for that too.
329  PetscVector<Number> ineq_constraints(cineq, sys.comm());
330 
331  // Clear the gradient prior to assembly
332  ineq_constraints.zero();
333 
335  solver->inequality_constraints_object->inequality_constraints(*(sys.current_local_solution), ineq_constraints, sys);
336  else
337  libmesh_error_msg("Constraints function not defined in __libmesh_tao_inequality_constraints");
338 
339  ineq_constraints.close();
340 
341  return ierr;
342  }
343 
344  //---------------------------------------------------------------
345  // This function is called by Tao to evaluate the Jacobian of the
346  // equality constraints at x
347  PetscErrorCode
348  __libmesh_tao_inequality_constraints_jacobian(Tao /*tao*/, Vec x, Mat J, Mat Jpre, void * ctx)
349  {
350  LOG_SCOPE("inequality_constraints_jacobian()", "TaoOptimizationSolver");
351 
352  PetscErrorCode ierr = 0;
353 
354  libmesh_assert(x);
355  libmesh_assert(J);
356  libmesh_assert(Jpre);
357 
358  // ctx should be a pointer to the solver (it was passed in as void *)
360  static_cast<TaoOptimizationSolver<Number> *> (ctx);
361 
362  OptimizationSystem & sys = solver->system();
363 
364  // We'll use current_local_solution below, so let's ensure that it's consistent
365  // with the vector x that was passed in.
366  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
367  PetscVector<Number> X(x, sys.comm());
368 
369  // Perform a swap so that sys.solution points to X
370  X.swap(X_sys);
371  // Impose constraints on X
372  sys.get_dof_map().enforce_constraints_exactly(sys);
373  // Update sys.current_local_solution based on X
374  sys.update();
375  // Swap back
376  X.swap(X_sys);
377 
378  // Let's also wrap J and Jpre in PetscMatrix objects for convenience
379  PetscMatrix<Number> J_petsc(J, sys.comm());
380  PetscMatrix<Number> Jpre_petsc(Jpre, sys.comm());
381 
383  solver->inequality_constraints_jacobian_object->inequality_constraints_jacobian(*(sys.current_local_solution), J_petsc, sys);
384  else
385  libmesh_error_msg("Constraints function not defined in __libmesh_tao_inequality_constraints_jacobian");
386 
387  J_petsc.close();
388  Jpre_petsc.close();
389 
390  return ierr;
391  }
392 
393 } // end extern "C"
394 //---------------------------------------------------------------------
395 
396 
397 
398 //---------------------------------------------------------------------
399 // TaoOptimizationSolver<> methods
400 template <typename T>
402  OptimizationSolver<T>(system_in),
403  _reason(TAO_CONVERGED_USER) // Arbitrary initial value...
404 {
405 }
406 
407 
408 
409 template <typename T>
411 {
412  this->clear ();
413 }
414 
415 
416 
417 template <typename T>
419 {
420  if (this->initialized())
421  {
422  this->_is_initialized = false;
423 
424  PetscErrorCode ierr=0;
425 
426  ierr = TaoDestroy(&_tao);
427  LIBMESH_CHKERR(ierr);
428  }
429 }
430 
431 
432 
433 template <typename T>
435 {
436  // Initialize the data structures if not done so already.
437  if (!this->initialized())
438  {
439  this->_is_initialized = true;
440 
441  PetscErrorCode ierr=0;
442 
443  ierr = TaoCreate(this->comm().get(),&_tao);
444  LIBMESH_CHKERR(ierr);
445  }
446 }
447 
448 template <typename T>
450 {
451  LOG_SCOPE("solve()", "TaoOptimizationSolver");
452 
453  this->init ();
454 
455  this->system().solution->zero();
456 
457  PetscMatrix<T> * hessian = cast_ptr<PetscMatrix<T> *>(this->system().matrix);
458  // PetscVector<T> * gradient = cast_ptr<PetscVector<T> *>(this->system().rhs);
459  PetscVector<T> * x = cast_ptr<PetscVector<T> *>(this->system().solution.get());
460  PetscVector<T> * ceq = cast_ptr<PetscVector<T> *>(this->system().C_eq.get());
461  PetscMatrix<T> * ceq_jac = cast_ptr<PetscMatrix<T> *>(this->system().C_eq_jac.get());
462  PetscVector<T> * cineq = cast_ptr<PetscVector<T> *>(this->system().C_ineq.get());
463  PetscMatrix<T> * cineq_jac = cast_ptr<PetscMatrix<T> *>(this->system().C_ineq_jac.get());
464  PetscVector<T> * lb = cast_ptr<PetscVector<T> *>(&this->system().get_vector("lower_bounds"));
465  PetscVector<T> * ub = cast_ptr<PetscVector<T> *>(&this->system().get_vector("upper_bounds"));
466 
467  // Set the starting guess to zero.
468  x->zero();
469 
470  PetscErrorCode ierr = 0;
471 
472  // Workaround for bug where TaoSetFromOptions *reset*
473  // programmatically set tolerance and max. function evaluation
474  // values when "-tao_type ipm" was specified on the command line: we
475  // call TaoSetFromOptions twice (both before and after setting
476  // custom options programmatically)
477  ierr = TaoSetFromOptions(_tao);
478  LIBMESH_CHKERR(ierr);
479 
480  // Set convergence tolerances
481  // f(X) - f(X*) (estimated) <= fatol
482  // |f(X) - f(X*)| (estimated) / |f(X)| <= frtol
483  // ||g(X)|| <= gatol
484  // ||g(X)|| / |f(X)| <= grtol
485  // ||g(X)|| / ||g(X0)|| <= gttol
486  // Command line equivalents: -tao_fatol, -tao_frtol, -tao_gatol, -tao_grtol, -tao_gttol
487  ierr = TaoSetTolerances(_tao,
488 #if PETSC_RELEASE_LESS_THAN(3,7,0)
489  // Releases up to 3.X.Y had fatol and frtol, after that they were removed.
490  // Hopefully we'll be able to know X and Y soon. Guessing at 3.7.0.
491  /*fatol=*/PETSC_DEFAULT,
492  /*frtol=*/PETSC_DEFAULT,
493 #endif
494  /*gatol=*/PETSC_DEFAULT,
496  /*gttol=*/PETSC_DEFAULT);
497  LIBMESH_CHKERR(ierr);
498 
499  // Set the max-allowed number of objective function evaluations
500  // Command line equivalent: -tao_max_funcs
501  ierr = TaoSetMaximumFunctionEvaluations(_tao, this->max_objective_function_evaluations);
502  LIBMESH_CHKERR(ierr);
503 
504  // Set the max-allowed number of optimization iterations.
505  // Command line equivalent: -tao_max_it
506  // Not implemented for now as it seems fairly similar to
507  // ierr = TaoSetMaximumIterations(_tao, 4);
508  // LIBMESH_CHKERR(ierr);
509 
510  // Set solution vec and an initial guess
511  ierr = TaoSetInitialVector(_tao, x->vec());
512  LIBMESH_CHKERR(ierr);
513 
514  // We have to have an objective function
516 
517  // Set routines for objective, gradient, hessian evaluation
518  ierr = TaoSetObjectiveRoutine(_tao, __libmesh_tao_objective, this);
519  LIBMESH_CHKERR(ierr);
520 
521  if (this->gradient_object)
522  {
523  ierr = TaoSetGradientRoutine(_tao, __libmesh_tao_gradient, this);
524  LIBMESH_CHKERR(ierr);
525  }
526 
527  if (this->hessian_object)
528  {
529  ierr = TaoSetHessianRoutine(_tao, hessian->mat(), hessian->mat(), __libmesh_tao_hessian, this);
530  LIBMESH_CHKERR(ierr);
531  }
532 
534  {
535  // Need to actually compute the bounds vectors first
537 
538  ierr = TaoSetVariableBounds(_tao,
539  lb->vec(),
540  ub->vec());
541  LIBMESH_CHKERR(ierr);
542  }
543 
544  if (this->equality_constraints_object)
545  {
546  ierr = TaoSetEqualityConstraintsRoutine(_tao, ceq->vec(), __libmesh_tao_equality_constraints, this);
547  LIBMESH_CHKERR(ierr);
548  }
549 
551  {
552  ierr = TaoSetJacobianEqualityRoutine(_tao,
553  ceq_jac->mat(),
554  ceq_jac->mat(),
556  this);
557  LIBMESH_CHKERR(ierr);
558  }
559 
560  // Optionally set inequality constraints
562  {
563  ierr = TaoSetInequalityConstraintsRoutine(_tao, cineq->vec(), __libmesh_tao_inequality_constraints, this);
564  LIBMESH_CHKERR(ierr);
565  }
566 
567  // Optionally set inequality constraints Jacobian
569  {
570  ierr = TaoSetJacobianInequalityRoutine(_tao,
571  cineq_jac->mat(),
572  cineq_jac->mat(),
574  this);
575  LIBMESH_CHKERR(ierr);
576  }
577 
578  // Check for Tao command line options
579  ierr = TaoSetFromOptions(_tao);
580  LIBMESH_CHKERR(ierr);
581 
582  // Perform the optimization
583  ierr = TaoSolve(_tao);
584  LIBMESH_CHKERR(ierr);
585 
586  // Store the convergence/divergence reason
587  ierr = TaoGetConvergedReason(_tao, &_reason);
588  LIBMESH_CHKERR(ierr);
589 }
590 
591 
592 template <typename T>
594 {
595  LOG_SCOPE("get_dual_variables()", "TaoOptimizationSolver");
596 
597  PetscVector<T> * lambda_eq_petsc =
598  cast_ptr<PetscVector<T> *>(this->system().lambda_eq.get());
599  PetscVector<T> * lambda_ineq_petsc =
600  cast_ptr<PetscVector<T> *>(this->system().lambda_ineq.get());
601 
602  Vec lambda_eq_petsc_vec = lambda_eq_petsc->vec();
603  Vec lambda_ineq_petsc_vec = lambda_ineq_petsc->vec();
604 
605  PetscErrorCode ierr = 0;
606  ierr = TaoGetDualVariables(_tao,
607  &lambda_eq_petsc_vec,
608  &lambda_ineq_petsc_vec);
609  LIBMESH_CHKERR(ierr);
610 }
611 
612 
613 template <typename T>
615 {
616  libMesh::out << "Tao optimization solver convergence/divergence reason: "
617  << TaoConvergedReasons[this->get_converged_reason()] << std::endl;
618 }
619 
620 
621 
622 template <typename T>
624 {
625  PetscErrorCode ierr=0;
626 
627  if (this->initialized())
628  {
629  ierr = TaoGetConvergedReason(_tao, &_reason);
630  LIBMESH_CHKERR(ierr);
631  }
632 
633  return static_cast<int>(_reason);
634 }
635 
636 
637 //------------------------------------------------------------------
638 // Explicit instantiations
639 template class TaoOptimizationSolver<Number>;
640 
641 } // namespace libMesh
642 
643 
644 
645 #endif // #if defined(LIBMESH_HAVE_PETSC_TAO) && !defined(LIBMESH_USE_COMPLEX_NUMBERS)
OptimizationSystem::ComputeObjective * objective_object
Object that computes the objective function f(X) at the input iterate X.
friend PetscErrorCode __libmesh_tao_inequality_constraints(Tao tao, Vec x, Vec cineq, void *ctx)
This class provides an interface to the Tao optimization solvers.
virtual void gradient(const NumericVector< Number > &X, NumericVector< Number > &grad_f, sys_type &S)=0
This function will be called to compute the gradient of the objective function, and must be implement...
friend PetscErrorCode __libmesh_tao_hessian(Tao tao, Vec x, Mat h, Mat pc, void *ctx)
UniquePtr< NumericVector< Number > > lambda_ineq
ImplicitSystem & sys
Real objective_function_relative_tolerance
Required change in objective function which signals convergence.
virtual void hessian(const NumericVector< Number > &X, SparseMatrix< Number > &H_f, sys_type &S)=0
This function will be called to compute the Hessian of the objective function, and must be implemente...
virtual void init() libmesh_override
Initialize data structures if not done so already.
virtual void inequality_constraints_jacobian(const NumericVector< Number > &X, SparseMatrix< Number > &C_ineq_jac, sys_type &S)=0
This function will be called to evaluate the Jacobian of C_ineq(X).
const class libmesh_nullptr_t libmesh_nullptr
friend PetscErrorCode __libmesh_tao_equality_constraints(Tao tao, Vec x, Vec ce, void *ctx)
OptimizationSystem::ComputeGradient * gradient_object
Object that computes the gradient grad_f(X) of the objective function at the input iterate X...
The libMesh namespace provides an interface to certain functionality in the library.
PetscErrorCode Vec Mat Mat pc
bool _is_initialized
Flag indicating if the data structures have been initialized.
const sys_type & system() const
UniquePtr< NumericVector< Number > > lambda_eq
Vectors to store the dual variables associated with equality and inequality constraints.
OptimizationSystem::ComputeHessian * hessian_object
Object that computes the Hessian H_f(X) of the objective function at the input iterate X...
virtual void clear() libmesh_override
Release all memory and clear data structures.
friend PetscErrorCode __libmesh_tao_equality_constraints_jacobian(Tao tao, Vec x, Mat J, Mat Jpre, void *ctx)
libmesh_assert(j)
PetscDiffSolver & solver
OptimizationSystem::ComputeInequalityConstraints * inequality_constraints_object
Object that computes the inequality constraints vector C_ineq(X).
PetscErrorCode __libmesh_tao_gradient(Tao tao, Vec x, Vec g, void *ctx)
PetscErrorCode __libmesh_tao_hessian(Tao tao, Vec x, Mat h, Mat pc, void *ctx)
PetscErrorCode Vec x
OptimizationSystem::ComputeEqualityConstraintsJacobian * equality_constraints_jacobian_object
Object that computes the Jacobian of C_eq(X).
Tao _tao
Optimization solver context.
TaoOptimizationSolver(sys_type &system)
Constructor.
friend PetscErrorCode __libmesh_tao_inequality_constraints_jacobian(Tao tao, Vec x, Mat J, Mat Jpre, void *ctx)
This System subclass enables us to assemble an objective function, gradient, Hessian and bounds for o...
friend PetscErrorCode __libmesh_tao_objective(Tao tao, Vec x, PetscReal *objective, void *ctx)
unsigned int max_objective_function_evaluations
Maximum number of objective function evaluations allowed.
UniquePtr< SparseMatrix< Number > > C_eq_jac
The sparse matrix that stores the Jacobian of C_eq.
friend PetscErrorCode __libmesh_tao_gradient(Tao tao, Vec x, Vec g, void *ctx)
PetscErrorCode Vec Mat Mat void * ctx
UniquePtr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1523
virtual int get_converged_reason() libmesh_override
const NumericVector< Number > & get_vector(const std::string &vec_name) const
Definition: system.C:794
virtual void equality_constraints(const NumericVector< Number > &X, NumericVector< Number > &C_eq, sys_type &S)=0
This function will be called to evaluate the equality constraints vector C_eq(X). ...
PetscErrorCode __libmesh_tao_objective(Tao tao, Vec x, PetscReal *objective, void *ctx)
virtual void lower_and_upper_bounds(sys_type &S)=0
This function should update the following two vectors: this->get_vector("lower_bounds"), this->get_vector("upper_bounds").
virtual void get_dual_variables() libmesh_override
Get the current values of dual variables associated with inequality and equality constraints.
virtual Number objective(const NumericVector< Number > &X, sys_type &S)=0
This function will be called to compute the objective function to be minimized, and must be implement...
virtual void equality_constraints_jacobian(const NumericVector< Number > &X, SparseMatrix< Number > &C_eq_jac, sys_type &S)=0
This function will be called to evaluate the Jacobian of C_eq(X).
TaoConvergedReason _reason
Store the reason for Tao convergence/divergence for use even after _tao has been cleared.
OptimizationSystem::ComputeLowerAndUpperBounds * lower_and_upper_bounds_object
Object that computes the lower and upper bounds vectors.
This base class can be inherited from to provide interfaces to optimization solvers from different pa...
PetscErrorCode ierr
const Parallel::Communicator & comm() const
OStreamProxy out
SparseMatrix< Number > * matrix
The system matrix.
virtual void solve() libmesh_override
Call the Tao solver.
PetscErrorCode __libmesh_tao_equality_constraints_jacobian(Tao tao, Vec x, Mat J, Mat Jpre, void *ctx)
virtual void print_converged_reason() libmesh_override
Prints a useful message about why the latest optimization solve con(di)verged.
PetscErrorCode __libmesh_tao_equality_constraints(Tao tao, Vec x, Vec ce, void *ctx)
PetscErrorCode __libmesh_tao_inequality_constraints_jacobian(Tao tao, Vec x, Mat J, Mat Jpre, void *ctx)
UniquePtr< NumericVector< Number > > C_eq
The vector that stores equality constraints.
PetscErrorCode __libmesh_tao_inequality_constraints(Tao tao, Vec x, Vec cineq, void *ctx)
virtual void inequality_constraints(const NumericVector< Number > &X, NumericVector< Number > &C_ineq, sys_type &S)=0
This function will be called to evaluate the equality constraints vector C_ineq(X).
OptimizationSystem::ComputeInequalityConstraintsJacobian * inequality_constraints_jacobian_object
Object that computes the Jacobian of C_ineq(X).
UniquePtr< SparseMatrix< Number > > C_ineq_jac
The sparse matrix that stores the Jacobian of C_ineq.
OptimizationSystem::ComputeEqualityConstraints * equality_constraints_object
Object that computes the equality constraints vector C_eq(X).
UniquePtr< NumericVector< Number > > C_ineq
The vector that stores inequality constraints.