22 #if defined(LIBMESH_HAVE_PETSC_TAO) && !defined(LIBMESH_USE_COMPLEX_NUMBERS) 51 LOG_SCOPE(
"objective()",
"TaoOptimizationSolver");
53 PetscErrorCode
ierr = 0;
56 libmesh_assert(objective);
73 sys.get_dof_map().enforce_constraints_exactly(sys);
82 libmesh_error_msg(
"Objective function not defined in __libmesh_tao_objective");
94 LOG_SCOPE(
"gradient()",
"TaoOptimizationSolver");
96 PetscErrorCode
ierr = 0;
116 sys.get_dof_map().enforce_constraints_exactly(sys);
132 libmesh_error_msg(
"Gradient function not defined in __libmesh_tao_gradient");
144 LOG_SCOPE(
"hessian()",
"TaoOptimizationSolver");
146 PetscErrorCode
ierr = 0;
167 sys.get_dof_map().enforce_constraints_exactly(sys);
177 hessian.attach_dof_map(sys.get_dof_map());
186 libmesh_error_msg(
"Hessian function not defined in __libmesh_tao_hessian");
200 LOG_SCOPE(
"equality_constraints()",
"TaoOptimizationSolver");
202 PetscErrorCode
ierr = 0;
222 sys.get_dof_map().enforce_constraints_exactly(sys);
233 eq_constraints.
zero();
238 libmesh_error_msg(
"Constraints function not defined in __libmesh_tao_equality_constraints");
240 eq_constraints.close();
251 LOG_SCOPE(
"equality_constraints_jacobian()",
"TaoOptimizationSolver");
253 PetscErrorCode
ierr = 0;
257 libmesh_assert(Jpre);
273 sys.get_dof_map().enforce_constraints_exactly(sys);
286 libmesh_error_msg(
"Constraints function not defined in __libmesh_tao_equality_constraints_jacobian");
299 LOG_SCOPE(
"inequality_constraints()",
"TaoOptimizationSolver");
301 PetscErrorCode
ierr = 0;
304 libmesh_assert(cineq);
321 sys.get_dof_map().enforce_constraints_exactly(sys);
332 ineq_constraints.
zero();
337 libmesh_error_msg(
"Constraints function not defined in __libmesh_tao_inequality_constraints");
339 ineq_constraints.close();
350 LOG_SCOPE(
"inequality_constraints_jacobian()",
"TaoOptimizationSolver");
352 PetscErrorCode
ierr = 0;
356 libmesh_assert(Jpre);
372 sys.get_dof_map().enforce_constraints_exactly(sys);
385 libmesh_error_msg(
"Constraints function not defined in __libmesh_tao_inequality_constraints_jacobian");
400 template <
typename T>
403 _reason(TAO_CONVERGED_USER)
409 template <
typename T>
417 template <
typename T>
424 PetscErrorCode
ierr=0;
426 ierr = TaoDestroy(&
_tao);
427 LIBMESH_CHKERR(ierr);
433 template <
typename T>
441 PetscErrorCode
ierr=0;
443 ierr = TaoCreate(this->
comm().
get(),&
_tao);
444 LIBMESH_CHKERR(ierr);
448 template <
typename T>
451 LOG_SCOPE(
"solve()",
"TaoOptimizationSolver");
470 PetscErrorCode
ierr = 0;
477 ierr = TaoSetFromOptions(
_tao);
478 LIBMESH_CHKERR(ierr);
487 ierr = TaoSetTolerances(
_tao,
488 #
if PETSC_RELEASE_LESS_THAN(3,7,0)
497 LIBMESH_CHKERR(ierr);
502 LIBMESH_CHKERR(ierr);
511 ierr = TaoSetInitialVector(
_tao, x->vec());
512 LIBMESH_CHKERR(ierr);
519 LIBMESH_CHKERR(ierr);
524 LIBMESH_CHKERR(ierr);
530 LIBMESH_CHKERR(ierr);
538 ierr = TaoSetVariableBounds(
_tao,
541 LIBMESH_CHKERR(ierr);
547 LIBMESH_CHKERR(ierr);
552 ierr = TaoSetJacobianEqualityRoutine(
_tao,
557 LIBMESH_CHKERR(ierr);
564 LIBMESH_CHKERR(ierr);
570 ierr = TaoSetJacobianInequalityRoutine(
_tao,
575 LIBMESH_CHKERR(ierr);
579 ierr = TaoSetFromOptions(
_tao);
580 LIBMESH_CHKERR(ierr);
583 ierr = TaoSolve(
_tao);
584 LIBMESH_CHKERR(ierr);
588 LIBMESH_CHKERR(ierr);
592 template <
typename T>
595 LOG_SCOPE(
"get_dual_variables()",
"TaoOptimizationSolver");
602 Vec lambda_eq_petsc_vec = lambda_eq_petsc->
vec();
603 Vec lambda_ineq_petsc_vec = lambda_ineq_petsc->vec();
605 PetscErrorCode
ierr = 0;
606 ierr = TaoGetDualVariables(
_tao,
607 &lambda_eq_petsc_vec,
608 &lambda_ineq_petsc_vec);
609 LIBMESH_CHKERR(ierr);
613 template <
typename T>
616 libMesh::out <<
"Tao optimization solver convergence/divergence reason: " 622 template <
typename T>
625 PetscErrorCode
ierr=0;
630 LIBMESH_CHKERR(ierr);
633 return static_cast<int>(
_reason);
645 #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
virtual void zero() libmesh_override
friend PetscErrorCode __libmesh_tao_inequality_constraints(Tao tao, Vec x, Vec cineq, void *ctx)
PetscErrorCode __libmesh_tao_equality_constraints_jacobian(Tao, Vec x, Mat J, Mat Jpre, void *ctx)
NumericVector interface to PETSc Vec.
virtual void gradient(const NumericVector< Number > &X, NumericVector< Number > &grad_f, sys_type &S)=0
friend PetscErrorCode __libmesh_tao_hessian(Tao tao, Vec x, Mat h, Mat pc, void *ctx)
std::unique_ptr< NumericVector< Number > > lambda_ineq
Real objective_function_relative_tolerance
virtual void hessian(const NumericVector< Number > &X, SparseMatrix< Number > &H_f, sys_type &S)=0
virtual void init() libmesh_override
virtual void inequality_constraints_jacobian(const NumericVector< Number > &X, SparseMatrix< Number > &C_ineq_jac, sys_type &S)=0
std::unique_ptr< NumericVector< Number > > C_eq
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
PetscErrorCode __libmesh_tao_gradient(Tao, Vec x, Vec g, void *ctx)
std::unique_ptr< SparseMatrix< Number > > C_ineq_jac
const sys_type & system() const
OptimizationSystem::ComputeHessian * hessian_object
virtual void clear() libmesh_override
friend PetscErrorCode __libmesh_tao_equality_constraints_jacobian(Tao tao, Vec x, Mat J, Mat Jpre, void *ctx)
OptimizationSystem::ComputeInequalityConstraints * inequality_constraints_object
OptimizationSystem::ComputeEqualityConstraintsJacobian * equality_constraints_jacobian_object
TaoOptimizationSolver(sys_type &system)
friend PetscErrorCode __libmesh_tao_inequality_constraints_jacobian(Tao tao, Vec x, Mat J, Mat Jpre, void *ctx)
std::unique_ptr< NumericVector< Number > > solution
friend PetscErrorCode __libmesh_tao_objective(Tao tao, Vec x, PetscReal *objective, void *ctx)
unsigned int max_objective_function_evaluations
friend PetscErrorCode __libmesh_tao_gradient(Tao tao, Vec x, Vec g, void *ctx)
virtual int get_converged_reason() libmesh_override
PetscErrorCode __libmesh_tao_inequality_constraints_jacobian(Tao, Vec x, Mat J, Mat Jpre, void *ctx)
const NumericVector< Number > & get_vector(const std::string &vec_name) const
std::unique_ptr< SparseMatrix< Number > > C_eq_jac
virtual void equality_constraints(const NumericVector< Number > &X, NumericVector< Number > &C_eq, sys_type &S)=0
std::unique_ptr< NumericVector< Number > > lambda_eq
virtual void lower_and_upper_bounds(sys_type &S)=0
virtual void get_dual_variables() libmesh_override
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
TaoConvergedReason _reason
OptimizationSystem::ComputeLowerAndUpperBounds * lower_and_upper_bounds_object
const Parallel::Communicator & comm() const
SparseMatrix< Number > * matrix
virtual void solve() libmesh_override
SparseMatrix interface to PETSc Mat.
virtual void print_converged_reason() libmesh_override
virtual void swap(NumericVector< T > &v) libmesh_override
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 inequality_constraints(const NumericVector< Number > &X, NumericVector< Number > &C_ineq, sys_type &S)=0
OptimizationSystem::ComputeInequalityConstraintsJacobian * inequality_constraints_jacobian_object
std::unique_ptr< NumericVector< Number > > C_ineq
OptimizationSystem::ComputeEqualityConstraints * equality_constraints_object