MoochoPack : Framework for Large-Scale Optimization Algorithms Version of the Day
MoochoPack_PreProcessBarrierLineSearch_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 // 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 "AbstractLinAlgPack_assert_print_nan_inf.hpp"
00035 #include "AbstractLinAlgPack_VectorAuxiliaryOps.hpp"
00036 #include "AbstractLinAlgPack_MatrixSymDiagStd.hpp"
00037 #include "AbstractLinAlgPack_VectorOut.hpp"
00038 #include "AbstractLinAlgPack_LinAlgOpPack.hpp"
00039 #include "NLPInterfacePack_NLPBarrier.hpp"
00040 #include "MoochoPack_PreProcessBarrierLineSearch_Step.hpp"
00041 #include "MoochoPack_IpState.hpp"
00042 #include "MoochoPack_moocho_algo_conversion.hpp"
00043 #include "IterationPack_print_algorithm_step.hpp"
00044 #include "Teuchos_dyn_cast.hpp"
00045 #include "Teuchos_TestForException.hpp"
00046 
00047 #define min(a,b) ( (a < b) ? a : b )
00048 #define max(a,b) ( (a > b) ? a : b )
00049 
00050 namespace MoochoPack {
00051 
00052 PreProcessBarrierLineSearch_Step::PreProcessBarrierLineSearch_Step(
00053   Teuchos::RCP<NLPInterfacePack::NLPBarrier> barrier_nlp,
00054   const value_type tau_boundary_frac
00055   )
00056   :
00057   tau_boundary_frac_(tau_boundary_frac),
00058   barrier_nlp_(barrier_nlp),
00059   filter_(FILTER_IQ_STRING)
00060 {
00061   TEST_FOR_EXCEPTION(
00062     !barrier_nlp_.get(),
00063     std::logic_error,
00064     "PreProcessBarrierLineSearch_Step given NULL NLPBarrier."
00065     );
00066 }
00067   
00068 
00069 bool PreProcessBarrierLineSearch_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     using AbstractLinAlgPack::assert_print_nan_inf;
00077   using AbstractLinAlgPack::fraction_to_boundary;
00078   using AbstractLinAlgPack::fraction_to_zero_boundary;
00079   using LinAlgOpPack::Vp_StV;
00080 
00081   NLPAlgo            &algo   = dyn_cast<NLPAlgo>(_algo);
00082   IpState             &s      = dyn_cast<IpState>(_algo.state());
00083   NLP                 &nlp    = algo.nlp();
00084 
00085   EJournalOutputLevel olevel = algo.algo_cntr().journal_output_level();
00086   std::ostream& out = algo.track().journal_out();
00087   
00088   // print step header.
00089   if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) 
00090   {
00091     using IterationPack::print_algorithm_step;
00092     print_algorithm_step( _algo, step_poss, type, assoc_step_poss, out );
00093   }
00094 
00095   const value_type& mu_k = s.barrier_parameter().get_k(0);
00096 
00097   // if using filter and u changed, clear filter
00098   if (filter_.exists_in(s))
00099   {
00100     if ( s.barrier_parameter().updated_k(-1) )
00101     {
00102       const value_type mu_km1 = s.barrier_parameter().get_k(-1);
00103       if (mu_k != mu_km1)
00104       {
00105         if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) 
00106         {
00107           out << "\nBarrier Parameter changed - resetting the filter ...\n";
00108         }
00109         // reset the filter
00110         MoochoPack::Filter_T &filter_k = filter_(s).set_k(0);
00111         filter_k.clear();
00112       }
00113     }
00114   }
00115 
00116   // Update the barrier parameter in the NLP
00117   barrier_nlp_->mu(s.barrier_parameter().get_k(0));
00118     
00119   // Calculate the barrier k terms
00120   barrier_nlp_->set_Gf( &(s.grad_barrier_obj().set_k(0)) );
00121   barrier_nlp_->calc_Gf(s.x().get_k(0), true);
00122 
00123   barrier_nlp_->set_f( &(s.barrier_obj().set_k(0)) );
00124   barrier_nlp_->calc_f(s.x().get_k(0), true);
00125 
00126 
00127   // Calculate the k+1 terms
00128   // Get iteration quantities...
00129   value_type& alpha_k = s.alpha().set_k(0);
00130   value_type& alpha_vl_k = s.alpha_vl().set_k(0);
00131   value_type& alpha_vu_k = s.alpha_vu().set_k(0);
00132 
00133   const Vector& x_k = s.x().get_k(0);
00134   VectorMutable& x_kp1 = s.x().set_k(+1);
00135 
00136   const Vector& d_k = s.d().get_k(0);
00137   const Vector& dvl_k = s.dvl().get_k(0);
00138   const Vector& dvu_k = s.dvu().get_k(0);
00139 
00140   const Vector& vl_k = s.Vl().get_k(0).diag();
00141   VectorMutable& vl_kp1 = s.Vl().set_k(+1).diag();
00142 
00143   const Vector& vu_k = s.Vu().get_k(0).diag();
00144   VectorMutable& vu_kp1 = s.Vu().set_k(+1).diag();
00145 
00146   alpha_k = fraction_to_boundary(
00147     tau_boundary_frac_, 
00148     x_k, 
00149     d_k,
00150     nlp.xl(),
00151     nlp.xu()
00152     );
00153 
00154   alpha_vl_k = fraction_to_zero_boundary(
00155     tau_boundary_frac_,
00156     vl_k,
00157     dvl_k
00158     );
00159 
00160   alpha_vu_k = fraction_to_zero_boundary(
00161     tau_boundary_frac_,
00162     vu_k,
00163     dvu_k
00164     );
00165 
00166   TEST_FOR_EXCEPT( !( alpha_k <= 1.0 && alpha_vl_k <= 1.0 && alpha_vu_k <= 1.0 ) );
00167   TEST_FOR_EXCEPT( !( alpha_k >= 0.0 && alpha_vl_k >= 0.0 && alpha_vu_k >= 0.0 ) );
00168 
00169   x_kp1 = x_k;
00170   Vp_StV( &x_kp1, alpha_k, d_k);
00171 
00172   alpha_vl_k = alpha_vu_k = min(alpha_vl_k, alpha_vu_k);
00173 
00174   vl_kp1 = vl_k;
00175   Vp_StV( &vl_kp1, alpha_vl_k, dvl_k);
00176 
00177   vu_kp1 = vu_k;
00178   Vp_StV( &vu_kp1, alpha_vu_k, dvu_k);
00179  
00180 
00181     IterQuantityAccess<VectorMutable>
00182     *c_iq   = nlp.m() > 0 ? &s.c() : NULL;
00183 
00184     if (assert_print_nan_inf(x_kp1, "x", true, NULL))
00185   {
00186     // Calcuate f and c at the new point.
00187     barrier_nlp_->unset_quantities();
00188     barrier_nlp_->set_f( &s.barrier_obj().set_k(+1) );
00189     if (c_iq) {
00190       barrier_nlp_->set_c( &c_iq->set_k(+1) );
00191       barrier_nlp_->calc_c( x_kp1, true );
00192     }
00193     barrier_nlp_->calc_f( x_kp1, false ); 
00194     barrier_nlp_->unset_quantities();
00195   }
00196   
00197   if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) )
00198   {
00199     out << "\nalpha_vl_k = " << alpha_vl_k
00200       << "\nalpha_vu_k = " << alpha_vu_k
00201       << "\nalpha_k    = " << alpha_k
00202       << std::endl;
00203   }
00204 
00205   if( static_cast<int>(olevel) >= static_cast<int>(PRINT_VECTORS) ) 
00206   {
00207     out << "\nvl_kp1 = \n" << vl_kp1
00208       << "\nvu_kp1 = \n" << vu_kp1
00209       << "\nx_kp1 = \n" << x_kp1;
00210   }
00211 
00212   return true;
00213 }
00214 
00215 void PreProcessBarrierLineSearch_Step::print_step(
00216   const Algorithm& _algo, poss_type step_poss, IterationPack::EDoStepType type
00217   ,poss_type assoc_step_poss, std::ostream& out, const std::string& L
00218   ) const
00219 {
00220   //const NLPAlgo   &algo = rsqp_algo(_algo);
00221   //const NLPAlgoState  &s    = algo.rsqp_state();
00222   out << L << "*** calculate alpha max by the fraction to boundary rule\n"
00223     << L << "ToDo: Complete documentation\n";
00224 }
00225 
00226 namespace {
00227 
00228 const int local_num_options = 1;
00229 
00230 enum local_EOptions 
00231 {
00232   TAU_BOUNDARY_FRAC
00233 };
00234 
00235 const char* local_SOptions[local_num_options] = 
00236 {
00237   "tau_boundary_frac"
00238 };
00239 
00240 }
00241 
00242  
00243 PreProcessBarrierLineSearch_StepSetOptions::PreProcessBarrierLineSearch_StepSetOptions(
00244   PreProcessBarrierLineSearch_Step* target
00245   , const char opt_grp_name[] )
00246   :
00247   OptionsFromStreamPack::SetOptionsFromStreamNode(
00248     opt_grp_name, local_num_options, local_SOptions ),
00249   OptionsFromStreamPack::SetOptionsToTargetBase< PreProcessBarrierLineSearch_Step >( target )
00250 {}
00251 
00252 void PreProcessBarrierLineSearch_StepSetOptions::setOption( 
00253   int option_num, const std::string& option_value )
00254 {
00255   typedef PreProcessBarrierLineSearch_Step target_t;
00256   switch( (local_EOptions)option_num ) 
00257   {
00258     case TAU_BOUNDARY_FRAC:
00259       target().tau_boundary_frac(std::atof(option_value.c_str()));
00260       break;
00261     default:
00262       TEST_FOR_EXCEPT(true);  // Local error only?
00263   }
00264 }
00265 
00266 } // end namespace MoochoPack 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends