MoochoPack_CheckSkipBFGSUpdateStd_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 <math.h>
00030 
00031 #include <ostream>
00032 
00033 #include "MoochoPack_CheckSkipBFGSUpdateStd_Step.hpp"
00034 #include "MoochoPack_moocho_algo_conversion.hpp"
00035 #include "IterationPack_print_algorithm_step.hpp"
00036 #include "AbstractLinAlgPack_MatrixSymOp.hpp"
00037 #include "AbstractLinAlgPack_MatrixOpOut.hpp"
00038 
00039 namespace MoochoPack {
00040 
00041 CheckSkipBFGSUpdateStd_Step::CheckSkipBFGSUpdateStd_Step(
00042   value_type  skip_bfgs_prop_const
00043   )
00044   :skip_bfgs_prop_const_(skip_bfgs_prop_const)
00045 {}
00046 
00047 bool CheckSkipBFGSUpdateStd_Step::do_step(
00048   Algorithm& _algo, poss_type step_poss, IterationPack::EDoStepType type
00049   , poss_type assoc_step_poss
00050   )
00051 {
00052   NLPAlgo &algo = rsqp_algo(_algo);
00053   NLPAlgoState  &s    = algo.rsqp_state();
00054 
00055   EJournalOutputLevel olevel = algo.algo_cntr().journal_output_level();
00056   std::ostream& out = algo.track().journal_out();
00057 
00058   // print step header.
00059   if( static_cast<int>(olevel) >= static_cast<int>(PRINT_ALGORITHM_STEPS) ) {
00060     using IterationPack::print_algorithm_step;
00061     print_algorithm_step( algo, step_poss, type, assoc_step_poss, out );
00062   }
00063 
00064   IterQuantityAccess<MatrixSymOp>
00065     &rHL = s.rHL();
00066   if( rHL.updated_k(-1) ) {
00067     bool skip_update = true;
00068     if( !s.Ypy().updated_k(-1) || !s.Zpz().updated_k(-1)
00069       || !s.rGL().updated_k(-1) || !s.c().updated_k(-1) )
00070     {
00071       if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00072         out
00073           << "\n*** Warning, rHL_km1 is updated but one or more of the quantities Ypy_km1, Zpz_km1\n"
00074             ",rGL_km1 and c_km1 are not updated.  Therefore, there is not sufficient information\n"
00075             "to determine if to skip the BFGS update or not.  Check storage requirements for\n"
00076             " the above quantities\n"
00077             "The BFGS wupdate will be skipped in this case ...\n";
00078       } 
00079       skip_update = true;
00080     }
00081     else {
00082       // The information exists for the update so determine
00083       // if we are in the region to perform the BFGS update.
00084       //
00085       // Check if we are to skip the update for this iteration
00086       const value_type
00087         nrm_rGL_km1 = s.rGL().get_k(-1).norm_2(),
00088         nrm_c_km1 = s.c().get_k(-1).norm_2(),
00089         nrm_Zpz_km1 = s.Zpz().get_k(-1).norm_2(),
00090         nrm_Ypy_km1 = s.Ypy().get_k(-1).norm_2();
00091       // ratio = (skip_bfgs_prop_const / sqrt(||rGL_km1|| + ||c_km1||)) * ( ||Zpz_km1|| / ||Ypy_km1|| )
00092       value_type
00093         ratio = ( skip_bfgs_prop_const() / ::sqrt( nrm_rGL_km1 + nrm_c_km1 ) )
00094           * ( nrm_Zpz_km1 / nrm_Ypy_km1 );
00095       // If ratio < 1.0 then skip the update
00096       skip_update = ratio < 1.0;
00097       if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00098         out
00099           << std::endl
00100           << "ratio = (skip_bfgs_prop_const/sqrt(||rGL_km1||+||c_km1||))*(||Zpz_km1||/||Ypy_km1||)\n"
00101           << "      = (" << skip_bfgs_prop_const() << "/sqrt("<<nrm_rGL_km1<<"+"<<nrm_c_km1<<"))\n"
00102           << "        * ("<<nrm_Zpz_km1<<"/"<<nrm_Ypy_km1<<")\n"
00103           << "      = " << ratio << std::endl
00104           << "ratio " << (skip_update ? '<' : '>' ) << " 1\n"
00105           << (skip_update
00106             ? "Skipping BFGS update ...\n"
00107             : "Perform BFGS update if you can ...\n"  );
00108       }
00109     }
00110     if( skip_update ) {
00111       if( (int)olevel >= (int)PRINT_ALGORITHM_STEPS ) {
00112         out
00113           << "\nrHL_k = rHL_km1\n";
00114       }
00115       const MatrixSymOp &rHL_km1 = rHL.get_k(-1);
00116       rHL.set_k(0) = rHL_km1;
00117       quasi_newton_stats_(s).set_k(0).set_updated_stats(
00118         QuasiNewtonStats::SKIPED );
00119       if( (int)olevel >= (int)PRINT_ITERATION_QUANTITIES ) {
00120         out << "\nrHL_k =\n" << rHL.get_k(0);
00121       }
00122     }
00123   }
00124   
00125   return true;
00126 }
00127 
00128 void CheckSkipBFGSUpdateStd_Step::print_step( const Algorithm& algo
00129   , poss_type step_poss, IterationPack::EDoStepType type, poss_type assoc_step_poss
00130   , std::ostream& out, const std::string& L ) const
00131 {
00132   out
00133     << L << "*** Check if we should do the BFGS update\n"
00134     << L << "if rHL_km1 is update then\n"
00135     << L << "    If Ypy_km1, Zpz_km1, rGL_km1, or c_km1 is not updated then\n"
00136     << L << "        *** Warning, insufficient information to determine if we should\n"
00137     << L << "        *** perform the update.  Check for sufficient backward storage.\n"
00138     << L << "        rHL_k = rHL_km1\n"
00139     << L << "    else\n"
00140     << L << "        *** Check if we are in the proper region\n"
00141     << L << "        ratio = (skip_bfgs_prop_const/sqrt(norm(rGL_km1,2)+norm(c_km1,2)))\n"
00142     << L << "                 * (norm(Zpz_km1,2)/norm(Ypy_km1,2) )\n"
00143     << L << "        if ratio < 1 then \n"
00144     << L << "            rHL_k = rHL_km1\n"
00145     << L << "        end\n"
00146     << L << "    end\n"
00147     << L << "end\n";
00148 }
00149 
00150 } // end namespace MoochoPack

Generated on Wed May 12 21:52:30 2010 for MOOCHO (Single Doxygen Collection) by  doxygen 1.4.7