tao_optimization_solver.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2018 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"
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 the input vector
71  // "x", update sys.current_local_solution based on "x", then swap
72  // back.
73  X.swap(X_sys);
74  sys.update();
75  X.swap(X_sys);
76 
77  // Enforce constraints (if any) exactly on the
78  // current_local_solution. This is the solution vector that is
79  // actually used in the computation of the objective function
80  // below, and is not locked by debug-enabled PETSc the way that
81  // the solution vector is.
82  sys.get_dof_map().enforce_constraints_exactly(sys, sys.current_local_solution.get());
83 
84  if (solver->objective_object != nullptr)
85  (*objective) = solver->objective_object->objective(*(sys.current_local_solution), sys);
86  else
87  libmesh_error_msg("Objective function not defined in __libmesh_tao_objective");
88 
89  return ierr;
90  }
91 
92 
93 
94  //---------------------------------------------------------------
95  // This function is called by Tao to evaluate the gradient at x
96  PetscErrorCode
97  __libmesh_tao_gradient(Tao /*tao*/, Vec x, Vec g, void * ctx)
98  {
99  LOG_SCOPE("gradient()", "TaoOptimizationSolver");
100 
101  PetscErrorCode ierr = 0;
102 
103  libmesh_assert(x);
104  libmesh_assert(g);
105  libmesh_assert(ctx);
106 
107  // ctx should be a pointer to the solver (it was passed in as void *)
109  static_cast<TaoOptimizationSolver<Number> *> (ctx);
110 
111  OptimizationSystem & sys = solver->system();
112 
113  // We'll use current_local_solution below, so let's ensure that it's consistent
114  // with the vector x that was passed in.
115  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
116  PetscVector<Number> X(x, sys.comm());
117 
118  // Perform a swap so that sys.solution points to the input vector
119  // "x", update sys.current_local_solution based on "x", then swap
120  // back.
121  X.swap(X_sys);
122  sys.update();
123  X.swap(X_sys);
124 
125  // We'll also pass the gradient in to the assembly routine
126  // so let's make a PETSc vector for that too.
127  PetscVector<Number> gradient(g, sys.comm());
128 
129  // Clear the gradient prior to assembly
130  gradient.zero();
131 
132  // Enforce constraints exactly on the current_local_solution.
133  sys.get_dof_map().enforce_constraints_exactly(sys, sys.current_local_solution.get());
134 
135  if (solver->gradient_object != nullptr)
136  solver->gradient_object->gradient(*(sys.current_local_solution), gradient, sys);
137  else
138  libmesh_error_msg("Gradient function not defined in __libmesh_tao_gradient");
139 
140  gradient.close();
141 
142  return ierr;
143  }
144 
145  //---------------------------------------------------------------
146  // This function is called by Tao to evaluate the Hessian at x
147  PetscErrorCode
148  __libmesh_tao_hessian(Tao /*tao*/, Vec x, Mat h, Mat pc, void * ctx)
149  {
150  LOG_SCOPE("hessian()", "TaoOptimizationSolver");
151 
152  PetscErrorCode ierr = 0;
153 
154  libmesh_assert(x);
155  libmesh_assert(h);
156  libmesh_assert(pc);
157  libmesh_assert(ctx);
158 
159  // ctx should be a pointer to the solver (it was passed in as void *)
161  static_cast<TaoOptimizationSolver<Number> *> (ctx);
162 
163  OptimizationSystem & sys = solver->system();
164 
165  // We'll use current_local_solution below, so let's ensure that it's consistent
166  // with the vector x that was passed in.
167  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
168  PetscVector<Number> X(x, sys.comm());
169 
170  // Perform a swap so that sys.solution points to the input vector
171  // "x", update sys.current_local_solution based on "x", then swap
172  // back.
173  X.swap(X_sys);
174  sys.update();
175  X.swap(X_sys);
176 
177  // Let's also wrap pc and h in PetscMatrix objects for convenience
178  PetscMatrix<Number> PC(pc, sys.comm());
179  PetscMatrix<Number> hessian(h, sys.comm());
180  PC.attach_dof_map(sys.get_dof_map());
181  hessian.attach_dof_map(sys.get_dof_map());
182 
183  // Enforce constraints exactly on the current_local_solution.
184  sys.get_dof_map().enforce_constraints_exactly(sys, sys.current_local_solution.get());
185 
186  if (solver->hessian_object != nullptr)
187  {
188  // Following PetscNonlinearSolver by passing in PC. It's not clear
189  // why we pass in PC and not hessian though?
190  solver->hessian_object->hessian(*(sys.current_local_solution), PC, sys);
191  }
192  else
193  libmesh_error_msg("Hessian function not defined in __libmesh_tao_hessian");
194 
195  PC.close();
196  hessian.close();
197 
198  return ierr;
199  }
200 
201 
202  //---------------------------------------------------------------
203  // This function is called by Tao to evaluate the equality constraints at x
204  PetscErrorCode
205  __libmesh_tao_equality_constraints(Tao /*tao*/, Vec x, Vec ce, void * ctx)
206  {
207  LOG_SCOPE("equality_constraints()", "TaoOptimizationSolver");
208 
209  PetscErrorCode ierr = 0;
210 
211  libmesh_assert(x);
212  libmesh_assert(ce);
213  libmesh_assert(ctx);
214 
215  // ctx should be a pointer to the solver (it was passed in as void *)
217  static_cast<TaoOptimizationSolver<Number> *> (ctx);
218 
219  OptimizationSystem & sys = solver->system();
220 
221  // We'll use current_local_solution below, so let's ensure that it's consistent
222  // with the vector x that was passed in.
223  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
224  PetscVector<Number> X(x, sys.comm());
225 
226  // Perform a swap so that sys.solution points to the input vector
227  // "x", update sys.current_local_solution based on "x", then swap
228  // back.
229  X.swap(X_sys);
230  sys.update();
231  X.swap(X_sys);
232 
233  // We'll also pass the constraints vector ce into the assembly routine
234  // so let's make a PETSc vector for that too.
235  PetscVector<Number> eq_constraints(ce, sys.comm());
236 
237  // Clear the gradient prior to assembly
238  eq_constraints.zero();
239 
240  // Enforce constraints exactly on the current_local_solution.
241  sys.get_dof_map().enforce_constraints_exactly(sys, sys.current_local_solution.get());
242 
243  if (solver->equality_constraints_object != nullptr)
244  solver->equality_constraints_object->equality_constraints(*(sys.current_local_solution), eq_constraints, sys);
245  else
246  libmesh_error_msg("Constraints function not defined in __libmesh_tao_equality_constraints");
247 
248  eq_constraints.close();
249 
250  return ierr;
251  }
252 
253  //---------------------------------------------------------------
254  // This function is called by Tao to evaluate the Jacobian of the
255  // equality constraints at x
256  PetscErrorCode
257  __libmesh_tao_equality_constraints_jacobian(Tao /*tao*/, Vec x, Mat J, Mat Jpre, void * ctx)
258  {
259  LOG_SCOPE("equality_constraints_jacobian()", "TaoOptimizationSolver");
260 
261  PetscErrorCode ierr = 0;
262 
263  libmesh_assert(x);
264  libmesh_assert(J);
265  libmesh_assert(Jpre);
266 
267  // ctx should be a pointer to the solver (it was passed in as void *)
269  static_cast<TaoOptimizationSolver<Number> *> (ctx);
270 
271  OptimizationSystem & sys = solver->system();
272 
273  // We'll use current_local_solution below, so let's ensure that it's consistent
274  // with the vector x that was passed in.
275  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
276  PetscVector<Number> X(x, sys.comm());
277 
278  // Perform a swap so that sys.solution points to the input vector
279  // "x", update sys.current_local_solution based on "x", then swap
280  // back.
281  X.swap(X_sys);
282  sys.update();
283  X.swap(X_sys);
284 
285  // Let's also wrap J and Jpre in PetscMatrix objects for convenience
286  PetscMatrix<Number> J_petsc(J, sys.comm());
287  PetscMatrix<Number> Jpre_petsc(Jpre, sys.comm());
288 
289  // Enforce constraints exactly on the current_local_solution.
290  sys.get_dof_map().enforce_constraints_exactly(sys, sys.current_local_solution.get());
291 
292  if (solver->equality_constraints_jacobian_object != nullptr)
293  solver->equality_constraints_jacobian_object->equality_constraints_jacobian(*(sys.current_local_solution), J_petsc, sys);
294  else
295  libmesh_error_msg("Constraints function not defined in __libmesh_tao_equality_constraints_jacobian");
296 
297  J_petsc.close();
298  Jpre_petsc.close();
299 
300  return ierr;
301  }
302 
303  //---------------------------------------------------------------
304  // This function is called by Tao to evaluate the inequality constraints at x
305  PetscErrorCode
306  __libmesh_tao_inequality_constraints(Tao /*tao*/, Vec x, Vec cineq, void * ctx)
307  {
308  LOG_SCOPE("inequality_constraints()", "TaoOptimizationSolver");
309 
310  PetscErrorCode ierr = 0;
311 
312  libmesh_assert(x);
313  libmesh_assert(cineq);
314  libmesh_assert(ctx);
315 
316  // ctx should be a pointer to the solver (it was passed in as void *)
318  static_cast<TaoOptimizationSolver<Number> *> (ctx);
319 
320  OptimizationSystem & sys = solver->system();
321 
322  // We'll use current_local_solution below, so let's ensure that it's consistent
323  // with the vector x that was passed in.
324  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
325  PetscVector<Number> X(x, sys.comm());
326 
327  // Perform a swap so that sys.solution points to the input vector
328  // "x", update sys.current_local_solution based on "x", then swap
329  // back.
330  X.swap(X_sys);
331  sys.update();
332  X.swap(X_sys);
333 
334  // We'll also pass the constraints vector ce into the assembly routine
335  // so let's make a PETSc vector for that too.
336  PetscVector<Number> ineq_constraints(cineq, sys.comm());
337 
338  // Clear the gradient prior to assembly
339  ineq_constraints.zero();
340 
341  // Enforce constraints exactly on the current_local_solution.
342  sys.get_dof_map().enforce_constraints_exactly(sys, sys.current_local_solution.get());
343 
344  if (solver->inequality_constraints_object != nullptr)
345  solver->inequality_constraints_object->inequality_constraints(*(sys.current_local_solution), ineq_constraints, sys);
346  else
347  libmesh_error_msg("Constraints function not defined in __libmesh_tao_inequality_constraints");
348 
349  ineq_constraints.close();
350 
351  return ierr;
352  }
353 
354  //---------------------------------------------------------------
355  // This function is called by Tao to evaluate the Jacobian of the
356  // equality constraints at x
357  PetscErrorCode
358  __libmesh_tao_inequality_constraints_jacobian(Tao /*tao*/, Vec x, Mat J, Mat Jpre, void * ctx)
359  {
360  LOG_SCOPE("inequality_constraints_jacobian()", "TaoOptimizationSolver");
361 
362  PetscErrorCode ierr = 0;
363 
364  libmesh_assert(x);
365  libmesh_assert(J);
366  libmesh_assert(Jpre);
367 
368  // ctx should be a pointer to the solver (it was passed in as void *)
370  static_cast<TaoOptimizationSolver<Number> *> (ctx);
371 
372  OptimizationSystem & sys = solver->system();
373 
374  // We'll use current_local_solution below, so let's ensure that it's consistent
375  // with the vector x that was passed in.
376  PetscVector<Number> & X_sys = *cast_ptr<PetscVector<Number> *>(sys.solution.get());
377  PetscVector<Number> X(x, sys.comm());
378 
379  // Perform a swap so that sys.solution points to the input vector
380  // "x", update sys.current_local_solution based on "x", then swap
381  // back.
382  X.swap(X_sys);
383  sys.update();
384  X.swap(X_sys);
385 
386  // Let's also wrap J and Jpre in PetscMatrix objects for convenience
387  PetscMatrix<Number> J_petsc(J, sys.comm());
388  PetscMatrix<Number> Jpre_petsc(Jpre, sys.comm());
389 
390  // Enforce constraints exactly on the current_local_solution.
391  sys.get_dof_map().enforce_constraints_exactly(sys, sys.current_local_solution.get());
392 
393  if (solver->inequality_constraints_jacobian_object != nullptr)
394  solver->inequality_constraints_jacobian_object->inequality_constraints_jacobian(*(sys.current_local_solution), J_petsc, sys);
395  else
396  libmesh_error_msg("Constraints function not defined in __libmesh_tao_inequality_constraints_jacobian");
397 
398  J_petsc.close();
399  Jpre_petsc.close();
400 
401  return ierr;
402  }
403 
404 } // end extern "C"
405 //---------------------------------------------------------------------
406 
407 
408 
409 //---------------------------------------------------------------------
410 // TaoOptimizationSolver<> methods
411 template <typename T>
413  OptimizationSolver<T>(system_in),
414  _reason(TAO_CONVERGED_USER) // Arbitrary initial value...
415 {
416 }
417 
418 
419 
420 template <typename T>
422 {
423  this->clear ();
424 }
425 
426 
427 
428 template <typename T>
430 {
431  if (this->initialized())
432  {
433  this->_is_initialized = false;
434 
435  PetscErrorCode ierr=0;
436 
437  ierr = TaoDestroy(&_tao);
438  LIBMESH_CHKERR(ierr);
439  }
440 }
441 
442 
443 
444 template <typename T>
446 {
447  // Initialize the data structures if not done so already.
448  if (!this->initialized())
449  {
450  this->_is_initialized = true;
451 
452  PetscErrorCode ierr=0;
453 
454  ierr = TaoCreate(this->comm().get(),&_tao);
455  LIBMESH_CHKERR(ierr);
456  }
457 }
458 
459 template <typename T>
461 {
462  LOG_SCOPE("solve()", "TaoOptimizationSolver");
463 
464  this->init ();
465 
466  this->system().solution->zero();
467 
468  PetscMatrix<T> * hessian = cast_ptr<PetscMatrix<T> *>(this->system().matrix);
469  // PetscVector<T> * gradient = cast_ptr<PetscVector<T> *>(this->system().rhs);
470  PetscVector<T> * x = cast_ptr<PetscVector<T> *>(this->system().solution.get());
471  PetscVector<T> * ceq = cast_ptr<PetscVector<T> *>(this->system().C_eq.get());
472  PetscMatrix<T> * ceq_jac = cast_ptr<PetscMatrix<T> *>(this->system().C_eq_jac.get());
473  PetscVector<T> * cineq = cast_ptr<PetscVector<T> *>(this->system().C_ineq.get());
474  PetscMatrix<T> * cineq_jac = cast_ptr<PetscMatrix<T> *>(this->system().C_ineq_jac.get());
475  PetscVector<T> * lb = cast_ptr<PetscVector<T> *>(&this->system().get_vector("lower_bounds"));
476  PetscVector<T> * ub = cast_ptr<PetscVector<T> *>(&this->system().get_vector("upper_bounds"));
477 
478  // Set the starting guess to zero.
479  x->zero();
480 
481  PetscErrorCode ierr = 0;
482 
483  // Workaround for bug where TaoSetFromOptions *reset*
484  // programmatically set tolerance and max. function evaluation
485  // values when "-tao_type ipm" was specified on the command line: we
486  // call TaoSetFromOptions twice (both before and after setting
487  // custom options programmatically)
488  ierr = TaoSetFromOptions(_tao);
489  LIBMESH_CHKERR(ierr);
490 
491  // Set convergence tolerances
492  // f(X) - f(X*) (estimated) <= fatol
493  // |f(X) - f(X*)| (estimated) / |f(X)| <= frtol
494  // ||g(X)|| <= gatol
495  // ||g(X)|| / |f(X)| <= grtol
496  // ||g(X)|| / ||g(X0)|| <= gttol
497  // Command line equivalents: -tao_fatol, -tao_frtol, -tao_gatol, -tao_grtol, -tao_gttol
498  ierr = TaoSetTolerances(_tao,
499 #if PETSC_RELEASE_LESS_THAN(3,7,0)
500  // Releases up to 3.X.Y had fatol and frtol, after that they were removed.
501  // Hopefully we'll be able to know X and Y soon. Guessing at 3.7.0.
502  /*fatol=*/PETSC_DEFAULT,
503  /*frtol=*/PETSC_DEFAULT,
504 #endif
505  /*gatol=*/PETSC_DEFAULT,
506  /*grtol=*/this->objective_function_relative_tolerance,
507  /*gttol=*/PETSC_DEFAULT);
508  LIBMESH_CHKERR(ierr);
509 
510  // Set the max-allowed number of objective function evaluations
511  // Command line equivalent: -tao_max_funcs
512  ierr = TaoSetMaximumFunctionEvaluations(_tao, this->max_objective_function_evaluations);
513  LIBMESH_CHKERR(ierr);
514 
515  // Set the max-allowed number of optimization iterations.
516  // Command line equivalent: -tao_max_it
517  // Not implemented for now as it seems fairly similar to
518  // ierr = TaoSetMaximumIterations(_tao, 4);
519  // LIBMESH_CHKERR(ierr);
520 
521  // Set solution vec and an initial guess
522  ierr = TaoSetInitialVector(_tao, x->vec());
523  LIBMESH_CHKERR(ierr);
524 
525  // We have to have an objective function
526  libmesh_assert( this->objective_object );
527 
528  // Set routines for objective, gradient, hessian evaluation
529  ierr = TaoSetObjectiveRoutine(_tao, __libmesh_tao_objective, this);
530  LIBMESH_CHKERR(ierr);
531 
532  if (this->gradient_object)
533  {
534  ierr = TaoSetGradientRoutine(_tao, __libmesh_tao_gradient, this);
535  LIBMESH_CHKERR(ierr);
536  }
537 
538  if (this->hessian_object)
539  {
540  ierr = TaoSetHessianRoutine(_tao, hessian->mat(), hessian->mat(), __libmesh_tao_hessian, this);
541  LIBMESH_CHKERR(ierr);
542  }
543 
544  if (this->lower_and_upper_bounds_object)
545  {
546  // Need to actually compute the bounds vectors first
547  this->lower_and_upper_bounds_object->lower_and_upper_bounds(this->system());
548 
549  ierr = TaoSetVariableBounds(_tao,
550  lb->vec(),
551  ub->vec());
552  LIBMESH_CHKERR(ierr);
553  }
554 
555  if (this->equality_constraints_object)
556  {
557  ierr = TaoSetEqualityConstraintsRoutine(_tao, ceq->vec(), __libmesh_tao_equality_constraints, this);
558  LIBMESH_CHKERR(ierr);
559  }
560 
561  if (this->equality_constraints_jacobian_object)
562  {
563  ierr = TaoSetJacobianEqualityRoutine(_tao,
564  ceq_jac->mat(),
565  ceq_jac->mat(),
567  this);
568  LIBMESH_CHKERR(ierr);
569  }
570 
571  // Optionally set inequality constraints
572  if (this->inequality_constraints_object)
573  {
574  ierr = TaoSetInequalityConstraintsRoutine(_tao, cineq->vec(), __libmesh_tao_inequality_constraints, this);
575  LIBMESH_CHKERR(ierr);
576  }
577 
578  // Optionally set inequality constraints Jacobian
579  if (this->inequality_constraints_jacobian_object)
580  {
581  ierr = TaoSetJacobianInequalityRoutine(_tao,
582  cineq_jac->mat(),
583  cineq_jac->mat(),
585  this);
586  LIBMESH_CHKERR(ierr);
587  }
588 
589  // Check for Tao command line options
590  ierr = TaoSetFromOptions(_tao);
591  LIBMESH_CHKERR(ierr);
592 
593  // Perform the optimization
594  ierr = TaoSolve(_tao);
595  LIBMESH_CHKERR(ierr);
596 
597  // Enforce constraints exactly now that the solve is done. We have
598  // been enforcing them on the current_local_solution during the
599  // solve, but now need to be sure they are enforced on the parallel
600  // solution vector as well.
601  this->system().get_dof_map().enforce_constraints_exactly(this->system());
602 
603  // Store the convergence/divergence reason
604  ierr = TaoGetConvergedReason(_tao, &_reason);
605  LIBMESH_CHKERR(ierr);
606 }
607 
608 
609 template <typename T>
611 {
612  LOG_SCOPE("get_dual_variables()", "TaoOptimizationSolver");
613 
614  PetscVector<T> * lambda_eq_petsc =
615  cast_ptr<PetscVector<T> *>(this->system().lambda_eq.get());
616  PetscVector<T> * lambda_ineq_petsc =
617  cast_ptr<PetscVector<T> *>(this->system().lambda_ineq.get());
618 
619  Vec lambda_eq_petsc_vec = lambda_eq_petsc->vec();
620  Vec lambda_ineq_petsc_vec = lambda_ineq_petsc->vec();
621 
622  PetscErrorCode ierr = 0;
623  ierr = TaoGetDualVariables(_tao,
624  &lambda_eq_petsc_vec,
625  &lambda_ineq_petsc_vec);
626  LIBMESH_CHKERR(ierr);
627 }
628 
629 
630 template <typename T>
632 {
633  libMesh::out << "Tao optimization solver convergence/divergence reason: "
634  << TaoConvergedReasons[this->get_converged_reason()] << std::endl;
635 }
636 
637 
638 
639 template <typename T>
641 {
642  PetscErrorCode ierr=0;
643 
644  if (this->initialized())
645  {
646  ierr = TaoGetConvergedReason(_tao, &_reason);
647  LIBMESH_CHKERR(ierr);
648  }
649 
650  return static_cast<int>(_reason);
651 }
652 
653 
654 //------------------------------------------------------------------
655 // Explicit instantiations
656 template class TaoOptimizationSolver<Number>;
657 
658 } // namespace libMesh
659 
660 
661 
662 #endif // #if defined(LIBMESH_HAVE_PETSC_TAO) && !defined(LIBMESH_USE_COMPLEX_NUMBERS)
PetscErrorCode __libmesh_tao_objective(Tao, Vec x, PetscReal *objective, void *ctx)
OptimizationSystem::ComputeObjective * objective_object
PetscErrorCode __libmesh_tao_equality_constraints_jacobian(Tao, Vec x, Mat J, Mat Jpre, void *ctx)
NumericVector interface to PETSc Vec.
Definition: petsc_vector.h:64
virtual void gradient(const NumericVector< Number > &X, NumericVector< Number > &grad_f, sys_type &S)=0
virtual void zero() override
Definition: petsc_vector.h:862
virtual void hessian(const NumericVector< Number > &X, SparseMatrix< Number > &H_f, sys_type &S)=0
virtual void inequality_constraints_jacobian(const NumericVector< Number > &X, SparseMatrix< Number > &C_ineq_jac, sys_type &S)=0
const Parallel::Communicator & comm() const
OptimizationSystem::ComputeGradient * gradient_object
virtual int get_converged_reason() override
PetscErrorCode __libmesh_tao_gradient(Tao, Vec x, Vec g, void *ctx)
virtual void print_converged_reason() override
OptimizationSystem::ComputeHessian * hessian_object
virtual void swap(NumericVector< T > &v) override
OptimizationSystem::ComputeInequalityConstraints * inequality_constraints_object
OptimizationSystem::ComputeEqualityConstraintsJacobian * equality_constraints_jacobian_object
void init(triangulateio &t)
PetscErrorCode __libmesh_tao_inequality_constraints_jacobian(Tao, Vec x, Mat J, Mat Jpre, void *ctx)
virtual void equality_constraints(const NumericVector< Number > &X, NumericVector< Number > &C_eq, sys_type &S)=0
const sys_type & system() const
void attach_dof_map(const DofMap &dof_map)
virtual Number objective(const NumericVector< Number > &X, sys_type &S)=0
virtual void equality_constraints_jacobian(const NumericVector< Number > &X, SparseMatrix< Number > &C_eq_jac, sys_type &S)=0
PetscErrorCode ierr
SparseMatrix interface to PETSc Mat.
bool initialized()
Definition: libmesh.C:258
PetscErrorCode __libmesh_tao_hessian(Tao, Vec x, Mat h, Mat pc, void *ctx)
PetscErrorCode __libmesh_tao_equality_constraints(Tao, Vec x, Vec ce, void *ctx)
OStreamProxy out(std::cout)
PetscErrorCode __libmesh_tao_inequality_constraints(Tao, Vec x, Vec cineq, void *ctx)
virtual void get_dual_variables() override
virtual void inequality_constraints(const NumericVector< Number > &X, NumericVector< Number > &C_ineq, sys_type &S)=0
OptimizationSystem::ComputeInequalityConstraintsJacobian * inequality_constraints_jacobian_object
OptimizationSystem::ComputeEqualityConstraints * equality_constraints_object