Base class for local/global surrogate-based optimization/least squares. More...
Protected Member Functions | |
SurrBasedMinimizer (Model &model) | |
constructor | |
~SurrBasedMinimizer () | |
destructor | |
void | initialize_graphics (bool graph_2d, bool tabular_data, const String &tabular_file) |
initialize graphics customized for surrogate-based iteration | |
void | run () |
run portion of run_iterator; implemented by all derived classes and may include pre/post steps in lieu of separate pre/post | |
void | print_results (std::ostream &s) |
virtual void | minimize_surrogates ()=0 |
Used for computing the optimal solution using a surrogate-based approach. Redefines the Iterator::run() virtual function. | |
void | update_lagrange_multipliers (const RealVector &fn_vals, const RealMatrix &fn_grads) |
initialize and update Lagrange multipliers for basic Lagrangian | |
void | update_augmented_lagrange_multipliers (const RealVector &fn_vals) |
initialize and update the Lagrange multipliers for augmented Lagrangian | |
bool | update_filter (const RealVector &fn_vals) |
update a filter from a set of function values | |
Real | lagrangian_merit (const RealVector &fn_vals, const RealVector &primary_wts, const RealVector &nln_ineq_l_bnds, const RealVector &nln_ineq_u_bnds, const RealVector &nln_eq_tgts) |
compute a Lagrangian function from a set of function values | |
void | lagrangian_gradient (const RealVector &fn_vals, const RealMatrix &fn_grads, const RealVector &primary_wts, const RealVector &nln_ineq_l_bnds, const RealVector &nln_ineq_u_bnds, const RealVector &nln_eq_tgts, RealVector &lag_grad) |
compute the gradient of the Lagrangian function | |
Real | augmented_lagrangian_merit (const RealVector &fn_vals, const RealVector &primary_wts, const RealVector &nln_ineq_l_bnds, const RealVector &nln_ineq_u_bnds, const RealVector &nln_eq_tgts) |
compute an augmented Lagrangian function from a set of function values | |
void | augmented_lagrangian_gradient (const RealVector &fn_vals, const RealMatrix &fn_grads, const RealVector &primary_wts, const RealVector &nln_ineq_l_bnds, const RealVector &nln_ineq_u_bnds, const RealVector &nln_eq_tgts, RealVector &alag_grad) |
compute the gradient of the augmented Lagrangian function | |
Real | penalty_merit (const RealVector &fn_vals, const RealVector &primary_wts) |
compute a penalty function from a set of function values | |
void | penalty_gradient (const RealVector &fn_vals, const RealMatrix &fn_grads, const RealVector &primary_wts, RealVector &pen_grad) |
compute the gradient of the penalty function | |
Real | constraint_violation (const RealVector &fn_vals, const Real &constraint_tol) |
compute the constraint violation from a set of function values | |
Protected Attributes | |
Iterator | approxSubProbMinimizer |
the minimizer used on the surrogate model to solve the approximate subproblem on each surrogate-based iteration | |
int | sbIterNum |
surrogate-based minimization iteration number | |
RealVectorArray | sbFilter |
Set of response function vectors defining a filter (objective vs. constraint violation) for iterate selection/rejection. | |
RealVector | lagrangeMult |
Lagrange multipliers for basic Lagrangian calculations. | |
RealVector | augLagrangeMult |
Lagrange multipliers for augmented Lagrangian calculations. | |
Real | penaltyParameter |
the penalization factor for violated constraints used in quadratic penalty calculations; increased in update_penalty() | |
RealVector | origNonlinIneqLowerBnds |
original nonlinear inequality constraint lower bounds (no relaxation) | |
RealVector | origNonlinIneqUpperBnds |
original nonlinear inequality constraint upper bounds (no relaxation) | |
RealVector | origNonlinEqTargets |
original nonlinear equality constraint targets (no relaxation) | |
Real | eta |
constant used in etaSequence updates | |
Real | alphaEta |
power for etaSequence updates when updating penalty | |
Real | betaEta |
power for etaSequence updates when updating multipliers | |
Real | etaSequence |
decreasing sequence of allowable constraint violation used in augmented Lagrangian updates (refer to Conn, Gould, and Toint, section 14.4) |
Base class for local/global surrogate-based optimization/least squares.
These minimizers use a SurrogateModel to perform optimization based either on local trust region methods or global updating methods.
void run | ( | ) | [inline, protected, virtual] |
run portion of run_iterator; implemented by all derived classes and may include pre/post steps in lieu of separate pre/post
Virtual run function for the iterator class hierarchy. All derived classes need to redefine it.
Reimplemented from Iterator.
References SurrBasedMinimizer::minimize_surrogates().
void print_results | ( | std::ostream & | s | ) | [protected, virtual] |
Redefines default iterator results printing to include optimization results (objective functions and constraints).
Reimplemented from Iterator.
References Dakota::abort_handler(), Iterator::activeSet, String::begins(), Iterator::bestResponseArray, Iterator::bestVariablesArray, Dakota::data_pairs, Model::interface_id(), Iterator::iteratedModel, Dakota::lookup_by_val(), Iterator::methodName, Iterator::numFunctions, Minimizer::numUserPrimaryFns, Minimizer::optimizationFlag, ActiveSet::request_values(), Model::truth_model(), and Dakota::write_data_partial().
void update_lagrange_multipliers | ( | const RealVector & | fn_vals, |
const RealMatrix & | fn_grads | ||
) | [protected] |
initialize and update Lagrange multipliers for basic Lagrangian
For the Rockafellar augmented Lagrangian, simple Lagrange multiplier updates are available which do not require the active constraint gradients. For the basic Lagrangian, Lagrange multipliers are estimated through solution of a nonnegative linear least squares problem.
References Dakota::abort_handler(), Minimizer::bigRealBoundSize, Minimizer::constraintTol, Iterator::iteratedModel, SurrBasedMinimizer::lagrangeMult, Iterator::numContinuousVars, Minimizer::numNonlinearEqConstraints, Minimizer::numNonlinearIneqConstraints, Minimizer::numUserPrimaryFns, Minimizer::objective_gradient(), SurrBasedMinimizer::origNonlinIneqLowerBnds, SurrBasedMinimizer::origNonlinIneqUpperBnds, and Model::primary_response_fn_weights().
Referenced by SurrBasedLocalMinimizer::hard_convergence_check().
void update_augmented_lagrange_multipliers | ( | const RealVector & | fn_vals | ) | [protected] |
initialize and update the Lagrange multipliers for augmented Lagrangian
For the Rockafellar augmented Lagrangian, simple Lagrange multiplier updates are available which do not require the active constraint gradients. For the basic Lagrangian, Lagrange multipliers are estimated through solution of a nonnegative linear least squares problem.
References SurrBasedMinimizer::augLagrangeMult, SurrBasedMinimizer::betaEta, Minimizer::bigRealBoundSize, SurrBasedMinimizer::etaSequence, Minimizer::numNonlinearEqConstraints, Minimizer::numNonlinearIneqConstraints, Minimizer::numUserPrimaryFns, SurrBasedMinimizer::origNonlinEqTargets, SurrBasedMinimizer::origNonlinIneqLowerBnds, SurrBasedMinimizer::origNonlinIneqUpperBnds, and SurrBasedMinimizer::penaltyParameter.
Referenced by SurrBasedLocalMinimizer::hard_convergence_check(), EffGlobalMinimizer::minimize_surrogates_on_model(), and SurrBasedLocalMinimizer::tr_ratio_check().
bool update_filter | ( | const RealVector & | fn_vals | ) | [protected] |
update a filter from a set of function values
Update the sbFilter with fn_vals if new iterate is non-dominated.
References SurrBasedMinimizer::constraint_violation(), Iterator::iteratedModel, Minimizer::numNonlinearConstraints, Minimizer::objective(), Model::primary_response_fn_weights(), and SurrBasedMinimizer::sbFilter.
Referenced by SurrBasedLocalMinimizer::hard_convergence_check(), and SurrBasedLocalMinimizer::tr_ratio_check().
Real lagrangian_merit | ( | const RealVector & | fn_vals, |
const RealVector & | primary_wts, | ||
const RealVector & | nln_ineq_l_bnds, | ||
const RealVector & | nln_ineq_u_bnds, | ||
const RealVector & | nln_eq_tgts | ||
) | [protected] |
compute a Lagrangian function from a set of function values
The Lagrangian function computation sums the objective function and the Lagrange multipler terms for inequality/equality constraints. This implementation follows the convention in Vanderplaats with g<=0 and h=0. The bounds/targets passed in may reflect the original constraints or the relaxed constraints.
References Minimizer::bigRealBoundSize, Minimizer::constraintTol, SurrBasedMinimizer::lagrangeMult, Minimizer::numNonlinearEqConstraints, Minimizer::numNonlinearIneqConstraints, Minimizer::numUserPrimaryFns, and Minimizer::objective().
Referenced by SurrBasedLocalMinimizer::approx_subprob_objective_eval(), and SurrBasedLocalMinimizer::tr_ratio_check().
Real augmented_lagrangian_merit | ( | const RealVector & | fn_vals, |
const RealVector & | primary_wts, | ||
const RealVector & | nln_ineq_l_bnds, | ||
const RealVector & | nln_ineq_u_bnds, | ||
const RealVector & | nln_eq_tgts | ||
) | [protected] |
compute an augmented Lagrangian function from a set of function values
The Rockafellar augmented Lagrangian function sums the objective function, Lagrange multipler terms for inequality/equality constraints, and quadratic penalty terms for inequality/equality constraints. This implementation follows the convention in Vanderplaats with g<=0 and h=0. The bounds/targets passed in may reflect the original constraints or the relaxed constraints.
References SurrBasedMinimizer::augLagrangeMult, Minimizer::bigRealBoundSize, Minimizer::numNonlinearEqConstraints, Minimizer::numNonlinearIneqConstraints, Minimizer::numUserPrimaryFns, Minimizer::objective(), and SurrBasedMinimizer::penaltyParameter.
Referenced by SurrBasedLocalMinimizer::approx_subprob_objective_eval(), EffGlobalMinimizer::get_best_sample(), EffGlobalMinimizer::minimize_surrogates_on_model(), and SurrBasedLocalMinimizer::tr_ratio_check().
Real penalty_merit | ( | const RealVector & | fn_vals, |
const RealVector & | primary_wts | ||
) | [protected] |
compute a penalty function from a set of function values
The penalty function computation applies a quadratic penalty to any constraint violations and adds this to the objective function(s) p = f + r_p cv.
References SurrBasedMinimizer::constraint_violation(), Minimizer::constraintTol, Minimizer::objective(), and SurrBasedMinimizer::penaltyParameter.
Referenced by SurrBasedLocalMinimizer::tr_ratio_check().
Real constraint_violation | ( | const RealVector & | fn_vals, |
const Real & | constraint_tol | ||
) | [protected] |
compute the constraint violation from a set of function values
Compute the quadratic constraint violation defined as cv = g+^T g+ + h+^T h+. This implementation supports equality constraints and 2-sided inequalities. The constraint_tol allows for a small constraint infeasibility (used for penalty methods, but not Lagrangian methods).
References Minimizer::bigRealBoundSize, Minimizer::numNonlinearEqConstraints, Minimizer::numNonlinearIneqConstraints, Minimizer::numUserPrimaryFns, SurrBasedMinimizer::origNonlinEqTargets, SurrBasedMinimizer::origNonlinIneqLowerBnds, and SurrBasedMinimizer::origNonlinIneqUpperBnds.
Referenced by SurrBasedLocalMinimizer::hard_convergence_check(), EffGlobalMinimizer::minimize_surrogates_on_model(), SurrBasedMinimizer::penalty_merit(), SurrBasedLocalMinimizer::relax_constraints(), SurrBasedLocalMinimizer::tr_ratio_check(), SurrBasedMinimizer::update_filter(), and SurrBasedLocalMinimizer::update_penalty().