42 #include "ConstrainedOptPack_MatrixDecompRangeOrthog.hpp"
43 #include "AbstractLinAlgPack_VectorSpace.hpp"
44 #include "AbstractLinAlgPack_VectorStdOps.hpp"
45 #include "AbstractLinAlgPack_MatrixSymOpNonsing.hpp"
46 #include "AbstractLinAlgPack_MatrixOpOut.hpp"
47 #include "AbstractLinAlgPack_AssertOp.hpp"
48 #include "AbstractLinAlgPack_LinAlgOpPack.hpp"
49 #include "Teuchos_Assert.hpp"
50 #include "Teuchos_FancyOStream.hpp"
52 namespace ConstrainedOptPack {
74 const char func_name[] =
"MatrixDecompRangeOrthog::initialize(...)";
76 C_ptr.get() == NULL, std::invalid_argument
77 ,func_name <<
" : Error!" );
79 D_ptr.get() == NULL, std::invalid_argument
80 ,func_name <<
" : Error!" );
82 S_ptr.get() == NULL, std::invalid_argument
83 ,func_name <<
" : Error!" );
84 #ifdef ABSTRACT_LIN_ALG_PACK_CHECK_VEC_SPCS
85 bool is_compatible = C_ptr->space_rows().is_compatible(D_ptr->space_cols());
87 !is_compatible, VectorSpace::IncompatibleVectorSpaces
88 ,func_name <<
" : Error, C_ptr->space_rows().is_compatible(D_ptr->space_cols()) == false!" );
89 is_compatible = S_ptr->space_cols().is_compatible(D_ptr->space_rows());
91 !is_compatible, VectorSpace::IncompatibleVectorSpaces
92 ,func_name <<
" : Error, S_ptr->space_cols().is_compatible(D_ptr->space_rows()) == false!" );
101 namespace rcp = MemMngPack;
102 C_ptr_ = Teuchos::null;
103 D_ptr_ = Teuchos::null;
104 S_ptr_ = Teuchos::null;
111 return C_ptr_.get() ? C_ptr_->rows() : 0;
116 return C_ptr_.get() ? C_ptr_->cols() : 0;
121 return C_ptr_->space_cols();
126 return C_ptr_->space_rows();
134 <<
"This is a " << this->
rows() <<
" x " << this->
cols()
135 <<
" nonsingular matrix (I + D'*D)*C with inverse inv(C')*(I + D*inv(S)*D') where C, D and S are:\n";
137 *out <<
"C =\n" << *
C_ptr();
138 *out <<
"D =\n" << *
D_ptr();
139 *out <<
"S =\n" << *
S_ptr();
145 ,
const Vector& x, value_type b
153 using LinAlgOpPack::V_MtV;
154 using LinAlgOpPack::Vp_MtV;
155 using LinAlgOpPack::V_StMtV;
157 assert_initialized(
"MatrixDecompRangeOrthog::Vp_StMtV(...)");
160 const MatrixOpNonsing &C = *C_ptr_;
161 const MatrixOp &D = *D_ptr_;
162 const MatrixSymOpNonsing &S = *S_ptr_;
166 VectorSpace::vec_mut_ptr_t
168 tD = D.space_cols().create_member();
179 V_MtV( tI.get(), D,
trans, x );
181 Vp_MtV( tD.get(), D,
no_trans, *tI );
182 Vp_StMtV( y, a, C, no_trans, *tD, b );
196 V_StMtV( tD.get(), a, C,
trans, x );
197 V_MtV( tI.get(), D,
trans, *tD );
216 using LinAlgOpPack::V_MtV;
217 using LinAlgOpPack::V_StMtV;
219 assert_initialized(
"MatrixDecompRangeOrthog::V_InvMtV(...)");
222 const MatrixOpNonsing &C = *C_ptr_;
223 const MatrixOp &D = *D_ptr_;
224 const MatrixSymOpNonsing &S = *S_ptr_;
228 VectorSpace::vec_mut_ptr_t
230 tIb = D.space_rows().create_member(),
231 tD = D.space_cols().create_member();
245 V_MtV( tIa.get(), D,
trans, *y );
247 Vp_StMtV( y, -1.0, D, no_trans, *tIb );
261 V_MtV( tIa.get(), D,
trans, x );
271 void MatrixDecompRangeOrthog::assert_initialized(
const char func_name[])
const
274 C_ptr_.get() == NULL, std::logic_error
275 ,func_name <<
" : Error, Must call initialize(...) first!" );
virtual const VectorSpace & space_rows() const =0
std::ostream & output(std::ostream &out) const
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
void Vp_StV(VectorMutable *v_lhs, const value_type &alpha, const Vector &v_rhs)
basic_OSTab< CharT, Traits > & incrTab(const int tabs=1)
void V_InvMtV(VectorMutable *v_lhs, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2) const
const D_ptr_t & D_ptr() const
void initialize(const C_ptr_t &C_ptr, const D_ptr_t &D_ptr, const S_ptr_t &S_ptr)
Initialize the matrix object.
const S_ptr_t & S_ptr() const
const VectorSpace & space_rows() const
TEUCHOS_DEPRECATED RCP< T > rcp(T *p, Dealloc_T dealloc, bool owns_mem)
void Vp_StMtV(VectorMutable *v_lhs, value_type alpha, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2, value_type beta) const
void V_InvMtV(VectorMutable *v_lhs, const MatrixNonsing &M_rhs1, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2)
void Vp_StMtV(VectorMutable *v_lhs, value_type alpha, const MatrixOp &M_rhs1, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2, value_type beta=1.0)
MatrixDecompRangeOrthog()
Constructs uninitialized.
Transp trans_not(Transp _trans)
void set_uninitialized()
Make uninitialized.
void Vp_MtV_assert_compatibility(VectorMutable *v_lhs, const MatrixOp &m_rhs1, BLAS_Cpp::Transp trans_rhs1, const Vector &v_rhs2)
virtual vec_mut_ptr_t create_member() const =0
const VectorSpace & space_cols() const
const C_ptr_t & C_ptr() const
void Vp_V(VectorMutable *v_lhs, const V &V_rhs)