MoochoPack_UpdateBarrierParameter_Step.cpp

Go to the documentation of this file.
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 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Roscoe A. Bartlett (rabartl@sandia.gov) 
00025 // 
00026 // ***********************************************************************
00027 // @HEADER
00028 
00029 #include <ostream>
00030 #include <typeinfo>
00031 #include <iostream>
00032 #include <math.h>
00033 
00034 #include "MoochoPack_UpdateBarrierParameter_Step.hpp"
00035 #include "MoochoPack_moocho_algo_conversion.hpp"
00036 #include "IterationPack_print_algorithm_step.hpp"
00037 #include "Teuchos_dyn_cast.hpp"
00038 #include "MoochoPack_IpState.hpp"
00039 #include "AbstractLinAlgPack_MatrixSymDiagStd.hpp"
00040 #include "AbstractLinAlgPack_VectorStdOps.hpp"
00041 #include "AbstractLinAlgPack_VectorAuxiliaryOps.hpp"
00042 #include "Teuchos_TestForException.hpp"
00043 
00044 #include "OptionsFromStreamPack_StringToBool.hpp"
00045 
00046 #define min(a,b) ( (a < b) ? a : b )
00047 #define max(a,b) ( (a > b) ? a : b )
00048 
00049 namespace MoochoPack {
00050 
00051 UpdateBarrierParameter_Step::UpdateBarrierParameter_Step(
00052   const value_type init_barrier_parameter,
00053   const value_type tau_mu,
00054   const value_type theta_mu,
00055   const value_type tau_epsilon,
00056   const value_type theta_epsilon,
00057   const value_type e_tol_max
00058   )
00059   :
00060   init_barrier_parameter_(init_barrier_parameter),
00061   tau_mu_(tau_mu),
00062   theta_mu_(theta_mu),
00063   tau_epsilon_(tau_epsilon),
00064   theta_epsilon_(theta_epsilon),
00065   e_tol_max_(e_tol_max)
00066   {}
00067   
00068 
00069 bool UpdateBarrierParameter_Step::do_step(
00070   Algorithm& _algo, poss_type step_poss, IterationPack::EDoStepType type
00071   ,poss_type assoc_step_poss
00072   )
00073   {
00074   using Teuchos::dyn_cast;
00075   using IterationPack::print_algorithm_step;
00076 
00077   NLPAlgo            &algo   = dyn_cast<NLPAlgo>(_algo);
00078   IpState             &s      = dyn_cast<IpState>(_algo.state());
00079   NLP                 &nlp    = algo.nlp();
00080 
00081   if (!nlp.is_initialized())
00082     {
00083     nlp.initialize(false);
00084     }
00085   
00086   EJournalOutputLevel olevel = algo.algo_cntr().journal_output_level();
00087   std::ostream& out = algo.track().journal_out();
00088   
00089   // print step header.
00090   if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) 
00091     {
00092     using IterationPack::print_algorithm_step;
00093     print_algorithm_step( _algo, step_poss, type, assoc_step_poss, out );
00094     }
00095 
00096   
00098   // Get iteration quantities
00100   IterQuantityAccess<value_type>& e_tol_iq = s.e_tol();
00101   IterQuantityAccess<value_type>& mu_iq = s.barrier_parameter();
00102 
00104   // Check values and initialize, if necessary
00106   /*  if (mu_iq.last_updated() == IterQuantity::NONE_UPDATED)
00107     {
00108     mu_iq.set_k(0) = init_barrier_parameter_;
00109     e_tol_iq.set_k(0) = Calculate_e_tol(mu_iq.get_k(0));
00110 
00111     if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) 
00112       {
00113       out << "\nInitializing barrier parameter (mu) and sub problem tolerance (e_tol) ...\n";
00114       }
00115     }
00116     else*/
00117     {
00119     // if tolerance is satisfied, calculate new barrier_parameter 
00120     //  and e_tol, otherwise update mu and e_tol from last iter
00122     const value_type opt_err = s.opt_kkt_err().get_k(0);
00123     const value_type feas_err = s.feas_kkt_err().get_k(0);
00124     const value_type comp_err_mu = s.comp_err_mu().get_k(0);
00125 
00126     const value_type mu_km1 = mu_iq.get_k(-1);
00127     if (e_tol_iq.last_updated() == IterQuantity::NONE_UPDATED)
00128       {
00129       // First time through, don't let mu change
00130       mu_iq.set_k(0,-1);
00131       e_tol_iq.set_k(0) = Calculate_e_tol(mu_iq.get_k(0));
00132           }
00133     else
00134       {
00135       const value_type e_tol_km1 = e_tol_iq.get_k(-1);
00136       bool sub_prob_converged = (opt_err < e_tol_km1 && feas_err < e_tol_km1 && comp_err_mu < e_tol_km1); 
00137       if (sub_prob_converged)
00138         {
00139         // Calculate new mu and e_tol
00140         value_type& mu_k = mu_iq.set_k(0);
00141         mu_k = min(tau_mu_*mu_km1, pow(mu_km1, theta_mu_));
00142         e_tol_iq.set_k(0) = Calculate_e_tol(mu_k);
00143         
00144         if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) 
00145           {
00146           out << "\nSub-problem converged!\n"
00147             << " Updating barrier parameter (mu) and sub problem tolerance (e_tol) ...\n";
00148           }
00149         }
00150       else
00151         {
00152         if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) 
00153           {
00154           out << "\nSub-problem not-converged!\n"
00155             << " Keeping existing barrier parameter (mu) and sub problem tolerance (e_tol) ...\n";
00156           }
00157         mu_iq.set_k(0,-1);
00158         e_tol_iq.set_k(0,-1);
00159         }
00160       }
00161 
00162     if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) 
00163       {
00164       out << "\nbarrier_parameter (mu) = " << mu_iq.get_k(0)
00165         << "\nsub problem tolerance (e_tol) = " << e_tol_iq.get_k(0)  << std::endl;
00166       }
00167     }
00168 
00169   return true;
00170   }
00171 
00172 
00173 void UpdateBarrierParameter_Step::print_step(
00174   const Algorithm& _algo, poss_type step_poss, IterationPack::EDoStepType type
00175   ,poss_type assoc_step_poss, std::ostream& out, const std::string& L
00176   ) const
00177   {
00178   //const NLPAlgo   &algo = rsqp_algo(_algo);
00179   //const NLPAlgoState  &s    = algo.rsqp_state();
00180   out << L << "# Update the interior point barrier parameter (mu)\n"
00181     << L << "if (KKTerror < e_tol) then\n"
00182     << L << "   mu_kp1 = min(tau_mu*mu_k, mu_k^theta_mu)\n"
00183     << L << "   e_tol_kp1 = min(e_tol_max, tau_epsilon*min(mu_k, mu_k^theta_epsilon))\n"
00184     << L << "else\n"
00185     << L << "   mu_kp1 = mu_k\n"
00186     << L << "   e_tol_kp1 = e_tol_k\n"
00187     << L << "end;\n";
00188   }
00189 
00190 value_type UpdateBarrierParameter_Step::Calculate_e_tol(value_type mu)
00191   { 
00192   value_type e_tol = tau_epsilon_*min(mu, pow(mu, theta_epsilon_));
00193   e_tol = min(e_tol_max_, e_tol);
00194 
00195   return e_tol;
00196   }
00197 
00198 
00199 namespace {
00200 
00201 const int local_num_options = 5;
00202 
00203 enum local_EOptions 
00204   {
00205     TAU_MU,
00206     THETA_MU,
00207     TAU_EPSILON,
00208     THETA_EPSILON,
00209     E_TOL_MAX
00210   };
00211 
00212 const char* local_SOptions[local_num_options] = 
00213   {
00214     "tau_mu",
00215     "theta_mu",
00216     "tau_epsilon",
00217     "theta_epsilon",
00218     "e_tol_max"
00219   };
00220 
00221 }
00222 
00223  
00224 UpdateBarrierParameter_StepSetOptions::UpdateBarrierParameter_StepSetOptions(
00225   UpdateBarrierParameter_Step* target
00226   , const char opt_grp_name[] )
00227   :
00228   OptionsFromStreamPack::SetOptionsFromStreamNode(
00229     opt_grp_name, local_num_options, local_SOptions ),
00230   OptionsFromStreamPack::SetOptionsToTargetBase< UpdateBarrierParameter_Step >( target )
00231   {
00232   }
00233 
00234 void UpdateBarrierParameter_StepSetOptions::setOption( 
00235   int option_num, const std::string& option_value )
00236   {
00237   using OptionsFromStreamPack::StringToBool;
00238   
00239   typedef UpdateBarrierParameter_Step target_t;
00240   switch( (local_EOptions)option_num ) 
00241     {
00242     case TAU_MU:
00243       target().tau_mu(std::atof(option_value.c_str()));
00244       break;
00245     case THETA_MU:
00246       target().theta_mu(std::atof(option_value.c_str()));
00247       break;
00248     case TAU_EPSILON:
00249       target().tau_epsilon(std::atof(option_value.c_str()));
00250       break;
00251     case THETA_EPSILON:
00252       target().theta_epsilon(std::atof(option_value.c_str()));
00253       break;
00254     case E_TOL_MAX:
00255       target().e_tol_max(std::atof(option_value.c_str()));
00256       break;
00257     default:
00258       TEST_FOR_EXCEPT(true);  // Local error only?
00259     }
00260   }
00261 
00262 } // end namespace MoochoPack 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
Generated on Wed Apr 13 10:11:00 2011 for MOOCHO (Single Doxygen Collection) by  doxygen 1.6.3