1 /* _______________________________________________________________________ 2 3 DAKOTA: Design Analysis Kit for Optimization and Terascale Applications 4 Copyright 2014-2020 National Technology & Engineering Solutions of Sandia, LLC (NTESS). 5 This software is distributed under the GNU Lesser General Public License. 6 For more information, see the README file in the top Dakota directory. 7 _______________________________________________________________________ */ 8 9 //- Class: SurrBasedMinimizer 10 //- Description: Base class for local and global surrogate-based optimization 11 //- and nonlinear least squares. 12 //- Owner: Mike Eldred 13 //- Checked by: 14 //- Version: $Id: SurrBasedMinimizer.hpp 4461 2007-08-28 17:40:08Z mseldre $ 15 16 #ifndef SURR_BASED_MINIMIZER_H 17 #define SURR_BASED_MINIMIZER_H 18 19 #include "DakotaMinimizer.hpp" 20 #include "DakotaModel.hpp" 21 22 namespace Dakota { 23 24 class SurrBasedLevelData; 25 26 27 /// Base class for local/global surrogate-based optimization/least squares. 28 29 /** These minimizers use a SurrogateModel to perform optimization based 30 either on local trust region methods or global updating methods. */ 31 32 class SurrBasedMinimizer: public Minimizer 33 { 34 protected: 35 36 // 37 //- Heading: Constructor and destructor 38 // 39 40 SurrBasedMinimizer(ProblemDescDB& problem_db, Model& model, std::shared_ptr<TraitsBase> traits); ///< constructor 41 ~SurrBasedMinimizer(); ///< destructor 42 43 // 44 //- Heading: Virtual function redefinitions 45 // 46 47 void derived_init_communicators(ParLevLIter pl_iter); 48 void derived_set_communicators(ParLevLIter pl_iter); 49 void derived_free_communicators(ParLevLIter pl_iter); 50 51 void print_results(std::ostream& s, short results_state = FINAL_RESULTS); 52 53 // 54 //- Heading: Utility member functions 55 // 56 57 /// initialize and update Lagrange multipliers for basic Lagrangian 58 void update_lagrange_multipliers(const RealVector& fn_vals, 59 const RealMatrix& fn_grads, SurrBasedLevelData& tr_data); 60 61 /// initialize and update the Lagrange multipliers for augmented Lagrangian 62 void update_augmented_lagrange_multipliers(const RealVector& fn_vals); 63 64 /// (re-)initialize filter from a set of function values 65 void initialize_filter(SurrBasedLevelData& tr_data, 66 const RealVector& fn_vals); 67 /// update filter using a new set of function values 68 bool update_filter(SurrBasedLevelData& tr_data, const RealVector& fn_vals); 69 70 // compute a filter merit function from a set of function values 71 //Real filter_merit(const RealVector& fns_center, const RealVector& fns_star); 72 73 /// compute a Lagrangian function from a set of function values 74 Real lagrangian_merit(const RealVector& fn_vals, 75 const BoolDeque& sense, const RealVector& primary_wts, 76 const RealVector& nln_ineq_l_bnds, 77 const RealVector& nln_ineq_u_bnds, 78 const RealVector& nln_eq_tgts); 79 /// compute the gradient of the Lagrangian function 80 void lagrangian_gradient(const RealVector& fn_vals, 81 const RealMatrix& fn_grads, 82 const BoolDeque& sense, 83 const RealVector& primary_wts, 84 const RealVector& nln_ineq_l_bnds, 85 const RealVector& nln_ineq_u_bnds, 86 const RealVector& nln_eq_tgts, 87 RealVector& lag_grad); 88 /// compute the Hessian of the Lagrangian function 89 void lagrangian_hessian(const RealVector& fn_vals, 90 const RealMatrix& fn_grads, 91 const RealSymMatrixArray& fn_hessians, 92 const BoolDeque& sense, 93 const RealVector& primary_wts, 94 const RealVector& nln_ineq_l_bnds, 95 const RealVector& nln_ineq_u_bnds, 96 const RealVector& nln_eq_tgts, 97 RealSymMatrix& lag_hess); 98 99 /// compute an augmented Lagrangian function from a set of function values 100 Real augmented_lagrangian_merit(const RealVector& fn_vals, 101 const BoolDeque& sense, 102 const RealVector& primary_wts, 103 const RealVector& nln_ineq_l_bnds, 104 const RealVector& nln_ineq_u_bnds, 105 const RealVector& nln_eq_tgts); 106 /// compute the gradient of the augmented Lagrangian function 107 void augmented_lagrangian_gradient(const RealVector& fn_vals, 108 const RealMatrix& fn_grads, 109 const BoolDeque& sense, 110 const RealVector& primary_wts, 111 const RealVector& nln_ineq_l_bnds, 112 const RealVector& nln_ineq_u_bnds, 113 const RealVector& nln_eq_tgts, 114 RealVector& alag_grad); 115 /// compute the Hessian of the augmented Lagrangian function 116 void augmented_lagrangian_hessian(const RealVector& fn_vals, 117 const RealMatrix& fn_grads, 118 const RealSymMatrixArray& fn_hessians, 119 const BoolDeque& sense, 120 const RealVector& primary_wts, 121 const RealVector& nln_ineq_l_bnds, 122 const RealVector& nln_ineq_u_bnds, 123 const RealVector& nln_eq_tgts, 124 RealSymMatrix& alag_hess); 125 126 /// compute a penalty function from a set of function values 127 Real penalty_merit(const RealVector& fn_vals, const BoolDeque& sense, 128 const RealVector& primary_wts); 129 130 /// compute the gradient of the penalty function 131 void penalty_gradient(const RealVector& fn_vals, const RealMatrix& fn_grads, 132 const BoolDeque& sense, const RealVector& primary_wts, 133 RealVector& pen_grad); 134 135 /// compute the constraint violation from a set of function values 136 Real constraint_violation(const RealVector& fn_vals, 137 const Real& constraint_tol); 138 139 // 140 //- Heading: Data members 141 // 142 143 /// the minimizer used on the surrogate model to solve the 144 /// approximate subproblem on each surrogate-based iteration 145 Iterator approxSubProbMinimizer; 146 147 /// global iteration counter corresponding to number of 148 /// surrogate-based minimizations 149 size_t globalIterCount; 150 151 /// Lagrange multipliers for basic Lagrangian calculations 152 RealVector lagrangeMult; 153 /// Lagrange multipliers for augmented Lagrangian calculations 154 RealVector augLagrangeMult; 155 /// the penalization factor for violated constraints used in quadratic 156 /// penalty calculations; increased in update_penalty() 157 Real penaltyParameter; 158 159 /// original nonlinear inequality constraint lower bounds (no relaxation) 160 RealVector origNonlinIneqLowerBnds; 161 /// original nonlinear inequality constraint upper bounds (no relaxation) 162 RealVector origNonlinIneqUpperBnds; 163 /// original nonlinear equality constraint targets (no relaxation) 164 RealVector origNonlinEqTargets; 165 166 /// constant used in etaSequence updates 167 Real eta; 168 /// power for etaSequence updates when updating penalty 169 Real alphaEta; 170 /// power for etaSequence updates when updating multipliers 171 Real betaEta; 172 /// decreasing sequence of allowable constraint violation used in augmented 173 /// Lagrangian updates (refer to Conn, Gould, and Toint, section 14.4) 174 Real etaSequence; 175 176 /// index for the active ParallelLevel within ParallelConfiguration::miPLIters 177 size_t miPLIndex; 178 }; 179 180 } // namespace Dakota 181 182 #endif 183