42 #ifndef MATRIX_SYM_POS_DEF_LBFGS_H
43 #define MATRIX_SYM_POS_DEF_LBFGS_H
55 namespace ConstrainedOptPack {
112 ,
public MatrixSymSecant
131 ,
bool auto_rescaling =
false
160 ,
bool auto_rescaling =
false
207 ,
bool auto_rescaling =
false
240 std::ostream&
output(std::ostream&
out)
const;
242 MatrixOp&
operator=(
const MatrixOp& mwo);
246 ,
const Vector& v_rhs2,
value_type beta )
const;
256 ,
const Vector& v_rhs2 )
const;
296 VectorSpace::space_ptr_t
411 #endif // MATRIX_SYM_POS_DEF_LBFGS_H
void V_invQtV(DVectorSlice *y, const DVectorSlice &x) const
y = inv(Q) * x
const DMatrixSliceTri Lb() const
Strictly lower triangular part of L.
void init_identity(const VectorSpace &space_diag, value_type alpha)
AbstractLinAlgPack::size_type size_type
bool maintain_original() const
Abstract base class for all polymorphic symmetrix nonsingular matrices that can be used to compute ma...
void secant_update(VectorMutable *s, VectorMutable *y, VectorMutable *Bs)
void initial_setup(size_type m=10, bool maintain_original=true, bool maintain_inverse=true, bool auto_rescaling=false)
Initial setup for the matrix.
PostMod(size_type m=10, bool maintain_original=true, bool maintain_inverse=true, bool auto_rescaling=false)
const VectorSpace & space_cols() const
bool original_is_updated_
VectorSpace::space_ptr_t vec_spc_
void V_InvMtV(VectorMutable *v_lhs, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2) const
void update_Q() const
Update Q.
MatrixSymPosDefLBFGS(size_type m=10, bool maintain_original=true, bool maintain_inverse=true, bool auto_rescaling=false)
Calls this->initial_setup()
void init_diagonal(const Vector &diag)
Actually this calls init_identity( diag.space(), diag.norm_inf() ).
size_type num_secant_updates_
PostMod class to use with MemMngPack::AbstractFactorStd.
const multi_vec_ptr_t S() const
STANDARD_MEMBER_COMPOSITION_MEMBERS(bool, auto_rescaling)
Set whether automatic rescaling is used or not.
Teuchos::RCP< const MultiVector > multi_vec_ptr_t
const DMatrixSliceTri R() const
const multi_vec_ptr_t Y() const
VectorSpace::multi_vec_mut_ptr_t multi_vec_mut_ptr_t
std::ostream & output(std::ostream &out) const
void assert_initialized() const
DenseLinAlgPack::VectorSliceTmpl< value_type > DVectorSlice
void Vp_StMtV(VectorMutable *v_lhs, value_type alpha, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2, value_type beta) const
Implementation of limited Memory BFGS matrix for arbitrary vector spaces.
void initialize(MatrixSymPosDefLBFGS *p) const
AbstractLinAlgPack::value_type value_type
DenseLinAlgPack::DMatrixSliceTri DMatrixSliceTri
bool maintain_inverse() const
size_type num_secant_updates() const
Returns the total number of successful secant updates performed since this->init_identity() was calle...
MatrixOp & operator=(const MatrixOp &mwo)
DenseLinAlgPack::DMatrixSlice DMatrixSlice
DenseLinAlgPack::DMatrixSliceSym DMatrixSliceSym
value_type gamma_k() const
void Vp_DtV(DVectorSlice *y, const DVectorSlice &x) const
y += D * x