ConstrainedOptPack::MatrixSymPosDefLBFGS Class Reference

Implementation of limited Memory BFGS matrix for arbitrary vector spaces. More...

#include <ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp>

Inheritance diagram for ConstrainedOptPack::MatrixSymPosDefLBFGS:

[legend]
List of all members.

Public types

typedef Teuchos::RefCountPtr<
const MultiVector > 
multi_vec_ptr_t
 

Constructors and initializers

 MatrixSymPosDefLBFGS (size_type m=10, bool maintain_original=true, bool maintain_inverse=true, bool auto_rescaling=false)
 Calls this->initial_setup().
void auto_rescaling (const bool &auto_rescaling)
 Set whether automatic rescaling is used or not.
const bool & auto_rescaling () const
void initial_setup (size_type m=10, bool maintain_original=true, bool maintain_inverse=true, bool auto_rescaling=false)
 Initial setup for the matrix.
size_type m () const
 
size_type m_bar () const
 
value_type gamma_k () const
 
const multi_vec_ptr_t S () const
 
const multi_vec_ptr_t Y () const
 
bool maintain_original () const
 
bool maintain_inverse () const
 
size_type num_secant_updates () const
 Returns the total number of successful secant updates performed since this->init_identity() was called.
bool auto_rescaling_

Overridden from MatrixOp

const VectorSpace & space_cols () const
 
std::ostream & output (std::ostream &out) const
 
MatrixOp & operator= (const MatrixOp &mwo)
 
void Vp_StMtV (VectorMutable *v_lhs, value_type alpha, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2, value_type beta) const
 

Overridden from MatrixOpNonsing

void V_InvMtV (VectorMutable *v_lhs, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2) const
 

Overridden from MatrixSymSecant

void init_identity (const VectorSpace &space_diag, value_type alpha)
 
void init_diagonal (const Vector &diag)
 Actually this calls init_identity( diag.space(), diag.norm_inf() ).
void secant_update (VectorMutable *s, VectorMutable *y, VectorMutable *Bs)
 

Overridden from MatrixOp

std::ostream & output (std::ostream &out) const
 
MatrixOp & operator= (const MatrixOp &m)
 
void Vp_StMtV (DVectorSlice *vs_lhs, value_type alpha, BLAS_Cpp::Transp trans_rhs1, const DVectorSlice &vs_rhs2, value_type beta) const
 

Overridden from MatrixWithOpFactorized

void V_InvMtV (DVectorSlice *v_lhs, BLAS_Cpp::Transp trans_rhs1, const DVectorSlice &vs_rhs2) const
 

Overridden from MatrixSymSecant

void init_identity (size_type n, value_type alpha)
 
void init_diagonal (const DVectorSlice &diag)
 Actually this calls init_identity( (&diag)->size(), norm_inf(diag) ).
void secant_update (DVectorSlice *s, DVectorSlice *y, DVectorSlice *Bs)
 

Overridden from MatrixSymAddDelUpdateble

void initialize (value_type alpha, size_type max_size)
 This is fine as long as alpha > 0.0.
void initialize (const DMatrixSliceSym &A, size_type max_size, bool force_factorization, Inertia inertia, PivotTolerances pivot_tols)
 Sorry, this will throw an exception!
size_type max_size () const
 
Inertia inertia () const
 Returns (0,0,rows()).
void set_uninitialized ()
 Will set rows() == 0.
void augment_update (const DVectorSlice *t, value_type alpha, bool force_refactorization, EEigenValType add_eigen_val, PivotTolerances pivot_tols)
 Augment the matrix to add a row and column.
void delete_update (size_type jd, bool force_refactorization, EEigenValType drop_eigen_val, PivotTolerances pivot_tols)
 Should always succeed unless user gives a wrong value for drop_eigen_val.

Public Member Functions

 MatrixSymPosDefLBFGS (size_type max_size=0, size_type m=10, bool maintain_original=true, bool maintain_inverse=true, bool auto_rescaling=false)
 Calls initial_setup(,,,).
void auto_rescaling (const bool &auto_rescaling)
 Set whether automatic rescaling is used or not.
const bool & auto_rescaling () const
void initial_setup (size_type max_size=0, size_type m=10, bool maintain_original=true, bool maintain_inverse=true, bool auto_rescaling=false)
 Initial setup for the matrix.
size_type m () const
 
size_type m_bar () const
 
size_type k_bar () const
 
value_type gamma_k () const
 
const DMatrixSlice S () const
 
const DMatrixSlice Y () const
 
bool maintain_original () const
 
bool maintain_inverse () const
 
size_type num_secant_updates () const
 Returns the total number of successful secant updates performed.
size_type rows () const
 

Private Types

typedef VectorSpace::multi_vec_mut_ptr_t multi_vec_mut_ptr_t

Private Member Functions

const DMatrixSliceTri R () const
 
const DMatrixSliceTri Lb () const
 Strictly lower triangular part of L.
DMatrixSlice STY ()
 
const DMatrixSlice STY () const
 
DMatrixSliceSym STS ()
 
const DMatrixSliceSym STS () const
 
DMatrixSliceSym YTY ()
 
const DMatrixSliceSym YTY () const
 
void V_invQtV (DVectorSlice *y, const DVectorSlice &x) const
 y = inv(Q) * x
void Vp_DtV (DVectorSlice *y, const DVectorSlice &x) const
 y += D * x
void update_Q () const
 Update Q.
void assert_initialized () const
 
const DMatrixSliceTri R () const
 
const DMatrixSliceTri Lb () const
 Strictly lower triangular part of L.
DMatrixSlice STY ()
 
const DMatrixSlice STY () const
 
DMatrixSliceSym STS ()
 
const DMatrixSliceSym STS () const
 
DMatrixSliceSym YTY ()
 
const DMatrixSliceSym YTY () const
 
void V_invQtV (DVectorSlice *y, const DVectorSlice &x) const
 y = inv(Q) * x
void Vp_DtV (DVectorSlice *y, const DVectorSlice &x) const
 y += D * x
void update_Q () const
 Update Q.
void assert_initialized () const
 

Private Attributes

bool maintain_original_
bool original_is_updated_
bool maintain_inverse_
bool inverse_is_updated_
VectorSpace::space_ptr_t vec_spc_
size_type n_
size_type m_
size_type m_bar_
size_type num_secant_updates_
value_type gamma_k_
multi_vec_mut_ptr_t S_
multi_vec_mut_ptr_t Y_
DMatrix STY_
DMatrix STSYTY_
bool Q_updated_
DMatrix QJ_
size_type n_max_
size_type k_bar_
DMatrix S_
DMatrix Y_
DVector work_

Detailed Description

Implementation of limited Memory BFGS matrix for arbitrary vector spaces.

The function set_num_updates_stored() must be called first to set the maximum number of the most recent updates that can be stored. The storage requirements for this class are O( n*m + m*m ) which is O(n*m) when n >> m which is expected (where n is the dimension of the vector space and m is the maximum number of updates stored).

This implementation is based on:

Byrd, Nocedal, and Schnabel, "Representations of quasi-Newton matrices and their use in limited memory methods", Mathematical Programming, 63 (1994)

Consider BFGS updates of the form:


 ( B^{k-1}, s^{k-1}, y^{k-1} ) -> B^{k}
 
 where:
 
 B^{k} = B^{k-1} - ( (B*s)*(B*s)' / (s'*B*s) )^{k-1} + ( (y*y') / (s'*y) )^{k-1}
 
 B <: R^(n x n)
 s <: R^(n)
 y <: R^(n)
 
Now let us consider limited memory updating. For this implementation we set:

 Bo = ( 1 / gamma_k ) * I
 
 where:
               / (s^{k-1}'*y^{k-1})/(y^{k-1}'*y^{k-1})           :  if auto_rescaling() == true
     gamma_k = |
             \  1/alpha from last call to init_identity(n,alpha) :  otherwise
 
 
Now let us define the matrices S and Y that store the update vectors s^{i} and y^{i} for i = 1 ... m_bar:

 S = [ s^{1}, s^{2},...,s^{m_bar} ] <: R^(n x m)
 Y = [ y^{1}, y^{2},...,y^{m_bar} ] <: R^(n x m)
  
Here we are only storing the m_bar <= m most recent update vectors and their ordering is significant. The columns S(:,m_bar) and Y(:,m_bar) contain the most recent update vectors. This is all client needs to know in order to reconstruct the updates themselves.

This class allows matrix-vector products x = B*y and the inverse matrix-vector products x = inv(B)*y to be performed at a cost of about O(n*m_bar^2).

Definition at line 97 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.


Member Typedef Documentation

typedef Teuchos::RefCountPtr<const MultiVector> ConstrainedOptPack::MatrixSymPosDefLBFGS::multi_vec_ptr_t
 

Definition at line 107 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

typedef VectorSpace::multi_vec_mut_ptr_t ConstrainedOptPack::MatrixSymPosDefLBFGS::multi_vec_mut_ptr_t [private]
 

Definition at line 273 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.


Constructor & Destructor Documentation

ConstrainedOptPack::MatrixSymPosDefLBFGS::MatrixSymPosDefLBFGS size_type  m = 10,
bool  maintain_original = true,
bool  maintain_inverse = true,
bool  auto_rescaling = false
 

Calls this->initial_setup().

Definition at line 148 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

ConstrainedOptPack::MatrixSymPosDefLBFGS::MatrixSymPosDefLBFGS size_type  max_size = 0,
size_type  m = 10,
bool  maintain_original = true,
bool  maintain_inverse = true,
bool  auto_rescaling = false
 

Calls initial_setup(,,,).

Definition at line 140 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.


Member Function Documentation

void ConstrainedOptPack::MatrixSymPosDefLBFGS::auto_rescaling const bool &  auto_rescaling  )  [inline]
 

Set whether automatic rescaling is used or not.

This function must be called before a BFGS update is performed in order for it to take effect for that update.

Definition at line 155 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

const bool& ConstrainedOptPack::MatrixSymPosDefLBFGS::auto_rescaling  )  const [inline]
 

Definition at line 155 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::initial_setup size_type  m = 10,
bool  maintain_original = true,
bool  maintain_inverse = true,
bool  auto_rescaling = false
 

Initial setup for the matrix.

This function must be called before init_identity(n) is called in order for these setting to have affect. When this function is called all current updates are lost and the matrix becomes uninitialized.

Parameters:
m [in] Max number of recent update vectors s and y stored.
maintain_original [in] If true then quantities needed to compute x = Bk*y will be maintained, otherwise they will not be unless needed. This is to save computational costs in case matrix-vector products will never be needed. However, if a matrix vector product is needed then these quantities will be computed on the fly in order to satisfy the request.
maintain_inverse [in] If true then quantities needed to compute x = inv(Bk)*y = x = Hk*y will be maintained, otherwise they will not be unless needed. This is to save computational costs in case inverse matrix-vector products will never be needed. However, if the inverse product is ever needed then the needed quantities will be computed on the fly in order to satisfiy the request. Because it takes so little extra work to maintain the quantities needed for Hk it is recommended to always set this to true.
auto_rescaling [in] See intro.

Definition at line 158 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::m  )  const [inline]
 

Definition at line 347 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::m_bar  )  const [inline]
 

Definition at line 353 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

value_type ConstrainedOptPack::MatrixSymPosDefLBFGS::gamma_k  )  const [inline]
 

Definition at line 359 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

const DMatrixSlice ConstrainedOptPack::MatrixSymPosDefLBFGS::S  )  const [inline]
 

Definition at line 366 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

const DMatrixSlice ConstrainedOptPack::MatrixSymPosDefLBFGS::Y  )  const [inline]
 

Definition at line 373 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

bool ConstrainedOptPack::MatrixSymPosDefLBFGS::maintain_original  )  const [inline]
 

Definition at line 379 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

bool ConstrainedOptPack::MatrixSymPosDefLBFGS::maintain_inverse  )  const [inline]
 

Definition at line 385 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::num_secant_updates  )  const [inline]
 

Returns the total number of successful secant updates performed since this->init_identity() was called.

Definition at line 391 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

const VectorSpace & ConstrainedOptPack::MatrixSymPosDefLBFGS::space_cols  )  const [virtual]
 

Implements AbstractLinAlgPack::MatrixBase.

Definition at line 184 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

std::ostream & ConstrainedOptPack::MatrixSymPosDefLBFGS::output std::ostream &  out  )  const [virtual]
 

Reimplemented from AbstractLinAlgPack::MatrixOp.

Definition at line 189 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

MatrixOp & ConstrainedOptPack::MatrixSymPosDefLBFGS::operator= const MatrixOp &  mwo  )  [virtual]
 

Reimplemented from AbstractLinAlgPack::MatrixOp.

Definition at line 213 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::Vp_StMtV VectorMutable *  v_lhs,
value_type  alpha,
BLAS_Cpp::Transp  trans_rhs1,
const Vector &  v_rhs2,
value_type  beta
const [virtual]
 

Implements AbstractLinAlgPack::MatrixOp.

Definition at line 249 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::V_InvMtV VectorMutable *  v_lhs,
BLAS_Cpp::Transp  trans_rhs1,
const Vector &  v_rhs2
const [virtual]
 

Implements AbstractLinAlgPack::MatrixNonsing.

Definition at line 339 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::init_identity const VectorSpace &  space_diag,
value_type  alpha
 

Definition at line 463 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::init_diagonal const Vector &  diag  ) 
 

Actually this calls init_identity( diag.space(), diag.norm_inf() ).

This initialization is not convienent for this implementation. Besides, when we are using automatric rescaling (auto_rescaling == true) then this will really not matter much anyway.

Definition at line 495 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::secant_update VectorMutable *  s,
VectorMutable *  y,
VectorMutable *  Bs
 

Definition at line 500 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

const DMatrixSliceTri ConstrainedOptPack::MatrixSymPosDefLBFGS::R  )  const [inline, private]
 

Definition at line 98 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

const DMatrixSliceTri ConstrainedOptPack::MatrixSymPosDefLBFGS::Lb  )  const [inline, private]
 

Strictly lower triangular part of L.

Definition at line 104 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

DMatrixSlice ConstrainedOptPack::MatrixSymPosDefLBFGS::STY  )  [inline, private]
 

Definition at line 110 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

const DMatrixSlice ConstrainedOptPack::MatrixSymPosDefLBFGS::STY  )  const [inline, private]
 

Definition at line 116 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

DMatrixSliceSym ConstrainedOptPack::MatrixSymPosDefLBFGS::STS  )  [inline, private]
 

Definition at line 122 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

const DMatrixSliceSym ConstrainedOptPack::MatrixSymPosDefLBFGS::STS  )  const [inline, private]
 

Definition at line 128 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

DMatrixSliceSym ConstrainedOptPack::MatrixSymPosDefLBFGS::YTY  )  [inline, private]
 

Definition at line 134 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

const DMatrixSliceSym ConstrainedOptPack::MatrixSymPosDefLBFGS::YTY  )  const [inline, private]
 

Definition at line 140 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::V_invQtV DVectorSlice y,
const DVectorSlice x
const [private]
 

y = inv(Q) * x

Definition at line 792 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::Vp_DtV DVectorSlice y,
const DVectorSlice x
const [private]
 

y += D * x

Definition at line 662 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::update_Q  )  const [private]
 

Update Q.

Definition at line 702 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::assert_initialized  )  const [private]
 

Definition at line 898 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::auto_rescaling const bool &  auto_rescaling  )  [inline]
 

Set whether automatic rescaling is used or not.

This function must be called before a BFGS update is performed in order for it to take effect for that update.

Definition at line 130 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

const bool& ConstrainedOptPack::MatrixSymPosDefLBFGS::auto_rescaling  )  const [inline]
 

Definition at line 130 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::initial_setup size_type  max_size = 0,
size_type  m = 10,
bool  maintain_original = true,
bool  maintain_inverse = true,
bool  auto_rescaling = false
 

Initial setup for the matrix.

This function must be called before init_identity(n) is called. When this function is called all current updates are lost and the matrix becomes uninitialized.

Parameters:
max_size [in] If max_size > 0 then this is the max size the matrix is allowed to become. If max_size == 0 then this size will be determined by one of the initialization methods.
m [in] Max number of recent update vectors stored
maintain_original [in] If true then quantities needed to compute x = Bk*y will be maintained, otherwise they will not be unless needed. This is to save computational costs in case matrix vector products will never be needed. However, if a matrix vector product is needed then these quantities will be computed on the fly in order to satisfy the request.
maintain_inverse [in] If true then quantities needed to compute x = inv(Bk)*y = x = Hk*y will be maintained , otherwise they will not be unless needed. This is to save computational costs in case inverse matrix vector products will never be needed. However, if the inverse product is ever needed then the needed quantities will be computed on the fly in order to satisfiy the request. Because it takes so little extra work to maintain the quantities needed for Hk it is recommended to always set this to true.
auto_rescaling [in] See intro.

Definition at line 151 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::m  )  const
 

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::m_bar  )  const
 

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::k_bar  )  const [inline]
 

Definition at line 385 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

value_type ConstrainedOptPack::MatrixSymPosDefLBFGS::gamma_k  )  const
 

const DMatrixSlice ConstrainedOptPack::MatrixSymPosDefLBFGS::S  )  const
 

const DMatrixSlice ConstrainedOptPack::MatrixSymPosDefLBFGS::Y  )  const
 

bool ConstrainedOptPack::MatrixSymPosDefLBFGS::maintain_original  )  const
 

bool ConstrainedOptPack::MatrixSymPosDefLBFGS::maintain_inverse  )  const
 

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::num_secant_updates  )  const
 

Returns the total number of successful secant updates performed.

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::rows  )  const [virtual]
 

Reimplemented from AbstractLinAlgPack::MatrixBase.

Definition at line 178 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

std::ostream& ConstrainedOptPack::MatrixSymPosDefLBFGS::output std::ostream &  out  )  const [virtual]
 

Reimplemented from AbstractLinAlgPack::MatrixOp.

MatrixOp& ConstrainedOptPack::MatrixSymPosDefLBFGS::operator= const MatrixOp &  m  )  [virtual]
 

Reimplemented from AbstractLinAlgPack::MatrixOp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::Vp_StMtV DVectorSlice vs_lhs,
value_type  alpha,
BLAS_Cpp::Transp  trans_rhs1,
const DVectorSlice vs_rhs2,
value_type  beta
const
 

Definition at line 239 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::V_InvMtV DVectorSlice v_lhs,
BLAS_Cpp::Transp  trans_rhs1,
const DVectorSlice vs_rhs2
const
 

Definition at line 336 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::init_identity size_type  n,
value_type  alpha
 

Definition at line 458 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::init_diagonal const DVectorSlice diag  ) 
 

Actually this calls init_identity( (&diag)->size(), norm_inf(diag) ).

This initialization is not convienent for this implementation. Besides, when we are using automatric rescaling (auto_rescaling == true) then this will really not matter much anyway.

Definition at line 498 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::secant_update DVectorSlice s,
DVectorSlice y,
DVectorSlice Bs
 

Definition at line 504 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::initialize value_type  alpha,
size_type  max_size
 

This is fine as long as alpha > 0.0.

Definition at line 661 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::initialize const DMatrixSliceSym A,
size_type  max_size,
bool  force_factorization,
Inertia  inertia,
PivotTolerances  pivot_tols
 

Sorry, this will throw an exception!

Definition at line 678 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::max_size  )  const
 

Definition at line 691 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

MatrixSymAddDelUpdateable::Inertia ConstrainedOptPack::MatrixSymPosDefLBFGS::inertia  )  const
 

Returns (0,0,rows()).

Definition at line 696 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::set_uninitialized  ) 
 

Will set rows() == 0.

Definition at line 701 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::augment_update const DVectorSlice t,
value_type  alpha,
bool  force_refactorization,
EEigenValType  add_eigen_val,
PivotTolerances  pivot_tols
 

Augment the matrix to add a row and column.

This function is very limited in what it will do. It will throw exceptions if alpha <= 0.0 or t != NULL or add_eigen_val == EIGEN_VAL_NEG or this->rows() == this->max_size(). The obvious postconditions for this function will only technically be satisfied if alpha == this->gamma_k().

Definition at line 706 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::delete_update size_type  jd,
bool  force_refactorization,
EEigenValType  drop_eigen_val,
PivotTolerances  pivot_tols
 

Should always succeed unless user gives a wrong value for drop_eigen_val.

Definition at line 760 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.cpp.

const DMatrixSliceTri ConstrainedOptPack::MatrixSymPosDefLBFGS::R  )  const [private]
 

const DMatrixSliceTri ConstrainedOptPack::MatrixSymPosDefLBFGS::Lb  )  const [private]
 

Strictly lower triangular part of L.

DMatrixSlice ConstrainedOptPack::MatrixSymPosDefLBFGS::STY  )  [private]
 

const DMatrixSlice ConstrainedOptPack::MatrixSymPosDefLBFGS::STY  )  const [private]
 

DMatrixSliceSym ConstrainedOptPack::MatrixSymPosDefLBFGS::STS  )  [private]
 

const DMatrixSliceSym ConstrainedOptPack::MatrixSymPosDefLBFGS::STS  )  const [private]
 

DMatrixSliceSym ConstrainedOptPack::MatrixSymPosDefLBFGS::YTY  )  [private]
 

const DMatrixSliceSym ConstrainedOptPack::MatrixSymPosDefLBFGS::YTY  )  const [private]
 

void ConstrainedOptPack::MatrixSymPosDefLBFGS::V_invQtV DVectorSlice y,
const DVectorSlice x
const [private]
 

y = inv(Q) * x

void ConstrainedOptPack::MatrixSymPosDefLBFGS::Vp_DtV DVectorSlice y,
const DVectorSlice x
const [private]
 

y += D * x

void ConstrainedOptPack::MatrixSymPosDefLBFGS::update_Q  )  const [private]
 

Update Q.

void ConstrainedOptPack::MatrixSymPosDefLBFGS::assert_initialized  )  const [private]
 


Member Data Documentation

bool ConstrainedOptPack::MatrixSymPosDefLBFGS::auto_rescaling_ [private]
 

Definition at line 130 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

bool ConstrainedOptPack::MatrixSymPosDefLBFGS::maintain_original_ [private]
 

Definition at line 303 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

bool ConstrainedOptPack::MatrixSymPosDefLBFGS::original_is_updated_ [private]
 

Definition at line 304 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

bool ConstrainedOptPack::MatrixSymPosDefLBFGS::maintain_inverse_ [private]
 

Definition at line 305 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

bool ConstrainedOptPack::MatrixSymPosDefLBFGS::inverse_is_updated_ [private]
 

Definition at line 306 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

VectorSpace::space_ptr_t ConstrainedOptPack::MatrixSymPosDefLBFGS::vec_spc_ [private]
 

Definition at line 284 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::n_ [private]
 

Definition at line 308 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::m_ [private]
 

Definition at line 308 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::m_bar_ [private]
 

Definition at line 308 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::num_secant_updates_ [private]
 

Definition at line 308 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

value_type ConstrainedOptPack::MatrixSymPosDefLBFGS::gamma_k_ [private]
 

Definition at line 316 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

multi_vec_mut_ptr_t ConstrainedOptPack::MatrixSymPosDefLBFGS::S_ [private]
 

Definition at line 294 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

multi_vec_mut_ptr_t ConstrainedOptPack::MatrixSymPosDefLBFGS::Y_ [private]
 

Definition at line 294 of file ConstrainedOptPack_MatrixSymPosDefLBFGS.hpp.

DMatrix ConstrainedOptPack::MatrixSymPosDefLBFGS::STY_ [private]
 

Definition at line 318 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

DMatrix ConstrainedOptPack::MatrixSymPosDefLBFGS::STSYTY_ [private]
 

Definition at line 318 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

bool ConstrainedOptPack::MatrixSymPosDefLBFGS::Q_updated_ [mutable, private]
 

Definition at line 328 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

DMatrix ConstrainedOptPack::MatrixSymPosDefLBFGS::QJ_ [mutable, private]
 

Definition at line 329 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::n_max_ [private]
 

Definition at line 308 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

size_type ConstrainedOptPack::MatrixSymPosDefLBFGS::k_bar_ [private]
 

Definition at line 308 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

DMatrix ConstrainedOptPack::MatrixSymPosDefLBFGS::S_ [private]
 

Definition at line 318 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

DMatrix ConstrainedOptPack::MatrixSymPosDefLBFGS::Y_ [private]
 

Definition at line 318 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.

DVector ConstrainedOptPack::MatrixSymPosDefLBFGS::work_ [mutable, private]
 

Definition at line 331 of file ConstrainedOptPack_MatrixSymPosDefLBFGSSerial.hpp.


The documentation for this class was generated from the following files:
Generated on Thu Sep 18 12:36:15 2008 for MOOCHO (Single Doxygen Collection) by doxygen 1.3.9.1