MoochoPack_LineSearchWatchDog_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 
00032 #include "MoochoPack_LineSearchWatchDog_Step.hpp"
00033 #include "MoochoPack_MoochoAlgorithmStepNames.hpp"
00034 #include "MoochoPack_moocho_algo_conversion.hpp"
00035 #include "IterationPack_print_algorithm_step.hpp"
00036 #include "ConstrainedOptPack_MeritFuncCalc1DQuadratic.hpp"
00037 #include "ConstrainedOptPack_MeritFuncCalcNLP.hpp"
00038 #include "ConstrainedOptPack_print_vector_change_stats.hpp"
00039 #include "ConstrainedOptPack/src/VectorWithNorms.h"
00040 #include "AbstractLinAlgPack/src/AbstractLinAlgPack_MatrixOp.hpp"
00041 #include "DenseLinAlgPack_DVectorClass.hpp"
00042 #include "DenseLinAlgPack_DVectorOp.hpp"
00043 #include "DenseLinAlgPack_DVectorOut.hpp"
00044 #include "DenseLinAlgPack_LinAlgOpPack.hpp"
00045 
00046 namespace {
00047   const int NORMAL_LINE_SEARCH = -1;
00048 }
00049 
00050 namespace LinAlgOpPack {
00051   using AbstractLinAlgPack::Vp_StMtV;
00052 }
00053 
00054 MoochoPack::LineSearchWatchDog_Step::LineSearchWatchDog_Step(
00055       const direct_line_search_ptr_t& direct_line_search
00056     , const merit_func_ptr_t&     merit_func
00057     , value_type            eta
00058     , value_type            opt_kkt_err_threshold
00059     , value_type            feas_kkt_err_threshold
00060     )
00061   :
00062       direct_line_search_(direct_line_search)
00063     , merit_func_(merit_func)
00064     , eta_(eta)
00065     , opt_kkt_err_threshold_(opt_kkt_err_threshold)
00066     , feas_kkt_err_threshold_(feas_kkt_err_threshold)
00067     , watch_k_(NORMAL_LINE_SEARCH)
00068 {}
00069 
00070 bool MoochoPack::LineSearchWatchDog_Step::do_step(Algorithm& _algo
00071   , poss_type step_poss, IterationPack::EDoStepType type, poss_type assoc_step_poss)
00072 {
00073   using DenseLinAlgPack::norm_inf;
00074   using DenseLinAlgPack::V_VpV;
00075   using DenseLinAlgPack::Vp_StV;
00076   using DenseLinAlgPack::Vt_S;
00077 
00078   using LinAlgOpPack::Vp_V;
00079   using LinAlgOpPack::V_MtV;
00080 
00081   using ConstrainedOptPack::print_vector_change_stats;
00082 
00083   NLPAlgo &algo = rsqp_algo(_algo);
00084   NLPAlgoState  &s    = algo.rsqp_state();
00085   NLP     &nlp  = algo.nlp();
00086 
00087   EJournalOutputLevel olevel = algo.algo_cntr().journal_output_level();
00088   std::ostream& out = algo.track().journal_out();
00089   out << std::boolalpha;
00090 
00091   // print step header.
00092   if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00093     using IterationPack::print_algorithm_step;
00094     print_algorithm_step( algo, step_poss, type, assoc_step_poss, out );
00095   }
00096 
00097   // /////////////////////////////////////////
00098   // Set references to iteration quantities
00099   //
00100   // Set k+1 first then go back to get k to ensure
00101   // we have backward storage.
00102   
00103   DVector
00104     &x_kp1 = s.x().set_k(+1).v();
00105   value_type
00106     &f_kp1 = s.f().set_k(+1);
00107   DVector
00108     &c_kp1 = s.c().set_k(+1).v();
00109 
00110   const value_type
00111     &f_k = s.f().get_k(0);
00112   const DVector
00113     &c_k = s.c().get_k(0).v();
00114   const DVector
00115     &x_k = s.x().get_k(0).v();
00116   const DVector
00117     &d_k = s.d().get_k(0).v();
00118   value_type
00119     &alpha_k = s.alpha().get_k(0);
00120 
00121   // /////////////////////////////////////
00122   // Compute Dphi_k, phi_kp1 and phi_k
00123 
00124   // Dphi_k
00125   const value_type
00126     Dphi_k = merit_func().deriv();
00127   if( Dphi_k >= 0 ) {
00128     throw LineSearchFailure( "LineSearch2ndOrderCorrect_Step::do_step(...) : " 
00129       "Error, d_k is not a descent direction for the merit function " );
00130   }
00131 
00132   // ph_kp1
00133   value_type
00134     &phi_kp1 = s.phi().set_k(+1) = merit_func().value( f_kp1, c_kp1 );
00135 
00136   // Must compute phi(x) at the base point x_k since the penalty parameter may have changed.
00137   const value_type
00138     &phi_k = s.phi().set_k(0) = merit_func().value( f_k, c_k );
00139 
00140   // //////////////////////////////////////
00141   // Setup the calculation merit function
00142 
00143   // Here f_kp1, and c_kp1 are updated at the same time the
00144   // line search is being performed.
00145   nlp.set_f( &f_kp1 );
00146   nlp.set_c( &c_kp1 );
00147   MeritFuncCalcNLP
00148     phi_calc( &merit_func(), &nlp );
00149 
00150   // ////////////////////////////////
00151   // Use Watchdog near the solution
00152 
00153   if( watch_k_ == NORMAL_LINE_SEARCH ) {
00154     const value_type
00155       opt_kkt_err_k = s.opt_kkt_err().get_k(0),
00156       feas_kkt_err_k  = s.feas_kkt_err().get_k(0);
00157     if( opt_kkt_err_k <= opt_kkt_err_threshold() && feas_kkt_err_k <= feas_kkt_err_threshold() ) {
00158       if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00159         out << "\nopt_kkt_err_k = " << opt_kkt_err_k << " <= opt_kkt_err_threshold = "
00160             << opt_kkt_err_threshold() << std::endl
00161           << "\nfeas_kkt_err_k = " << feas_kkt_err_k << " <= feas_kkt_err_threshold = "
00162             << feas_kkt_err_threshold() << std::endl
00163           << "\nSwitching to watchdog linesearch ...\n";
00164       }
00165       watch_k_ = 0;
00166     }
00167   }
00168 
00169   if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00170     out << "\nTrial point:\n"
00171       << "phi_k   = " << phi_k << std::endl
00172       << "Dphi_k  = " << Dphi_k << std::endl
00173       << "phi_kp1 = " << phi_kp1 << std::endl;
00174   }
00175 
00176   bool  ls_success = true,
00177       step_return = true;
00178 
00179   switch( watch_k_ ) {
00180     case 0:
00181     {
00182       // Take  a full step
00183       const value_type phi_cord = phi_k + eta() * Dphi_k;
00184       const bool accept_step = phi_kp1 <= phi_cord;
00185 
00186       if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00187         out << "\n*** Zeroth watchdog iteration:\n"
00188           << "\nphi_kp1 = " << phi_kp1 << ( accept_step ? " <= " : " > " )
00189             << "phi_k + eta * Dphi_k = " << phi_cord << std::endl;
00190       }
00191 
00192       if( phi_kp1 > phi_cord ) {
00193         if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00194           out << "\nAccept this increase for now but watch out next iteration!\n";
00195         }
00196         // Save this initial point
00197         xo_   = x_k;
00198         fo_   = f_k;
00199         nrm_co_ = norm_inf( c_k );
00200         do_   = d_k;
00201         phio_ = phi_k;
00202         Dphio_  = Dphi_k;
00203         phiop1_ = phi_kp1;
00204         // Slip the update of the penalty parameter
00205         const value_type mu_k = s.mu().get_k(0);
00206         s.mu().set_k(+1) = mu_k;
00207         // Move on to the next step in the watchdog procedure
00208         watch_k_ = 1;
00209       }
00210       else {
00211         if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00212           out << "\nAll is good!\n";
00213         }
00214         // watch_k_ stays 0
00215       }
00216       step_return = true;
00217       break;
00218     }
00219     case 1:
00220     {
00221       if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00222         out << "\n*** First watchdog iteration:\n"
00223           << "\nDo a line search to determine x_kp1 = x_k + alpha_k * d_k ...\n";
00224       }
00225       // Now do a line search but and we require some type of reduction
00226       const DVectorSlice xd[2] = { x_k(), d_k() };
00227       MeritFuncCalc1DQuadratic phi_calc_1d( phi_calc, 1, xd, &x_kp1() );
00228       ls_success = direct_line_search().do_line_search( phi_calc_1d, phi_k
00229         , &alpha_k, &phi_kp1
00230         , (int)olevel >= (int)PRINT_ALGORITHM_STEPS ?
00231           &out : static_cast<std::ostream*>(0)  );
00232 
00233       // If the linesearch failed then the rest of the tests will catch this.
00234 
00235       value_type phi_cord = 0;
00236       bool test1, test2;
00237 
00238       if(   ( test1 = ( phi_k <= phio_ ) )
00239         || ( test2 = phi_kp1 <= ( phi_cord = phio_ + eta() * Dphio_ ) )   )
00240       {
00241         // We will accept this step and and move on.
00242         if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00243           out
00244             << "\nphi_k = " << phi_k << ( test1 ? " <= " : " > " )
00245               << "phi_km1 = " << phio_ << std::endl
00246             << "phi_kp1 = " << phi_kp1 << ( test2 ? " <= " : " > " )
00247               << "phi_km1 + eta * Dphi_km1 = " << phi_cord << std::endl
00248             << "This is a sufficent reduction so reset watchdog.\n";
00249         }
00250         watch_k_ = 0;
00251         step_return = true;
00252       }
00253       else if ( ! ( test1 = ( phi_kp1 <= phio_ ) ) ) {
00254         // Even this reduction is no good!
00255         // Go back to original point and do a linesearch from there.
00256         if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00257           out
00258             << "\nphi_kp1 = " << phi_kp1 << " > phi_km1 = " << phio_ << std::endl
00259             << "This is not good reduction in phi so do linesearch from x_km1\n"
00260             << "\n* Go back to x_km1: x_kp1 = x_k - alpha_k * d_k ...\n";
00261         }
00262 
00263         // Go back from x_k to x_km1 for iteration k:
00264         //
00265         // x_kp1 = x_km1
00266         // x_kp1 = x_k - alpha_km1 * d_km1
00267         //
00268         // A negative sign for alpha is an indication that we are backtracking.
00269         //
00270         s.alpha().set_k(0)      = -1.0;
00271         s.d().set_k(0).v()      = do_;
00272         s.f().set_k(+1)       = fo_;
00273 
00274         if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00275           out << "Output iteration k ...\n"
00276             << "k = k+1\n";
00277         }
00278 
00279         // Output these iteration quantities
00280         algo.track().output_iteration( algo );  // k
00281         // Transition to iteration k+1
00282         s.next_iteration();
00283 
00284         // Take the step from x_k = x_km2 to x_kp1 for iteration k (k+1):
00285         //
00286         // x_kp1 = x_km2 + alpha_n * d_km2
00287         // x_kp1 = x_k   + alpha_n * d_km1
00288         // x_kp1 = x_n
00289         //        
00290         if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00291           out << "\n* Take the step from x_k = x_km2 to x_kp1 for iteration k (k+1)\n"
00292             << "Find: x_kp1 = x_k + alpha_k * d_k = x_km2 + alpha_k * d_km2\n ...\n";
00293         }
00294 
00295         // alpha_k = 1.0
00296         value_type &alpha_k = s.alpha().set_k(0) = 1.0;
00297         
00298         // /////////////////////////////////////
00299         // Compute Dphi_k and phi_k
00300 
00301         // x_k
00302         const DVector &x_k                = xo_;
00303 
00304         // d_k
00305         const DVector &d_k = s.d().set_k(0).v()     = do_;
00306 
00307         // Dphi_k
00308         const value_type &Dphi_k            = Dphio_;
00309 
00310         // phi_k
00311         const value_type &phi_k = s.phi().set_k(0)    = phio_;
00312 
00313         // Here f_kp1, and c_kp1 are updated at the same time the
00314         // line search is being performed.
00315         algo.nlp().set_f( &s.f().set_k(+1) );
00316         algo.nlp().set_c( &s.c().set_k(+1).v() );
00317         phi_calc.set_nlp( algo.get_nlp() );
00318 
00319         // ////////////////////////////////////////
00320         // Compute x_xp1 and ph_kp1 for full step
00321 
00322         // x_kp1 = x_k + alpha_k * d_k
00323         DVector &x_kp1 = s.x().set_k(+1).v();
00324         V_VpV( &x_kp1, x_k, d_k );
00325 
00326         // phi_kp1
00327         value_type &phi_kp1 = s.phi().set_k(+1)     = phiop1_;
00328 
00329         const DVectorSlice xd[2] = { x_k(), d_k() };
00330         MeritFuncCalc1DQuadratic phi_calc_1d( phi_calc, 1, xd, &x_kp1() );
00331         ls_success = direct_line_search().do_line_search(
00332             phi_calc_1d, phi_k
00333           , &alpha_k, &phi_kp1
00334           , (int)olevel >= (int)PRINT_ALGORITHM_STEPS ?
00335             &out : static_cast<std::ostream*>(0)  );
00336 
00337         if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00338           out << "\nOutput iteration k (k+1) ...\n"
00339             << "k = k+1 (k+2)\n"
00340             << "Reinitialize watchdog algorithm\n";
00341         }
00342 
00343         // Output these iteration quantities
00344         algo.track().output_iteration( algo );  // (k+1)
00345         // Transition to iteration k+1 (k+2)
00346         s.next_iteration();
00347 
00348         watch_k_ = 0; // Reinitialize the watchdog
00349 
00350         // Any update for k (k+2) should use the last updated value
00351         // which was for k-2 (k) since there is not much info for k-1 (k+1).
00352         // Be careful here and make sure this is square with other steps.
00353 
00354         algo.do_step_next( EvalNewPoint_name );
00355         step_return = false;  // Redirect control
00356       }
00357       else {
00358         if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00359           out
00360             << "phi_kp1 = " << phi_kp1 << " <= phi_km1 = " << phio_ << std::endl
00361             << "\nAccept this step but do a linesearch next iteration!\n";
00362         }
00363         // Slip the update of the penalty parameter
00364         const value_type mu_k = s.mu().get_k(0);
00365         s.mu().set_k(+1) = mu_k;
00366         // Do the last stage of the watchdog procedure next iteration.
00367         watch_k_ = 2;
00368         step_return = true;
00369       }
00370       break;
00371     }
00372     case NORMAL_LINE_SEARCH:
00373     case 2:
00374     {
00375       if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00376         if( watch_k_ == 2 ) {
00377           out << "\n*** Second watchdog iteration:\n"
00378             << "Do a line search to determine x_kp1 = x_k + alpha_k * d_k ...\n";
00379         }
00380         else {
00381           out << "\n*** Normal linesearch:\n"
00382             << "Do a line search to determine x_kp1 = x_k + alpha_k * d_k ...\n";
00383         }
00384       }
00385 
00386       const DVectorSlice xd[2] = { x_k(), d_k() };
00387       MeritFuncCalc1DQuadratic phi_calc_1d( phi_calc, 1, xd, &x_kp1() );
00388       ls_success = direct_line_search().do_line_search( phi_calc_1d, phi_k
00389         , &alpha_k, &phi_kp1
00390         , (int)olevel >= (int)PRINT_ALGORITHM_STEPS ?
00391           &out : static_cast<std::ostream*>(0)  );
00392 
00393       if( watch_k_ == 2 )
00394         watch_k_ = 0;
00395 
00396       step_return = true;
00397       break;
00398     }
00399     default:
00400       TEST_FOR_EXCEPT(true);  // Only local programming error
00401   }
00402 
00403   if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) {
00404     out << "\nalpha    = " << s.alpha().get_k(0) << "\n";
00405     out << "\nphi_kp1 = " << s.phi().get_k(+1) << "\n";
00406   }
00407 
00408   if( static_cast<int>(olevel) >= static_cast<int>(PRINT_VECTORS) ) {
00409     out << "\nd_k = \n" << s.d().get_k(0)();
00410     out << "\nx_kp1 = \n" << s.x().get_k(+1)();
00411   }
00412 
00413   if( !ls_success )
00414     throw LineSearchFailure("LineSearchWatchDog_Step::do_step(): Line search failure");
00415 
00416   return step_return;
00417 
00418 }
00419 
00420 void MoochoPack::LineSearchWatchDog_Step::print_step( const Algorithm& algo
00421   , poss_type step_poss, IterationPack::EDoStepType type, poss_type assoc_step_poss
00422   , std::ostream& out, const std::string& L ) const
00423 {
00424   out << L << "*** Use the Watchdog linesearch when near solution.\n"
00425     << L << "default: opt_kkt_err_threshold = 0.0\n"
00426     << L << "         feas_kkt_err_threshold = 0.0\n"
00427     << L << "         eta = 1.0e-4\n"
00428     << L << "         watch_k = NORMAL_LINE_SEARCH\n"
00429     << L << "begin definition of NLP merit function phi.value(f(x),c(x)):\n";
00430 
00431   merit_func().print_merit_func( out, L + "    " );
00432   
00433   out << L << "end definition\n"
00434     << L << "Dphi_k = phi.deriv()\n"
00435     << L << "if Dphi_k >= 0 then\n"
00436     << L << "    throw line_search_failure\n"
00437     << L << "end\n"
00438     << L << "phi_kp1 = phi_k.value(f_kp1,c_kp1)\n"
00439     << L << "phi_k = phi.value(f_k,c_k)\n"
00440     << L << "if watch_k == NORMAL_LINE_SEARCH then\n"
00441     << L << "    if opt_kkt_err <= opt_kkt_err_threshold\n"
00442     << L << "      and feas_kkt_err <= feas_kkt_err_threshold then\n"
00443     << L << "        *** Start using watchdog from now on!\n"
00444     << L << "        watch_k = 0\n"
00445     << L << "    end\n"
00446     << L << "end\n"
00447     << L << "if watch_k == 0 then\n"
00448     << L << "    *** Zeroth watchdog iteration\n"
00449     << L << "    if phi_kp1 >= phi_k + eta * Dphi_k then\n"
00450     << L << "        *** Accept this increase for now but watch out next iteration!\n"
00451     << L << "        *** Save the first point\n"
00452     << L << "        xo     = x_k\n"
00453     << L << "        fo     = f_k\n"
00454     << L << "        nrm_co = norm_inf_c_k\n"
00455     << L << "        do     = d_k\n"
00456     << L << "        phio   = phi_k\n"
00457     << L << "        Dphio  = Dphi_k\n"
00458     << L << "        phiop1 = phi_kp1\n"
00459     << L << "        *** Skip the update of the penalty parameter next iteration.\n"
00460     << L << "        mu_kp1 = mu_k\n"
00461     << L << "        *** Continue with next step in watchdog\n"
00462     << L << "        watch_k = 1\n"
00463     << L << "    else\n"
00464     << L << "        *** This is a good step so take it!\n"
00465     << L << "    end\n"
00466     << L << "else if watch_k == 1 then\n"
00467     << L << "    *** First watchdog iteration\n"
00468     << L << "    Do line search for: x_kp1 = x_k + alpha_k + d_k\n"
00469     << L << "        -> alpha_k, x_kp1, f_kp1, c_kp1, phi_kp1\n"
00470     << L << "    if ( phi_k <= phio ) or ( phi_kp1 <= phio + eta * Dphio ) then\n"
00471     << L << "        *** We will accept this step and reinitialize the watchdog\n"
00472     << L << "        watch_k = 0\n"
00473     << L << "    else if ( phi_kp1 > phio ) then\n"
00474     << L << "        *** This reduction is no good!\n"
00475     << L << "        *** Go back from x_k to x_km1 for this iteration (k)\n"
00476     << L << "        alpha_k        = -1.0\n"
00477     << L << "        d_k            = do\n"
00478     << L << "        f_kp1          = fo\n"
00479     << L << "        Output this iteration (k)\n"
00480     << L << "        k = k+1\n"
00481     << L << "        *** Go from x_k = x_km2 to x_kp1 for this iteration (k+1)\n"
00482     << L << "        alpha_k        = 1\n"
00483     << L << "        x_k            = xo\n"
00484     << L << "        d_k            = do\n"
00485     << L << "        Dphi_k         = Dphio\n"
00486     << L << "        phi_k          = phio\n"
00487     << L << "        x_kp1          = x_k + d_k\n"
00488     << L << "        phi_kp1        = phiop1\n"
00489     << L << "        Do line search for: x_kp1 = x_k + alpha_k + d_k\n"
00490     << L << "            -> alpha_k, x_kp1, f_kp1, c_kp1, phi_kp1\n"
00491     << L << "        Output this iteration (k+1)\n"
00492     << L << "        k = k+1\n"
00493     << L << "        *** Any updates for k (k+2) should use the last updated value\n"
00494     << L << "        *** which was for k-2 (k) since there is not much info for k-1 (k+1).\n"
00495     << L << "        *** Be careful here and make sure this works with other steps.\n"
00496     << L << "        goto EvalNewPoint\n"
00497     << L << "    else\n"
00498     << L << "        *** Accept this reduction but do a linesearch next iteration!\n"
00499     << L << "        *** Skip the update of the penalty parameter next iteration.\n"
00500     << L << "        mu_kp1 = mu_k\n"
00501     << L << "        *** Continue with next step in watchdog\n"
00502     << L << "        watch_k = 2\n"
00503     << L << "    end\n"
00504     << L << "else if ( watch_k == 2 ) then\n"
00505     << L << "    *** Second watchdog iteration\n"
00506     << L << "    Do line search for: x_kp1 = x_k + alpha_k + d_k\n"
00507     << L << "        -> alpha_k, x_kp1, f_kp1, c_kp1, phi_kp1\n"
00508     << L << "    *** Reset the watchdog algorithm\n"
00509     << L << "    watch_k = 1\n"
00510     << L << "else if ( watch_k == NORMAL_LINE_SEARCH ) then\n"
00511     << L << "    Do line search for: x_kp1 = x_k + alpha_k + d_k\n"
00512     << L << "        -> alpha_k, x_kp1, f_kp1, c_kp1, phi_kp1\n"
00513     << L << "    begin direct line search : \""
00514         << typeName(direct_line_search()) << "\"\n";
00515 
00516   direct_line_search().print_algorithm( out, L + "    " );
00517 
00518   out
00519     << L << "    end direct line search\n"
00520     << L << "end\n"
00521     << L << "if maximum number of linesearch iterations are exceeded then\n"
00522     << L << "    throw line_search_failure\n"
00523     << L << "end\n";
00524 }

Generated on Tue Oct 20 12:51:47 2009 for MOOCHO (Single Doxygen Collection) by doxygen 1.4.7