MOOCHO (Single Doxygen Collection) Version of the Day
ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp
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 #ifndef MATRIX_SYM_POS_DEF_LBFGS_H
00030 #define MATRIX_SYM_POS_DEF_LBFGS_H
00031 
00032 #include <vector>
00033 
00034 #include "ConstrainedOptPack_Types.hpp"
00035 #include "AbstractLinAlgPack_MatrixSymSecant.hpp"
00036 #include "AbstractLinAlgPack_MatrixSymOpNonsing.hpp"
00037 #include "AbstractLinAlgPack_MultiVectorMutable.hpp"
00038 #include "AbstractLinAlgPack_VectorSpace.hpp"
00039 #include "DenseLinAlgPack_DMatrixAsTriSym.hpp"
00040 #include "Teuchos_StandardMemberCompositionMacros.hpp"
00041 
00042 namespace ConstrainedOptPack {
00043 
00097 class MatrixSymPosDefLBFGS
00098   : public AbstractLinAlgPack::MatrixSymOpNonsing // doxygen needs full path
00099   , public MatrixSymSecant
00100 {
00101 public:
00102 
00105 
00107   typedef Teuchos::RCP<const MultiVector>   multi_vec_ptr_t;
00108 
00111   class PostMod {
00112   public:
00114     PostMod(
00115       size_type     m                   = 10
00116       ,bool         maintain_original   = true
00117       ,bool         maintain_inverse    = true
00118       ,bool         auto_rescaling      = false
00119       )
00120       :m_(m)
00121       ,maintain_original_(maintain_original)
00122       ,maintain_inverse_(maintain_inverse)
00123       ,auto_rescaling_(auto_rescaling)
00124     {}
00126     void initialize(MatrixSymPosDefLBFGS* p) const
00127     {
00128       p->initial_setup(m_,maintain_original_,maintain_inverse_,auto_rescaling_);
00129     }
00130   private:
00131     size_type   m_;
00132     bool        maintain_original_;
00133     bool        maintain_inverse_;
00134     bool        auto_rescaling_;
00135   }; // end PostMod
00136 
00138 
00141 
00143   MatrixSymPosDefLBFGS(
00144       size_type     m                   = 10
00145     ,bool         maintain_original   = true
00146     ,bool         maintain_inverse    = true
00147     ,bool         auto_rescaling      = false
00148     );
00149 
00155   STANDARD_MEMBER_COMPOSITION_MEMBERS( bool, auto_rescaling );
00156 
00190    void initial_setup(
00191       size_type     m                   = 10
00192     ,bool         maintain_original   = true
00193     ,bool         maintain_inverse    = true
00194     ,bool         auto_rescaling      = false
00195      );
00196 
00197   // //////////////////////////////////
00198   // Representation access
00199 
00201   size_type m() const;
00203   size_type m_bar() const;
00205   value_type gamma_k() const;
00207   const multi_vec_ptr_t S() const;
00209   const multi_vec_ptr_t Y() const;
00211   bool maintain_original() const;
00213   bool maintain_inverse() const;
00217   size_type num_secant_updates() const;
00218 
00220 
00223 
00225   const VectorSpace& space_cols() const;
00227   std::ostream& output(std::ostream& out) const;
00229   MatrixOp& operator=(const MatrixOp& mwo);
00231   void Vp_StMtV(
00232     VectorMutable* v_lhs, value_type alpha, BLAS_Cpp::Transp trans_rhs1
00233     , const Vector& v_rhs2, value_type beta ) const;
00234 
00236 
00239 
00241   void V_InvMtV(
00242     VectorMutable* v_lhs, BLAS_Cpp::Transp trans_rhs1
00243     , const Vector& v_rhs2 ) const;
00244 
00246 
00249 
00251   void init_identity( const VectorSpace& space_diag, value_type alpha );
00258   void init_diagonal( const Vector& diag );
00260   void secant_update(
00261     VectorMutable     *s
00262     ,VectorMutable    *y
00263     ,VectorMutable    *Bs
00264     );
00265 
00267 
00268 private:
00269 
00270   // //////////////////////////////////
00271   // Private types
00272 
00273   typedef VectorSpace::multi_vec_mut_ptr_t   multi_vec_mut_ptr_t;
00274 
00275   // //////////////////////////////////
00276   // Private data members
00277 
00278   bool        maintain_original_;  // If true, qualities needed for Bk will be maintained
00279   bool        original_is_updated_;// If true, qualities needed for Bk are already updated
00280   bool        maintain_inverse_;   // If true, quantities needed for Hk will be maintained
00281   bool        inverse_is_updated_; // If true, quantities needed for Hk are already updated
00282 
00283   VectorSpace::space_ptr_t
00284             vec_spc_; // The vector space that everything is based on.
00285 
00286   size_type   n_,   // Size of the matrix.  If 0 then is uninitialized
00287               m_,   // Maximum number of update vectors that can be stored.
00288                 m_bar_, // Current number of update vectors being stored.
00289                     // 0 <= m_bar <= m
00290                 num_secant_updates_; // Records the number of secant updates performed
00291   value_type  gamma_k_;// Scaling factor for Bo = (1/gamma_k) * I.
00292 
00293   multi_vec_mut_ptr_t
00294                 S_,     // (n x m) Matrix of stored update vectors = [ s1, ..., sm ]
00295                     // S(:,m_bar) is the most recently stored s update vector
00296               Y_;     // (n_max x m) Matrix of stored update vectors = [ y1, ..., ym ]
00297             // Y(:,k_bar) is the most recently stored y update vector
00298   DMatrix   STY_,   // (m x m) The matrix S'Y
00299               STSYTY_;// ((m+1) x (m+1)) The strictly upper triangular part stores the
00300                       // upper triangular part Y'Y and the strictly lower triangular
00301                       // part stores the lower triangular part of S'S.  The diagonal
00302                       // can be used for workspace.
00303 
00304   mutable bool        Q_updated_; // True if Q has been updated for the most current update.
00305   mutable DMatrix   QJ_;        // Used to store factorization of the schur complement of Q.
00306 
00307   // //////////////////////////////////
00308   // Private member functions
00309 
00310   // Access to important matrices.
00311 
00313   const DMatrixSliceTri R() const;
00315   const DMatrixSliceTri Lb() const;
00317   DMatrixSlice STY();
00319   const DMatrixSlice STY() const;
00321   DMatrixSliceSym STS();
00323   const DMatrixSliceSym STS() const;
00325   DMatrixSliceSym YTY();
00327   const DMatrixSliceSym YTY() const;
00329   void V_invQtV( DVectorSlice* y, const DVectorSlice& x ) const;
00331   void Vp_DtV( DVectorSlice* y, const DVectorSlice& x ) const;
00332 
00333   // Updates
00334 
00336   void update_Q() const;
00337 
00339   void assert_initialized() const;
00340 
00341 };  // end class MatrixSymPosDefLBFGS
00342 
00343 // //////////////////////////////////////////////
00344 // Inline member functions
00345 
00346 inline
00347 size_type MatrixSymPosDefLBFGS::m() const
00348 {
00349   return m_;
00350 }
00351 
00352 inline
00353 size_type MatrixSymPosDefLBFGS::m_bar() const
00354 {
00355   return m_bar_;
00356 }
00357 
00358 inline
00359 value_type MatrixSymPosDefLBFGS::gamma_k() const
00360 {
00361   return gamma_k_;
00362 }
00363 
00364 inline
00365 const MatrixSymPosDefLBFGS::multi_vec_ptr_t
00366 MatrixSymPosDefLBFGS::S() const
00367 {
00368   return S_->mv_sub_view(1,n_,1,m_bar_);
00369 }
00370 
00371 inline
00372 const MatrixSymPosDefLBFGS::multi_vec_ptr_t
00373 MatrixSymPosDefLBFGS::Y() const
00374 {
00375   return Y_->mv_sub_view(1,n_,1,m_bar_);
00376 }
00377 
00378 inline
00379 bool MatrixSymPosDefLBFGS::maintain_original() const
00380 {
00381   return maintain_original_;
00382 }
00383 
00384 inline
00385 bool MatrixSymPosDefLBFGS::maintain_inverse() const
00386 {
00387   return maintain_inverse_;
00388 }
00389 
00390 inline
00391 size_type MatrixSymPosDefLBFGS::num_secant_updates() const
00392 {
00393   return num_secant_updates_;
00394 }
00395 
00396 } // end namespace ConstrainedOptPack 
00397 
00398 #endif  // MATRIX_SYM_POS_DEF_LBFGS_H
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines