MoochoPack : Framework for Large-Scale Optimization Algorithms Version of the Day
MoochoPack_UpdateBarrierParameter_Step.cpp
00001 // @HEADER
00002 // ***********************************************************************
00003 // 
00004 // Moocho: Multi-functional Object-Oriented arCHitecture for Optimization
00005 //                  Copyright (2003) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov) 
00038 // 
00039 // ***********************************************************************
00040 // @HEADER
00041 
00042 #include <ostream>
00043 #include <typeinfo>
00044 #include <iostream>
00045 #include <math.h>
00046 
00047 #include "MoochoPack_UpdateBarrierParameter_Step.hpp"
00048 #include "MoochoPack_moocho_algo_conversion.hpp"
00049 #include "IterationPack_print_algorithm_step.hpp"
00050 #include "Teuchos_dyn_cast.hpp"
00051 #include "MoochoPack_IpState.hpp"
00052 #include "AbstractLinAlgPack_MatrixSymDiagStd.hpp"
00053 #include "AbstractLinAlgPack_VectorStdOps.hpp"
00054 #include "AbstractLinAlgPack_VectorAuxiliaryOps.hpp"
00055 #include "Teuchos_Assert.hpp"
00056 
00057 #include "OptionsFromStreamPack_StringToBool.hpp"
00058 
00059 #define min(a,b) ( (a < b) ? a : b )
00060 #define max(a,b) ( (a > b) ? a : b )
00061 
00062 namespace MoochoPack {
00063 
00064 UpdateBarrierParameter_Step::UpdateBarrierParameter_Step(
00065   const value_type init_barrier_parameter,
00066   const value_type tau_mu,
00067   const value_type theta_mu,
00068   const value_type tau_epsilon,
00069   const value_type theta_epsilon,
00070   const value_type e_tol_max
00071   )
00072   :
00073   init_barrier_parameter_(init_barrier_parameter),
00074   tau_mu_(tau_mu),
00075   theta_mu_(theta_mu),
00076   tau_epsilon_(tau_epsilon),
00077   theta_epsilon_(theta_epsilon),
00078   e_tol_max_(e_tol_max)
00079   {}
00080   
00081 
00082 bool UpdateBarrierParameter_Step::do_step(
00083   Algorithm& _algo, poss_type step_poss, IterationPack::EDoStepType type
00084   ,poss_type assoc_step_poss
00085   )
00086   {
00087   using Teuchos::dyn_cast;
00088   using IterationPack::print_algorithm_step;
00089 
00090   NLPAlgo            &algo   = dyn_cast<NLPAlgo>(_algo);
00091   IpState             &s      = dyn_cast<IpState>(_algo.state());
00092   NLP                 &nlp    = algo.nlp();
00093 
00094   if (!nlp.is_initialized())
00095     {
00096     nlp.initialize(false);
00097     }
00098   
00099   EJournalOutputLevel olevel = algo.algo_cntr().journal_output_level();
00100   std::ostream& out = algo.track().journal_out();
00101   
00102   // print step header.
00103   if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) 
00104     {
00105     using IterationPack::print_algorithm_step;
00106     print_algorithm_step( _algo, step_poss, type, assoc_step_poss, out );
00107     }
00108 
00109   
00111   // Get iteration quantities
00113   IterQuantityAccess<value_type>& e_tol_iq = s.e_tol();
00114   IterQuantityAccess<value_type>& mu_iq = s.barrier_parameter();
00115 
00117   // Check values and initialize, if necessary
00119   /*  if (mu_iq.last_updated() == IterQuantity::NONE_UPDATED)
00120     {
00121     mu_iq.set_k(0) = init_barrier_parameter_;
00122     e_tol_iq.set_k(0) = Calculate_e_tol(mu_iq.get_k(0));
00123 
00124     if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) 
00125       {
00126       out << "\nInitializing barrier parameter (mu) and sub problem tolerance (e_tol) ...\n";
00127       }
00128     }
00129     else*/
00130     {
00132     // if tolerance is satisfied, calculate new barrier_parameter 
00133     //  and e_tol, otherwise update mu and e_tol from last iter
00135     const value_type opt_err = s.opt_kkt_err().get_k(0);
00136     const value_type feas_err = s.feas_kkt_err().get_k(0);
00137     const value_type comp_err_mu = s.comp_err_mu().get_k(0);
00138 
00139     const value_type mu_km1 = mu_iq.get_k(-1);
00140     if (e_tol_iq.last_updated() == IterQuantity::NONE_UPDATED)
00141       {
00142       // First time through, don't let mu change
00143       mu_iq.set_k(0,-1);
00144       e_tol_iq.set_k(0) = Calculate_e_tol(mu_iq.get_k(0));
00145           }
00146     else
00147       {
00148       const value_type e_tol_km1 = e_tol_iq.get_k(-1);
00149       bool sub_prob_converged = (opt_err < e_tol_km1 && feas_err < e_tol_km1 && comp_err_mu < e_tol_km1); 
00150       if (sub_prob_converged)
00151         {
00152         // Calculate new mu and e_tol
00153         value_type& mu_k = mu_iq.set_k(0);
00154         mu_k = min(tau_mu_*mu_km1, pow(mu_km1, theta_mu_));
00155         e_tol_iq.set_k(0) = Calculate_e_tol(mu_k);
00156         
00157         if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) 
00158           {
00159           out << "\nSub-problem converged!\n"
00160             << " Updating barrier parameter (mu) and sub problem tolerance (e_tol) ...\n";
00161           }
00162         }
00163       else
00164         {
00165         if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) 
00166           {
00167           out << "\nSub-problem not-converged!\n"
00168             << " Keeping existing barrier parameter (mu) and sub problem tolerance (e_tol) ...\n";
00169           }
00170         mu_iq.set_k(0,-1);
00171         e_tol_iq.set_k(0,-1);
00172         }
00173       }
00174 
00175     if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) 
00176       {
00177       out << "\nbarrier_parameter (mu) = " << mu_iq.get_k(0)
00178         << "\nsub problem tolerance (e_tol) = " << e_tol_iq.get_k(0)  << std::endl;
00179       }
00180     }
00181 
00182   return true;
00183   }
00184 
00185 
00186 void UpdateBarrierParameter_Step::print_step(
00187   const Algorithm& _algo, poss_type step_poss, IterationPack::EDoStepType type
00188   ,poss_type assoc_step_poss, std::ostream& out, const std::string& L
00189   ) const
00190   {
00191   //const NLPAlgo   &algo = rsqp_algo(_algo);
00192   //const NLPAlgoState  &s    = algo.rsqp_state();
00193   out << L << "# Update the interior point barrier parameter (mu)\n"
00194     << L << "if (KKTerror < e_tol) then\n"
00195     << L << "   mu_kp1 = min(tau_mu*mu_k, mu_k^theta_mu)\n"
00196     << L << "   e_tol_kp1 = min(e_tol_max, tau_epsilon*min(mu_k, mu_k^theta_epsilon))\n"
00197     << L << "else\n"
00198     << L << "   mu_kp1 = mu_k\n"
00199     << L << "   e_tol_kp1 = e_tol_k\n"
00200     << L << "end;\n";
00201   }
00202 
00203 value_type UpdateBarrierParameter_Step::Calculate_e_tol(value_type mu)
00204   { 
00205   value_type e_tol = tau_epsilon_*min(mu, pow(mu, theta_epsilon_));
00206   e_tol = min(e_tol_max_, e_tol);
00207 
00208   return e_tol;
00209   }
00210 
00211 
00212 namespace {
00213 
00214 const int local_num_options = 5;
00215 
00216 enum local_EOptions 
00217   {
00218     TAU_MU,
00219     THETA_MU,
00220     TAU_EPSILON,
00221     THETA_EPSILON,
00222     E_TOL_MAX
00223   };
00224 
00225 const char* local_SOptions[local_num_options] = 
00226   {
00227     "tau_mu",
00228     "theta_mu",
00229     "tau_epsilon",
00230     "theta_epsilon",
00231     "e_tol_max"
00232   };
00233 
00234 }
00235 
00236  
00237 UpdateBarrierParameter_StepSetOptions::UpdateBarrierParameter_StepSetOptions(
00238   UpdateBarrierParameter_Step* target
00239   , const char opt_grp_name[] )
00240   :
00241   OptionsFromStreamPack::SetOptionsFromStreamNode(
00242     opt_grp_name, local_num_options, local_SOptions ),
00243   OptionsFromStreamPack::SetOptionsToTargetBase< UpdateBarrierParameter_Step >( target )
00244   {
00245   }
00246 
00247 void UpdateBarrierParameter_StepSetOptions::setOption( 
00248   int option_num, const std::string& option_value )
00249   {
00250   using OptionsFromStreamPack::StringToBool;
00251   
00252   typedef UpdateBarrierParameter_Step target_t;
00253   switch( (local_EOptions)option_num ) 
00254     {
00255     case TAU_MU:
00256       target().tau_mu(std::atof(option_value.c_str()));
00257       break;
00258     case THETA_MU:
00259       target().theta_mu(std::atof(option_value.c_str()));
00260       break;
00261     case TAU_EPSILON:
00262       target().tau_epsilon(std::atof(option_value.c_str()));
00263       break;
00264     case THETA_EPSILON:
00265       target().theta_epsilon(std::atof(option_value.c_str()));
00266       break;
00267     case E_TOL_MAX:
00268       target().e_tol_max(std::atof(option_value.c_str()));
00269       break;
00270     default:
00271       TEUCHOS_TEST_FOR_EXCEPT(true);  // Local error only?
00272     }
00273   }
00274 
00275 } // end namespace MoochoPack 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends