42 #ifndef MATRIX_SYM_HESSIAN_RELAX_NON_SING_H
43 #define MATRIX_SYM_HESSIAN_RELAX_NON_SING_H
45 #include "ConstrainedOptPack_Types.hpp"
46 #include "AbstractLinAlgPack_MatrixSymDiagStd.hpp"
49 namespace ConstrainedOptPack {
142 const MatrixSymOpNonsing&
G()
const;
156 MatrixOp* mwo_lhs, value_type alpha
161 ,
const Vector& v_rhs2, value_type beta)
const;
165 ,
const SpVectorSlice& sv_rhs2, value_type beta)
const;
168 VectorMutable* v_lhs, value_type alpha
171 ,
const Vector& v_rhs3, value_type beta)
const;
174 VectorMutable* v_lhs, value_type alpha
177 ,
const SpVectorSlice& sv_rhs3, value_type beta)
const;
186 MatrixSymOp* sym_lhs, value_type alpha
200 ,
const Vector& v_rhs2)
const;
204 ,
const SpVectorSlice& sv_rhs2)
const;
220 void assert_initialized()
const;
238 return M_.diag_ptr();
242 const MatrixSymOpNonsing&
245 assert_initialized();
250 const MatrixSymDiagStd&
253 assert_initialized();
259 #endif // MATRIX_SYM_HESSIAN_RELAX_NON_SING_H
void V_InvMtV(VectorMutable *v_lhs, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2) const
Teuchos::RCP< VectorMutable > vec_mut_ptr_t
const vec_mut_ptr_t & M_diag_ptr() const
Teuchos::RCP< const VectorSpace > space_ptr_t
void initialize(const G_ptr_t &G_ptr, const vec_mut_ptr_t &M_diag_ptr, const space_ptr_t &space=Teuchos::null)
Initialize the Hessian and the relaxation terms.
void Mp_StPtMtP(MatrixSymOp *sym_lhs, value_type alpha, EMatRhsPlaceHolder dummy_place_holder, const GenPermMatrixSlice &gpms_rhs, BLAS_Cpp::Transp gpms_rhs_trans, value_type beta) const
Matrix class for non-singular Hessian matrix augmented with a terms for "Big M" relaxation variables...
const VectorSpace & space_cols() const
const G_ptr_t & G_ptr() const
void Vp_StMtV(VectorMutable *v_lhs, value_type alpha, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2, value_type beta) const
const AbstractLinAlgPack::MatrixSymDiagStd & M() const
Teuchos::RCP< const MatrixSymOpNonsing > G_ptr_t
void set_uninitialized()
Set uninitalized.
void Vp_StPtMtV(VectorMutable *v_lhs, value_type alpha, const GenPermMatrixSlice &P_rhs1, BLAS_Cpp::Transp P_rhs1_trans, BLAS_Cpp::Transp M_rhs2_trans, const Vector &v_rhs3, value_type beta) const
bool Mp_StM(MatrixOp *mwo_lhs, value_type alpha, BLAS_Cpp::Transp trans_rhs) const
const MatrixSymOpNonsing & G() const
MatrixSymHessianRelaxNonSing()
Construct uninitialized.