Changeset 416


Ignore:
Timestamp:
Jul 29, 2005 3:11:41 PM (14 years ago)
Author:
andreasw
Message:
  • revised handling of "acceptable level of accuracy" (now in ConvCheck?)
  • fixed uncaught evaluation error exceptions
Location:
branches/dev/Algorithm
Files:
18 edited

Legend:

Unmodified
Added
Removed
  • branches/dev/Algorithm/IpAdaptiveMuUpdate.cpp

    r414 r416  
    5757      "This option allows to specify a lower bound on the adaptively chosen "
    5858      "barrier parameter.  By default, it is set to "
    59       "0.1*min(\"tol\",\"compl_inf_tol\"), which should be a very reasonable "
    60       "value.");
     59      "min(\"tol\",\"compl_inf_tol\")/(\"barrier_tol_factor\"+1), which "
     60      "should be a very reasonable value.");
    6161    roptions->AddLowerBoundedNumberOption(
    6262      "adaptive_mu_safeguard_factor",
     
    143143      "in QualityFunctionMuOracle]");
    144144
    145     // ToDo move the last two options here to QualityFunctionMuOracle?
    146     roptions->AddStringOption4(
    147       "adaptive_mu_kkt_centrality",
    148       "Determines whether a penalty term for centrality is added to KKT error.",
    149       "none",
    150       "none", "no penalty term is added",
    151       "log", "complementarity * the log of the centrality measure",
    152       "reciprocal", "complementarity * the reciprocal of the centrality measure",
    153       "cubed-reciprocal", "complementarity * the reciprocal of the centrality measure cubed",
    154       "This determines whether a term penalizing deviation from centrality "
    155       "with respect to complementarity is added to the KKT "
    156       "error measure used in the globalization atrategies.  The "
    157       "complementarity measure here is the xi in the Loqo update rule. "
    158       "(This originated from including such a term in the quality "
    159       "function for the QualityFunctionMuOracle and needs to be used "
    160       "here for consistency.)");
    161 
    162     roptions->AddStringOption2(
    163       "adaptive_mu_kkt_balancing_term",
    164       "Determines whether a balancing term for centrality is added to KKT error.",
    165       "none",
    166       "none", "no balancing term is added",
    167       "cubic", "Max(0,Max(dual_ing,primal_inf)-compl)^3",
    168       "This determines whether a term penalizing stuations there the "
    169       "complementality is much smaller than dual and primal "
    170       "infeasibilities is added to the KKT error measure used in the "
    171       "globalization strategies. (This originated from including such "
    172       "a term in the quality function for the QualityFunctionMuOracle "
    173       "and needs to be used here for consistency.)");
    174 
    175145  }
    176146
     
    179149  {
    180150    options.GetNumericValue("mu_max", mu_max_, prefix);
    181     if (!options.GetNumericValue("mu_min", mu_min_, prefix)) {
    182       mu_min_ = 0.1*Min(IpData().tol(), IpData().compl_inf_tol());
    183     }
    184151    options.GetNumericValue("tau_min", tau_min_, prefix);
    185152    options.GetNumericValue("adaptive_mu_safeguard_factor", adaptive_mu_safeguard_factor_, prefix);
     
    218185    options.GetEnumValue("quality_function_balancing_term", enum_int, prefix);
    219186    adaptive_mu_kkt_balancing_term_ = QualityFunctionMuOracle::BalancingTermEnum(enum_int);
     187    options.GetNumericValue("compl_inf_tol", compl_inf_tol_, prefix);
     188    if (!options.GetNumericValue("mu_min", mu_min_, prefix)) {
     189      // Defer computation of the default until the scaling of the NLP
     190      // is known
     191      mu_min_ = -1.;
     192    }
    220193
    221194    init_dual_inf_ = -1.;
     
    240213  void AdaptiveMuUpdate::UpdateBarrierParameter()
    241214  {
     215    DBG_START_METH("AdaptiveMuUpdate::UpdateBarrierParameter",
     216                   dbg_verbosity);
     217
     218    // if min_mu_ has not been given, we now set the default (can't do
     219    // that earlier, because during call of InitializeImpl, the
     220    // scaling in the NLP is not yet determined
     221    if (mu_min_ < 0.) {
     222      mu_min_ = Min(IpData().tol(),
     223                    IpNLP().NLP_scaling()->apply_obj_scaling(compl_inf_tol_))/
     224                (barrier_tol_factor_+1.);
     225    }
     226
    242227    // if there are not bounds, we always return the minimum MU value
    243228    if (!check_if_no_bounds_) {
     
    282267          // ToDo combine this code with MonotoneMuUpdate
    283268          Number tol = IpData().tol();
    284           Number compl_inf_tol = IpData().compl_inf_tol();
     269          Number compl_inf_tol =
     270            IpNLP().NLP_scaling()->apply_obj_scaling(compl_inf_tol_);
    285271
    286272          Number new_mu = Min( mu_linear_decrease_factor_*mu,
    287273                               pow(mu, mu_superlinear_decrease_power_) );
    288           new_mu = Max(new_mu, Min(compl_inf_tol, tol)/10.);
     274          DBG_PRINT((1,"new_mu = %e, compl_inf_tol = %e tol = %e\n", new_mu, compl_inf_tol, tol));
     275          new_mu = Max(new_mu,
     276                       Min(compl_inf_tol, tol)/(barrier_tol_factor_+1.));
    289277          if (tiny_step_flag && new_mu == mu) {
    290278            THROW_EXCEPTION(TINY_STEP_DETECTED,
     
    688676  AdaptiveMuUpdate::lower_mu_safeguard()
    689677  {
     678    DBG_START_METH("AdaptiveMuUpdate::lower_mu_safeguard",
     679                   dbg_verbosity);
    690680    if (adaptive_mu_safeguard_factor_ == 0.)
    691681      return 0.;
     
    713703      Max(adaptive_mu_safeguard_factor_ * (dual_inf/init_dual_inf_),
    714704          adaptive_mu_safeguard_factor_ * (primal_inf/init_primal_inf_));
     705    DBG_PRINT((1,"dual_inf=%e init_dual_inf_=%e primal_inf=%e init_primal_inf_=%e\n", dual_inf, init_dual_inf_, primal_inf, init_primal_inf_));
    715706
    716707    if (adaptive_mu_globalization_==KKT_ERROR) {
  • branches/dev/Algorithm/IpAdaptiveMuUpdate.hpp

    r409 r416  
    100100    /** Factor for filter margin */
    101101    Number filter_margin_fact_;
     102    /** Unscaled tolerance for complementarity */
     103    Number compl_inf_tol_;
    102104    //@}
    103105
  • branches/dev/Algorithm/IpAlgBuilder.cpp

    r378 r416  
    185185      new PDFullSpaceSolver(*resto_AugSolver);
    186186
     187    // Convergence check in the restoration phase
     188    SmartPtr<RestoFilterConvergenceCheck> resto_convCheck =
     189      new RestoFilterConvergenceCheck();
     190
    187191    // Line search method for the restoration phase
    188192    SmartPtr<RestoRestorationPhase> resto_resto =
    189193      new RestoRestorationPhase();
    190194    SmartPtr<FilterLineSearch> resto_LineSearch =
    191       new FilterLineSearch(GetRawPtr(resto_resto), GetRawPtr(resto_PDSolver));
     195      new FilterLineSearch(GetRawPtr(resto_resto), GetRawPtr(resto_PDSolver),
     196                           GetRawPtr(resto_convCheck));
    192197
    193198    // Create the mu update that will be used by the restoration phase
     
    235240                             resto_MuOracle, resto_FixMuOracle);
    236241    }
    237 
    238     // Convergence check in the restoration phase
    239     SmartPtr<RestoFilterConvergenceCheck> resto_convCheck =
    240       new RestoFilterConvergenceCheck();
    241242
    242243    // Initialization of the iterates for the restoration phase
     
    267268    // Create the line search to be used by the main algorithm
    268269    SmartPtr<FilterLineSearch> lineSearch =
    269       new FilterLineSearch(GetRawPtr(resto_phase), GetRawPtr(PDSolver));
     270      new FilterLineSearch(GetRawPtr(resto_phase), GetRawPtr(PDSolver),
     271                           convCheck);
    270272
    271273    // The following cross reference is not good: We have to store a
  • branches/dev/Algorithm/IpConvCheck.hpp

    r2 r416  
    3737      CONTINUE,
    3838      CONVERGED,
     39      CONVERGED_TO_ACCEPTABLE_POINT,
    3940      MAXITER_EXCEEDED,
    4041      FAILED
     
    4748    /** Pure virtual method for performing the convergence test */
    4849    virtual ConvergenceStatus CheckConvergence()=0;
     50
     51    /** Method for testing if the current iterate is considered to
     52     *  satisfy the "accptable level" of accuracy.  The idea is that
     53     *  if the desired convergence tolerance cannot be achieved, the
     54     *  algorithm might stop after a number of acceptable points have
     55     *  been encountered. */
     56    virtual bool CurrentIsAcceptable()=0;
    4957
    5058  private:
  • branches/dev/Algorithm/IpFilterLineSearch.cpp

    r412 r416  
    2828
    2929  FilterLineSearch::FilterLineSearch(const SmartPtr<RestorationPhase>& resto_phase,
    30                                      const SmartPtr<PDSystemSolver>& pd_solver)
     30                                     const SmartPtr<PDSystemSolver>& pd_solver,
     31                                     const SmartPtr<ConvergenceCheck>& conv_check)
    3132      :
    3233      LineSearch(),
    3334      resto_phase_(resto_phase),
    3435      pd_solver_(pd_solver),
     36      conv_check_(conv_check),
    3537      theta_min_(-1.0),
    3638      theta_max_(-1.0),
     
    254256      "Determines the number of trial iterations before the watchdog "
    255257      "procedure is aborted and the algorithm returns to the stored point.");
    256     roptions->AddLowerBoundedNumberOption(
    257       "acceptable_tol",
    258       "Threshold for NLP error to consider iterate as acceptable.",
    259       0.0, false, 1e-6,
    260       "Determines tolerance for which an iterate is considered as accepable "
    261       "solution.  If the algorithm would trigger the restoration phase at "
    262       "such a point, it instead terminates, returning this acceptable "
    263       "point.  Further, if the algorithm encounters \"acceptable_iter_max\" "
    264       "successive points satisfying this NLP tolerance, the algorithm "
    265       "terminates.");
    266     roptions->AddLowerBoundedIntegerOption(
    267       "acceptable_iter_max",
    268       "Number of acceptable iterates to trigger termination.",
    269       0, 15,
    270       "if the algorithm encounters so many successive acceptable iterates "
    271       "(see \"acceptable_tol\"), it terminates, assuming that the problem "
    272       "has been solved to best possible accuracy given round-off.");
    273258  }
    274259
     
    321306    options.GetIntegerValue("watchdog_trial_iter_max", watchdog_trial_iter_max_, prefix);
    322307    options.GetIntegerValue("watchdog_shortened_iter_trigger", watchdog_shortened_iter_trigger_, prefix);
    323     options.GetNumericValue("acceptable_tol", acceptable_tol_, prefix);
    324     options.GetIntegerValue("acceptable_iter_max", acceptable_iter_max_, prefix);
    325308
    326309    // ToDo decide if also the PDSystemSolver should be initialized here...
     
    333316
    334317    count_successive_shortened_steps_ = 0;
    335     count_acceptable_iter_ = 0;
    336318
    337319    acceptable_iterate_ = NULL;
     
    358340    // Store current iterate if the optimality error is on acceptable
    359341    // level to restored if things fail later
    360     if (IpCq().curr_nlp_error()<=acceptable_tol_) {
     342    if (CurrentIsAcceptable()) {
    361343      Jnlst().Printf(J_DETAILED, J_LINE_SEARCH,
    362344                     "Storing current iterate as backup acceptable point.\n");
    363       IpData().Append_info_string("A");
    364345      StoreAcceptablePoint();
    365346    }
     
    551532            AugmentFilter();
    552533          }
    553           if (IpCq().curr_nlp_error()<=acceptable_tol_) {
     534          if (CurrentIsAcceptable()) {
    554535            THROW_EXCEPTION(ACCEPTABLE_POINT_REACHED,
    555536                            "Restoration phase called at acceptable point.");
     
    559540            THROW_EXCEPTION(IpoptException, "No Restoration Phase given to this Filter Line Search Object!");
    560541          }
     542          // ToDo make the 1e-2 below a parameter?
    561543          if (IpCq().curr_constraint_violation()<=
    562               1e-2*Min(IpData().tol(),IpData().primal_inf_tol())) {
     544              1e-2*IpData().tol()) {
    563545            bool found_acceptable = RestoreAcceptablePoint();
    564546            if (found_acceptable) {
    565547              THROW_EXCEPTION(ACCEPTABLE_POINT_REACHED,
    566                               "Restoration phase called at almost feasible point, but acceptable point could be restore.\n");
     548                              "Restoration phase called at almost feasible point, but acceptable point could be restored.\n");
    567549            }
    568550            else {
     551              // ToDo does that happen too often?
    569552              THROW_EXCEPTION(RESTORATION_FAILED,
    570553                              "Restoration phase called, but point is almost feasible.");
     
    593576          }
    594577          count_successive_shortened_steps_ = 0;
    595           count_acceptable_iter_ = 0;
    596578          if (expect_infeasible_problem_) {
    597579            expect_infeasible_problem_ = false;
     
    611593
    612594      PerformDualStep(alpha_primal, alpha_dual_max, actual_delta);
    613 
    614       if (acceptable_iter_max_>0) {
    615         if (IpCq().curr_nlp_error()<=acceptable_tol_) {
    616           count_acceptable_iter_++;
    617           if (count_acceptable_iter_>=acceptable_iter_max_) {
    618             IpData().AcceptTrialPoint();
    619             THROW_EXCEPTION(ACCEPTABLE_POINT_REACHED,
    620                             "Algorithm seems stuck at acceptable level.");
    621           }
    622         }
    623         else {
    624           count_acceptable_iter_=0;
    625         }
    626       }
    627595
    628596      if (n_steps==0) {
     
    12601228        Jnlst().Printf(J_WARNING, J_MAIN, "Warning: SOC step rejected due to evaluation error\n");
    12611229        accept = false;
     1230        // There is no point in continuing SOC procedure
     1231        break;
    12621232      }
    12631233
     
    16241594  }
    16251595
     1596  bool FilterLineSearch::CurrentIsAcceptable()
     1597  {
     1598    return (IsValid(conv_check_) &&
     1599            conv_check_->CurrentIsAcceptable());
     1600  }
     1601
    16261602  void FilterLineSearch::StoreAcceptablePoint()
    16271603  {
  • branches/dev/Algorithm/IpFilterLineSearch.hpp

    r395 r416  
    1616#include "IpPDSystemSolver.hpp"
    1717#include "IpIpoptType.hpp"
     18#include "IpConvCheck.hpp"
    1819
    1920namespace Ipopt
     
    3233    /** Constructor.  The PDSystemSolver object only needs to be
    3334     *  provided (i.e. not NULL) if second order correction is to be
    34      *  used. */
     35     *  used.  The ConvergenceCheck object is used to determine
     36     *  whether the current iterate is acceptable (for example, the
     37     *  restoration phase is not started if the acceptability level
     38     *  has been reached).  If conv_check is NULL, we assume that the
     39     *  current iterate is not acceptable. */
    3540    FilterLineSearch(const SmartPtr<RestorationPhase>& resto_phase,
    36                      const SmartPtr<PDSystemSolver>& pd_solver
     41                     const SmartPtr<PDSystemSolver>& pd_solver,
     42                     const SmartPtr<ConvergenceCheck>& conv_check
    3743                    );
    3844
     
    242248     *  found. Returns true if such as point is available. */
    243249    bool RestoreAcceptablePoint();
     250
     251    /** Method for determining if the current iterate is acceptable.
     252     *  This is a wrapper for same method from ConvergenceCheck, but
     253     *  returns false, if no ConvergenceCheck object is provided. */
     254    bool CurrentIsAcceptable();
    244255
    245256    /** @name Parameters for the filter algorithm.  Names as in the paper */
     
    335346    Index watchdog_shortened_iter_trigger_;
    336347
    337     /** Acceptable tolerance for the problem to terminate earlier if
    338      *  algorithm seems stuck or cycling */
    339     Number acceptable_tol_;
    340     /** Maximum number of iterations with acceptable level of accuracy
    341      *  and full steps, after which the algorithm terminates.  If 0,
    342      *  this heuristic is disabled. */
    343     Index acceptable_iter_max_;
    344 
    345348    /** Indicates whether the algorithm should start directly with the
    346349     *  restoratin phase */
     
    408411    Index count_successive_shortened_steps_;
    409412
    410     /** Counter for the number of successive iterations in which the
    411      *  nlp error was below the acceptable tolerance and a full step
    412      *  was accepted. */
    413     Index count_acceptable_iter_;
    414 
    415413    /** Flag indicating if a tiny step was detected in previous
    416414     *  iteration */
    417415    bool tiny_step_last_iteration_;
    418416
     417    /** @name Strategy objective that are used */
     418    //@{
    419419    SmartPtr<RestorationPhase> resto_phase_;
    420420    SmartPtr<PDSystemSolver> pd_solver_;
     421    SmartPtr<ConvergenceCheck> conv_check_;
     422    //@}
    421423  };
    422424
  • branches/dev/Algorithm/IpIpoptAlg.cpp

    r410 r416  
    171171      if (conv_status == ConvergenceCheck::CONVERGED) {
    172172        return SUCCESS;
     173      }
     174      if (conv_status == ConvergenceCheck::CONVERGED_TO_ACCEPTABLE_POINT) {
     175        return STOP_AT_ACCEPTABLE_POINT;
    173176      }
    174177      else if (conv_status == ConvergenceCheck::MAXITER_EXCEEDED) {
  • branches/dev/Algorithm/IpIpoptData.cpp

    r414 r416  
    5353      "implementation paper).  [Some other algorithmic features also use "
    5454      "this quantity.]");
    55     reg_options->AddLowerBoundedNumberOption(
    56       "dual_inf_tol",
    57       "Acceptance threshold for the unscaled dual infeasibility.",
    58       0.0, true, 1e-2,
    59       "Absolute tolerance on the dual infesaibility.  Successful termination "
    60       "requires that the (unscaled) dual infeasibility is less than this "
    61       "threshold.");
    62     reg_options->AddLowerBoundedNumberOption(
    63       "primal_inf_tol",
    64       "Acceptance threshold for the unscaled primal infeasibility.",
    65       0.0, true, 1e-2,
    66       "Absolute tolerance on the dual infesaibility.  Successful termination "
    67       "requires that the (unscaled) primal infeasibility is less than this "
    68       "threshold.");
    69     reg_options->AddLowerBoundedNumberOption(
    70       "compl_inf_tol",
    71       "Acceptance threshold for the complementarity conditions.",
    72       0.0, true, 1e-2,
    73       "Absolute tolerance on the complementarity.  Successful termination "
    74       "requires that the (unscaled) complementarity is less than this "
    75       "threshold.");
    7655  }
    7756
     
    8160  {
    8261    options.GetNumericValue("tol", tol_, prefix);
    83     options.GetNumericValue("dual_inf_tol", dual_inf_tol_, prefix);
    84     options.GetNumericValue("primal_inf_tol", primal_inf_tol_, prefix);
    85     options.GetNumericValue("compl_inf_tol", compl_inf_tol_, prefix);
    8662
    8763    iter_count_=0;
  • branches/dev/Algorithm/IpIpoptData.hpp

    r414 r416  
    293293      return tol_;
    294294    }
    295     Number dual_inf_tol() const
    296     {
    297       DBG_ASSERT(initialize_called_);
    298       return dual_inf_tol_;
    299     }
    300     Number primal_inf_tol() const
    301     {
    302       DBG_ASSERT(initialize_called_);
    303       return primal_inf_tol_;
    304     }
    305     Number compl_inf_tol() const
    306     {
    307       DBG_ASSERT(initialize_called_);
    308       return compl_inf_tol_;
    309     }
    310295    //@}
    311296
     
    457442    /** Overall convergence tolerance */
    458443    Number tol_;
    459     /** Tolerance on unscaled dual infeasibility */
    460     Number dual_inf_tol_;
    461     /** Tolerance on unscaled primal infeasibility */
    462     Number primal_inf_tol_;
    463     /** Tolerance on unscaled complementarity */
    464     Number compl_inf_tol_;
    465444    //@}
    466445
  • branches/dev/Algorithm/IpMonotoneMuUpdate.cpp

    r378 r416  
    9696    options.GetNumericValue("mu_superlinear_decrease_power", mu_superlinear_decrease_power_, prefix);
    9797    options.GetNumericValue("tau_min", tau_min_, prefix);
     98    options.GetNumericValue("compl_inf_tol", compl_inf_tol_, prefix);
    9899
    99100    IpData().Set_mu(mu_init_);
     
    182183    Number mu = IpData().curr_mu();
    183184    Number tol = IpData().tol();
    184     Number compl_inf_tol = IpData().compl_inf_tol();
     185
     186    // Here we need the complementarity tolerance that is posed to the
     187    // scaled problem
     188    Number compl_inf_tol =
     189      IpNLP().NLP_scaling()->apply_obj_scaling(compl_inf_tol_);
    185190
    186191    new_mu = Min( mu_linear_decrease_factor_*mu,
    187192                  pow(mu, mu_superlinear_decrease_power_) );
    188     new_mu = Max(new_mu, Min(tol, compl_inf_tol)/10.);
     193    new_mu = Max(new_mu, Min(tol, compl_inf_tol)/(barrier_tol_factor_+1.));
    189194
    190195    // update the fraction to the boundary parameter
  • branches/dev/Algorithm/IpMonotoneMuUpdate.hpp

    r375 r416  
    8181    /** Tau_min for fraction to boundary rule */
    8282    Number tau_min_;
     83    Number compl_inf_tol_;
    8384    //@}
    8485
  • branches/dev/Algorithm/IpOptErrorConvCheck.cpp

    r381 r416  
    2929      "The algorithm terminates with an error message if the number of "
    3030      "iterations exceeded this number. [Also used in RestoFilterConvCheck]");
     31    roptions->AddLowerBoundedNumberOption(
     32      "dual_inf_tol",
     33      "Acceptance threshold for the dual infeasibility.",
     34      0.0, true, 1e-4,
     35      "Absolute tolerance on the dual infesaibility. Successful termination "
     36      "requires that the (unscaled) dual infeasibility is less than this "
     37      "threshold.");
     38    roptions->AddLowerBoundedNumberOption(
     39      "constr_viol_tol",
     40      "Acceptance threshold for the constraint violation.",
     41      0.0, true, 1e-4,
     42      "Absolute tolerance on the constraint violation. Successful termination "
     43      "requires that the (unscaled) constraint violation is less than this "
     44      "threshold.");
     45    roptions->AddLowerBoundedNumberOption(
     46      "compl_inf_tol",
     47      "Acceptance threshold for the complementarity conditions.",
     48      0.0, true, 1e-4,
     49      "Absolute tolerance on the complementarity. Successful termination "
     50      "requires that the (unscaled) complementarity is less than this "
     51      "threshold.");
     52    roptions->AddLowerBoundedIntegerOption(
     53      "acceptable_iter",
     54      "Number of acceptable iterates to trigger termination.",
     55      0, 15,
     56      "If the algorithm encounters so many successive acceptable iterates "
     57      "(see \"acceptable_tol\"), it terminates, assuming that the problem "
     58      "has been solved to best possible accuracy given round-off.  If it is "
     59      "set to zero, this heuristic is disabled.");
     60    roptions->AddLowerBoundedNumberOption(
     61      "acceptable_tol",
     62      "Acceptable convergence tolerance (relative).",
     63      0.0, true,  1e-6,
     64      "Determines which (scaled) overall optimality error is considered to be "
     65      "acceptable. If the algorithm encounters \"acceptable_iter\" iterations "
     66      "in a row that are considered \"acceptable\", it will terminate before "
     67      "the desired convergence tolerance (\"tol\", \"dual_inf_tol\", etc) is "
     68      "met.");
     69    roptions->AddLowerBoundedNumberOption(
     70      "acceptable_dual_inf_tol",
     71      "Acceptance threshold for the dual infeasibility.",
     72      0.0, true, 1e-2,
     73      "Absolute tolerance on the dual infesaibility. Acceptable termination "
     74      "requires that the (unscaled) dual infeasibility is less than this "
     75      "threshold.");
     76    roptions->AddLowerBoundedNumberOption(
     77      "acceptable_constr_viol_tol",
     78      "Acceptance threshold for the constraint violation.",
     79      0.0, true, 1e-2,
     80      "Absolute tolerance on the constraint violation. Acceptable termination "
     81      "requires that the (unscaled) constraint violation is less than this "
     82      "threshold.");
     83    roptions->AddLowerBoundedNumberOption(
     84      "acceptable_compl_inf_tol",
     85      "Acceptance threshold for the complementarity conditions.",
     86      0.0, true, 1e-2,
     87      "Absolute tolerance on the complementarity. Acceptable termination "
     88      "requires that the (unscaled) complementarity is less than this "
     89      "threshold.");
    3190  }
    3291
     
    3695  {
    3796    options.GetIntegerValue("max_iter", max_iterations_, prefix);
     97    options.GetNumericValue("dual_inf_tol", dual_inf_tol_, prefix);
     98    options.GetNumericValue("constr_viol_tol", constr_viol_tol_, prefix);
     99    options.GetNumericValue("compl_inf_tol", compl_inf_tol_, prefix);
     100    options.GetIntegerValue("acceptable_iter", acceptable_iter_, prefix);
     101    options.GetNumericValue("acceptable_tol", acceptable_tol_, prefix);
     102    options.GetNumericValue("acceptable_dual_inf_tol", acceptable_dual_inf_tol_, prefix);
     103    options.GetNumericValue("acceptable_constr_viol_tol", acceptable_constr_viol_tol_, prefix);
     104    options.GetNumericValue("acceptable_compl_inf_tol", acceptable_compl_inf_tol_, prefix);
     105    acceptable_counter_ = 0;
    38106
    39107    return true;
    40108  }
    41109
    42   ConvergenceCheck::ConvergenceStatus OptimalityErrorConvergenceCheck::CheckConvergence()
     110  ConvergenceCheck::ConvergenceStatus
     111  OptimalityErrorConvergenceCheck::CheckConvergence()
    43112  {
    44113    DBG_START_METH("OptimalityErrorConvergenceCheck::CheckConvergence", dbg_verbosity);
    45     // maybe we should throw exceptions here instead?
    46114
    47115    if (IpData().iter_count() >= max_iterations_) {
     
    50118
    51119    Number overall_error = IpCq().curr_nlp_error();
    52     Number dual_inf = IpCq().curr_dual_infeasibility(NORM_MAX);
    53     Number primal_inf = IpCq().curr_primal_infeasibility(NORM_MAX);
    54     Number compl_inf = IpCq().curr_complementarity(0., NORM_MAX);
    55     DBG_PRINT((1,"overall_error = %8.2e\n",overall_error));
     120    Number dual_inf = IpCq().unscaled_curr_dual_infeasibility(NORM_MAX);
     121    Number constr_viol = IpCq().unscaled_curr_nlp_constraint_violation(NORM_MAX);
     122    Number compl_inf = IpCq().unscaled_curr_complementarity(0., NORM_MAX);
     123
    56124    if (overall_error <= IpData().tol() &&
    57         dual_inf <= IpData().dual_inf_tol() &&
    58         primal_inf <= IpData().primal_inf_tol() &&
    59         compl_inf <= IpData().compl_inf_tol()) {
     125        dual_inf <= dual_inf_tol_ &&
     126        constr_viol <= constr_viol_tol_ &&
     127        compl_inf <= compl_inf_tol_) {
    60128      return ConvergenceCheck::CONVERGED;
     129    }
     130
     131    if (acceptable_iter_>0 && CurrentIsAcceptable()) {
     132      IpData().Append_info_string("A");
     133      acceptable_counter_++;
     134      if (acceptable_counter_ >= acceptable_iter_) {
     135        return ConvergenceCheck::CONVERGED_TO_ACCEPTABLE_POINT;
     136      }
     137    }
     138    else {
     139      acceptable_counter_ = 0;
    61140    }
    62141
    63142    return ConvergenceCheck::CONTINUE;
    64143  }
     144
     145  bool OptimalityErrorConvergenceCheck::CurrentIsAcceptable()
     146  {
     147    DBG_START_METH("OptimalityErrorConvergenceCheck::CurrentIsAcceptable",
     148                   dbg_verbosity);
     149
     150    Number overall_error = IpCq().curr_nlp_error();
     151    Number dual_inf = IpCq().unscaled_curr_dual_infeasibility(NORM_MAX);
     152    Number constr_viol = IpCq().unscaled_curr_nlp_constraint_violation(NORM_MAX);
     153    Number compl_inf = IpCq().unscaled_curr_complementarity(0., NORM_MAX);
     154
     155    DBG_PRINT((1, "overall_error = %e\n", overall_error));
     156    DBG_PRINT((1, "dual_inf = %e\n", dual_inf));
     157    DBG_PRINT((1, "constr_viol = %e\n", constr_viol));
     158    DBG_PRINT((1, "compl_inf = %e\n", compl_inf));
     159
     160    DBG_PRINT((1, "acceptable_tol_ = %e\n", acceptable_tol_));
     161    DBG_PRINT((1, "acceptable_dual_inf_tol_ = %e\n", acceptable_dual_inf_tol_));
     162    DBG_PRINT((1, "acceptable_constr_viol_tol_ = %e\n", acceptable_constr_viol_tol_));
     163    DBG_PRINT((1, "acceptable_compl_inf_tol_ = %e\n", acceptable_compl_inf_tol_));
     164
     165    return (overall_error <= acceptable_tol_ &&
     166            dual_inf <= acceptable_dual_inf_tol_ &&
     167            constr_viol <= acceptable_constr_viol_tol_ &&
     168            compl_inf <= acceptable_compl_inf_tol_);
     169  }
     170
     171
    65172} // namespace Ipopt
  • branches/dev/Algorithm/IpOptErrorConvCheck.hpp

    r310 r416  
    4040    virtual ConvergenceStatus CheckConvergence();
    4141
     42    /** Auxilliary function for testing whether current iterate
     43     *  satisfies the acceptable level of optimality */
     44    virtual bool CurrentIsAcceptable();
     45
    4246    /** Methods for IpoptType */
    4347    //@{
    4448    static void RegisterOptions(SmartPtr<RegisteredOptions> roptions);
     49    //@}
     50
     51  protected:
     52    /** @name Algorithmic parameters */
     53    //@{
     54    /** Maximal number of iterations */
     55    Index max_iterations_;
     56    /** Tolerance on unscaled dual infeasibility */
     57    Number dual_inf_tol_;
     58    /** Tolerance on unscaled constraint violation */
     59    Number constr_viol_tol_;
     60    /** Tolerance on unscaled complementarity */
     61    Number compl_inf_tol_;
     62    /** Number of iterations with acceptable level of accuracy, after
     63     *  which the algorithm terminates.  If 0, this heuristic is
     64     *  disabled. */
     65    Index acceptable_iter_;
     66    /** Acceptable tolerance for the problem to terminate earlier if
     67     *  algorithm seems stuck or cycling */
     68    Number acceptable_tol_;
     69    /** Acceptable tolerance on unscaled dual infeasibility */
     70    Number acceptable_dual_inf_tol_;
     71    /** Acceptable tolerance on unscaled constraint violation */
     72    Number acceptable_constr_viol_tol_;
     73    /** Acceptable tolerance on unscaled complementarity */
     74    Number acceptable_compl_inf_tol_;
    4575    //@}
    4676
     
    5989    //@}
    6090
    61     /** @name Algorithmic parameters */
    62     //@{
    63     /** Maximal number of iterations */
    64     Index max_iterations_;
    65     //@}
     91    /** Counter for successive iterations in which acceptability
     92     *  criteria are met. */
     93    Index acceptable_counter_;
    6694  };
    6795
  • branches/dev/Algorithm/IpOrigIpoptNLP.cpp

    r414 r416  
    248248    DBG_START_METH("OrigIpoptNLP::f", dbg_verbosity);
    249249    Number ret = 0.0;
    250     DBG_PRINT((1, "x.Tag = %d\n", x.GetTag()));
     250    DBG_PRINT((2, "x.Tag = %d\n", x.GetTag()));
    251251    if (!f_cache_.GetCachedResult1Dep(ret, &x)) {
    252252      f_evals_++;
    253253      SmartPtr<const Vector> unscaled_x = NLP_scaling()->unapply_vector_scaling_x(&x);
    254254      bool success = nlp_->Eval_f(*unscaled_x, ret);
     255      DBG_PRINT((1, "success = %d ret = %e\n", success, ret));
    255256      ASSERT_EXCEPTION(success && FiniteNumber(ret), Eval_Error,
    256257                       "Error evaluating the objective function");
     
    284285  SmartPtr<const Vector> OrigIpoptNLP::c(const Vector& x)
    285286  {
    286     SmartPtr<Vector> unscaled_c;
    287287    SmartPtr<const Vector> retValue;
    288     if (!c_cache_.GetCachedResult1Dep(retValue, &x)) {
     288    SmartPtr<const Vector> dep;
     289    if (c_space_->Dim()==0) {
     290      // We do this caching of an empty vector so that the returned
     291      // Vector has always the same tag (this might make a difference
     292      // in cases where only the constraints are supposed to change...
     293      dep = NULL;
     294      if (!c_cache_.GetCachedResult1Dep(retValue, GetRawPtr(dep))) {
     295        retValue = c_space_->MakeNew();
     296        c_cache_.AddCachedResult1Dep(retValue, GetRawPtr(dep));
     297      }
     298    }
     299    else {
     300      dep = &x;
     301    }
     302    if (!c_cache_.GetCachedResult1Dep(retValue, GetRawPtr(dep))) {
     303      SmartPtr<Vector> unscaled_c = c_space_->MakeNew();
    289304      c_evals_++;
    290       unscaled_c = c_space_->MakeNew();
    291305      SmartPtr<const Vector> unscaled_x = NLP_scaling()->unapply_vector_scaling_x(&x);
    292306      bool success = nlp_->Eval_c(*unscaled_x, *unscaled_c);
     
    294308                       Eval_Error, "Error evaluating the equality constraints");
    295309      retValue = NLP_scaling()->apply_vector_scaling_c(ConstPtr(unscaled_c));
    296       c_cache_.AddCachedResult1Dep(retValue, &x);
     310      c_cache_.AddCachedResult1Dep(retValue, GetRawPtr(dep));
    297311    }
    298312
     
    303317  {
    304318    DBG_START_METH("OrigIpoptNLP::d", dbg_verbosity);
    305     SmartPtr<Vector> unscaled_d;
    306319    SmartPtr<const Vector> retValue;
    307     if (!d_cache_.GetCachedResult1Dep(retValue, &x)) {
     320    SmartPtr<const Vector> dep;
     321    if (d_space_->Dim()==0) {
     322      // We do this caching of an empty vector so that the returned
     323      // Vector has always the same tag (this might make a difference
     324      // in cases where only the constraints are supposed to change...
     325      dep = NULL;
     326      if (!d_cache_.GetCachedResult1Dep(retValue, GetRawPtr(dep))) {
     327        retValue = d_space_->MakeNew();
     328        d_cache_.AddCachedResult1Dep(retValue, GetRawPtr(dep));
     329      }
     330    }
     331    else {
     332      dep = &x;
     333    }
     334    if (!d_cache_.GetCachedResult1Dep(retValue, GetRawPtr(dep))) {
    308335      d_evals_++;
    309       unscaled_d = d_space_->MakeNew();
     336      SmartPtr<Vector> unscaled_d = d_space_->MakeNew();
    310337
    311338      DBG_PRINT_VECTOR(2, "scaled_x", x);
     
    316343                       Eval_Error, "Error evaluating the inequality constraints");
    317344      retValue = NLP_scaling()->apply_vector_scaling_d(ConstPtr(unscaled_d));
    318       d_cache_.AddCachedResult1Dep(retValue, &x);
     345      d_cache_.AddCachedResult1Dep(retValue, GetRawPtr(dep));
    319346    }
    320347
     
    324351  SmartPtr<const Matrix> OrigIpoptNLP::jac_c(const Vector& x)
    325352  {
    326     SmartPtr<Matrix> unscaled_jac_c;
    327353    SmartPtr<const Matrix> retValue;
    328     if (!jac_c_cache_.GetCachedResult1Dep(retValue, &x)) {
     354    SmartPtr<const Vector> dep;
     355    if (c_space_->Dim()==0) {
     356      // We do this caching of an empty vector so that the returned
     357      // Matrix has always the same tag
     358      dep = NULL;
     359      if (!jac_c_cache_.GetCachedResult1Dep(retValue, GetRawPtr(dep))) {
     360        retValue = jac_c_space_->MakeNew();
     361        jac_c_cache_.AddCachedResult1Dep(retValue, GetRawPtr(dep));
     362      }
     363    }
     364    else {
     365      dep = &x;
     366    }
     367    if (!jac_c_cache_.GetCachedResult1Dep(retValue, GetRawPtr(dep))) {
    329368      jac_c_evals_++;
    330       unscaled_jac_c = jac_c_space_->MakeNew();
     369      SmartPtr<Matrix> unscaled_jac_c = jac_c_space_->MakeNew();
    331370
    332371      SmartPtr<const Vector> unscaled_x = NLP_scaling()->unapply_vector_scaling_x(&x);
     
    334373      ASSERT_EXCEPTION(success, Eval_Error, "Error evaluating the jacobian of the equality constraints");
    335374      retValue = NLP_scaling()->apply_jac_c_scaling(ConstPtr(unscaled_jac_c));
    336       jac_c_cache_.AddCachedResult1Dep(retValue, &x);
     375      jac_c_cache_.AddCachedResult1Dep(retValue, GetRawPtr(dep));
    337376    }
    338377
     
    342381  SmartPtr<const Matrix> OrigIpoptNLP::jac_d(const Vector& x)
    343382  {
    344     SmartPtr<Matrix> unscaled_jac_d;
    345383    SmartPtr<const Matrix> retValue;
    346     if (!jac_d_cache_.GetCachedResult1Dep(retValue, &x)) {
     384    SmartPtr<const Vector> dep;
     385    if (d_space_->Dim()==0) {
     386      // We do this caching of an empty vector so that the returned
     387      // Matrix has always the same tag
     388      dep = NULL;
     389      if (!jac_d_cache_.GetCachedResult1Dep(retValue, GetRawPtr(dep))) {
     390        retValue = jac_d_space_->MakeNew();
     391        jac_d_cache_.AddCachedResult1Dep(retValue, GetRawPtr(dep));
     392      }
     393    }
     394    else {
     395      dep = &x;
     396    }
     397    if (!jac_d_cache_.GetCachedResult1Dep(retValue, GetRawPtr(dep))) {
    347398      jac_d_evals_++;
    348       unscaled_jac_d = jac_d_space_->MakeNew();
     399      SmartPtr<Matrix> unscaled_jac_d = jac_d_space_->MakeNew();
    349400
    350401      SmartPtr<const Vector> unscaled_x = NLP_scaling()->unapply_vector_scaling_x(&x);
     
    352403      ASSERT_EXCEPTION(success, Eval_Error, "Error evaluating the jacobian of the inequality constraints");
    353404      retValue = NLP_scaling()->apply_jac_d_scaling(ConstPtr(unscaled_jac_d));
    354       jac_d_cache_.AddCachedResult1Dep(retValue, &x);
     405      jac_d_cache_.AddCachedResult1Dep(retValue, GetRawPtr(dep));
    355406    }
    356407
  • branches/dev/Algorithm/IpRestoFilterConvCheck.cpp

    r411 r416  
    4242  void RestoFilterConvergenceCheck::RegisterOptions(SmartPtr<RegisteredOptions> roptions)
    4343  {
    44     roptions->AddBoundedNumberOption("kappa_resto", "???",
    45                                      0.0, false, 1.0, true, 0.9);
     44    roptions->AddBoundedNumberOption(
     45      "required_infeasibility_reduction",
     46      "Required reduction of infeasibility before restoration phase is left.",
     47      0.0, false, 1.0, true,
     48      0.9,
     49      "The restoration phase algorithm is performed, until a point is found "
     50      "that is acceptable to the filter and the infeasibility has been "
     51      "reduced by at least the fraction given by this option.");
    4652  }
    4753
     
    5056  {
    5157    DBG_ASSERT(orig_filter_line_search_ && "Need to call RestoFilterConvergenceCheck::SetOrigFilterLineSearch before Initialize");
    52     options.GetNumericValue("kappa_resto", kappa_resto_, prefix);
     58    options.GetNumericValue("required_infeasibility_reduction", kappa_resto_, prefix);
    5359    options.GetIntegerValue("max_iter", maximum_iters_, prefix);
    5460
     
    100106    Number orig_theta_max = Max(kappa_resto_*orig_curr_theta,
    101107                                1.e1*Min(orig_ip_data->tol(),
    102                                          orig_ip_data->primal_inf_tol()));
     108                                         constr_viol_tol_));
    103109
    104110    if (first_resto_iter_) {
     
    148154
    149155      status = OptimalityErrorConvergenceCheck::CheckConvergence();
    150       if (status == CONVERGED) {
     156      if (status == CONVERGED || status == CONVERGED_TO_ACCEPTABLE_POINT) {
    151157        Number orig_trial_primal_inf =
    152158          orig_ip_cq->trial_primal_infeasibility(NORM_MAX);
  • branches/dev/Algorithm/IpRestoIpoptNLP.cpp

    r414 r416  
    2222{
    2323  DBG_SET_VERBOSITY(0);
     24
     25  DefineIpoptType(RestoIpoptNLP);
    2426
    2527  RestoIpoptNLP::RestoIpoptNLP(IpoptNLP& orig_ip_nlp,
     
    3941  {}
    4042
     43  void RestoIpoptNLP::RegisterOptions(SmartPtr<RegisteredOptions> roptions)
     44  {
     45    roptions->AddStringOption2(
     46      "evaluate_orig_obj_at_resto_trial",
     47      "Determines if the original objective function should be evalutated at restoration phase trial points.",
     48      "yes",
     49      "no", "skip evaluation",
     50      "yes", "evaluate at every trial point",
     51      "Setting this option to true makes the restoration phase algorithm "
     52      "evaluate the objective function of the original problem at every trial "
     53      "point encountered during the restoration phase, even if this value is "
     54      "not required.  In this way, it is guaranteed that the original "
     55      "objective function can be evaluated without error at all accepted "
     56      "iterates; otherwise the algorithm might fail at a point where the "
     57      "restoration phase accepts an iterate that is good for the restoration "
     58      "phase problem, but not the original problem.  On the other hand, if "
     59      "the evaluation of the original objective is expensive, this might be "
     60      "costly.");
     61  }
     62
    4163  bool RestoIpoptNLP::Initialize(const Journalist& jnlst,
    4264                                 const OptionsList& options,
    4365                                 const std::string& prefix)
    4466  {
     67    options.GetBoolValue("evaluate_orig_obj_at_resto_trial",
     68                         evaluate_orig_obj_at_resto_trial_, prefix);
     69
    4570    initialized_ = true;
    46     return true;
     71    return IpoptNLP::Initialize(jnlst, options, prefix);
    4772  }
    4873
     
    281306    // All remaining blocks are zero'ed out
    282307
     308    SmartPtr<const MatrixSpace> scaled_jac_c_space;
     309    SmartPtr<const MatrixSpace> scaled_jac_d_space;
     310    SmartPtr<const SymMatrixSpace> scaled_h_space;
     311    NLP_scaling()->DetermineScaling(GetRawPtr(x_space_),
     312                                    c_space_, d_space_,
     313                                    GetRawPtr(jac_c_space_),
     314                                    GetRawPtr(jac_d_space_),
     315                                    GetRawPtr(h_space_),
     316                                    scaled_jac_c_space, scaled_jac_d_space,
     317                                    scaled_h_space);
     318    // For now we assume that no scaling is done inside the NLP_Scaling
     319    DBG_ASSERT(scaled_jac_c_space == jac_c_space_);
     320    DBG_ASSERT(scaled_jac_d_space == jac_d_space_);
     321    DBG_ASSERT(scaled_h_space == h_space_);
     322
    283323    ///////////////////////////
    284324    // Create the bound data //
     
    411451
    412452    ret += ret2;
     453
     454    // We evaluate also the objective function for the original
     455    // problem here.  This might be wasteful, but it will detect if
     456    // the original objective function cannot be evaluated at the
     457    // trial point in the restoration phase
     458    if (evaluate_orig_obj_at_resto_trial_) {
     459      Number orig_f = orig_ip_nlp_->f(*x_only);
     460    }
     461
    413462    return ret;
    414463  }
  • branches/dev/Algorithm/IpRestoIpoptNLP.hpp

    r414 r416  
    2222namespace Ipopt
    2323{
     24
     25  DeclareIpoptType(RestoIpoptNLP);
    2426
    2527  /** This class maps the traditional NLP into
     
    262264    //@}
    263265
     266    /** Methods for IpoptType */
     267    //@{
     268    /** Called by IpoptType to register the options */
     269    static void RegisterOptions(SmartPtr<RegisteredOptions> roptions);
     270    //@}
     271
    264272  private:
    265273    /** @name Pointers for the original NLP information. */
     
    370378    //@}
    371379
     380    /** Flag indicating if evalution of the objective should be
     381     *  performed for every restoration phase objective function
     382     *  evaluation. */
     383    bool evaluate_orig_obj_at_resto_trial_;
     384
    372385    /** Flag indicating if initialization method has been called */
    373386    bool initialized_;
  • branches/dev/Algorithm/IpRestoMinC_1Nrm.cpp

    r413 r416  
    120120      // that we do not return from the restoration phase is the
    121121      // problem is infeasible
    122       actual_resto_options->SetNumericValue("resto.kappa_resto", 1e-3);
     122      actual_resto_options->SetNumericValue("resto.required_infeasibility_reduction", 1e-3);
    123123    }
    124124    else if(square_problem) {
     
    126126      // If this is a square problem, the want the restoration phase
    127127      // never to be left until the problem is converged
    128       actual_resto_options->SetNumericValue("resto.kappa_resto", 0.);
     128      actual_resto_options->SetNumericValue("resto.required_infeasibility_reduction", 0.);
    129129    }
    130130
Note: See TracChangeset for help on using the changeset viewer.