MOOCHO (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Public Types | Public Member Functions | Private Member Functions | Private Attributes | List of all members
ConstrainedOptPack::MatrixSymPosDefBandedChol Class Reference

Matrix subclass for banded symmetric positive definite matrices and their Cholesky factors. More...

#include <ConstrainedOptPack_MatrixSymPosDefBandedChol.hpp>

Inherits MatrixSymWithOpFactorized.

Public Types

typedef Teuchos::RCP
< MemMngPack::ReleaseResource
release_resource_ptr_t
 

Public Member Functions

 MatrixSymPosDefBandedChol (size_type n=0, size_type kd=0, DMatrixSlice *MB=NULL, const release_resource_ptr_t &MB_release_resource_ptr=NULL, BLAS_Cpp::Uplo MB_uplo=BLAS_Cpp::lower, DMatrixSlice *UB=NULL, const release_resource_ptr_t &UB_release_resource_ptr=NULL, BLAS_Cpp::Uplo UB_uplo=BLAS_Cpp::lower, bool update_factor=false, value_type scale=1.0)
 Construct and Initialize. More...
 
void initialize (size_type n=0, size_type kd=0, DMatrixSlice *MB=NULL, const release_resource_ptr_t &MB_release_resource_ptr=NULL, BLAS_Cpp::Uplo MB_uplo=BLAS_Cpp::lower, DMatrixSlice *UB=NULL, const release_resource_ptr_t &UB_release_resource_ptr=NULL, BLAS_Cpp::Uplo UB_uplo=BLAS_Cpp::lower, bool update_factor=false, value_type scale=1.0)
 Initialize. More...
 
size_type kd () const
 
DMatrixSlice & MB ()
 Get view of MB. More...
 
const DMatrixSlice & MB () const
 
BLAS_Cpp::Uplo MB_uplo () const
 
DMatrixSlice & UB ()
 Get view of UB. More...
 
const DMatrixSlice & UB () const
 
BLAS_Cpp::Uplo UB_uplo () const
 
size_type rows () const
 
size_type nz () const
 
std::ostream & output (std::ostream &out) const
 
void Vp_StMtV (DVectorSlice *vs_lhs, value_type alpha, BLAS_Cpp::Transp trans_rhs1, const DVectorSlice &vs_rhs2, value_type beta) const
 
void Vp_StMtV (DVectorSlice *vs_lhs, value_type alpha, BLAS_Cpp::Transp trans_rhs1, const SpVectorSlice &sv_rhs2, value_type beta) const
 
void Vp_StPtMtV (DVectorSlice *vs_lhs, value_type alpha, const GenPermMatrixSlice &P_rhs1, BLAS_Cpp::Transp P_rhs1_trans, BLAS_Cpp::Transp M_rhs2_trans, const DVectorSlice &vs_rhs3, value_type beta) const
 
void Vp_StPtMtV (DVectorSlice *vs_lhs, value_type alpha, const GenPermMatrixSlice &P_rhs1, BLAS_Cpp::Transp P_rhs1_trans, BLAS_Cpp::Transp M_rhs2_trans, const SpVectorSlice &sv_rhs3, value_type beta) const
 
void V_InvMtV (DVectorSlice *vs_lhs, BLAS_Cpp::Transp trans_rhs1, const DVectorSlice &vs_rhs2) const
 With throw exception if factorization is not allowed. More...
 

Private Member Functions

void assert_initialized () const
 
void update_factorization () const
 

Private Attributes

size_type n_
 
size_type kd_
 
DMatrixSlice MB_
 
release_resource_ptr_t MB_release_resource_ptr_
 
BLAS_Cpp::Uplo MB_uplo_
 
DMatrixSlice UB_
 
release_resource_ptr_t UB_release_resource_ptr_
 
BLAS_Cpp::Uplo UB_uplo_
 
bool factor_updated_
 
value_type scale_
 

Detailed Description

Matrix subclass for banded symmetric positive definite matrices and their Cholesky factors.

This class is designed to support the LAPACK routines for banded symmetric positive definite matrices. The banded matrix and/or its cholesky factor are stored in simple flat rectangular matrices compatible with the LAPACK routines.

For example, for n = 8, kd = 3# the original matrix M# (if set) is stored in the following format MB#: {verbatim}

    M                                   MB

[ x x x x ] [ x x x x x x x x ] [ x x x x x ] lower triangle [ x x x x x x x o ] [ x x x x x x ] => [ x x x x x x o o ] [ x x x x x x x ] [ x x x x x o o o ] [ x x x x x x x ] [ x x x x x x ] [ o o o x x x x x ] [ x x x x x ] upper triangle [ o o x x x x x x ] [ x x x x ] => [ o x x x x x x x ] [ x x x x x x x x ] {verbatim} The Cholesky factor #U# is sorted in a similar format UB#. Technically, the matrix is M = scale * U'*U# so that M# may be negative definite as well.

Definition at line 78 of file ConstrainedOptPack_MatrixSymPosDefBandedChol.hpp.

Member Typedef Documentation

Constructor & Destructor Documentation

ConstrainedOptPack::MatrixSymPosDefBandedChol::MatrixSymPosDefBandedChol ( size_type  n = 0,
size_type  kd = 0,
DMatrixSlice *  MB = NULL,
const release_resource_ptr_t MB_release_resource_ptr = NULL,
BLAS_Cpp::Uplo  MB_uplo = BLAS_Cpp::lower,
DMatrixSlice *  UB = NULL,
const release_resource_ptr_t UB_release_resource_ptr = NULL,
BLAS_Cpp::Uplo  UB_uplo = BLAS_Cpp::lower,
bool  update_factor = false,
value_type  scale = 1.0 
)

Construct and Initialize.

This constructor just calls #this->initialize(...)#.

Member Function Documentation

void ConstrainedOptPack::MatrixSymPosDefBandedChol::initialize ( size_type  n = 0,
size_type  kd = 0,
DMatrixSlice *  MB = NULL,
const release_resource_ptr_t MB_release_resource_ptr = NULL,
BLAS_Cpp::Uplo  MB_uplo = BLAS_Cpp::lower,
DMatrixSlice *  UB = NULL,
const release_resource_ptr_t UB_release_resource_ptr = NULL,
BLAS_Cpp::Uplo  UB_uplo = BLAS_Cpp::lower,
bool  update_factor = false,
value_type  scale = 1.0 
)

Initialize.

If called with all of the default arguments then #this# will become uninitialized.

ToDo: Finish pre and post conditions!

Parameters
n[in] Determines the size of the banded matrix (n x n). If n == 0# then all of the following arguments should be left at their defaults and #this# will become uninitialized.
kd[in] Determines the band width of the matrix as defined by xPBTRF(...).
MB[in/state] If MB != NULL# then this matrix (size (kd+1) x n) is used to store the original banded matrix MB# in the format of xPBTRF(...). This matrix must be initialized on input.
MB_release_resource_ptr[in] Only significant if MB != NULL#. Points to a resource to be released when MB# is no longer needed.
MB_uplo[in] Determines if MB# is stores the upper or lower triangular elements.
UB[in/state] If UB != NULL# then this matrix (size (kd+1) x n) is used to store the Cholesky factor of the banded matrix in the format of xPBTRF(...). This matrix may or may not be initialized on input. If #update_factor == false# this this matrix must already be initialized. If #update_factor == true# then this matrix will be computed. 2 * If UB == NULL# then storage for the Cholesky factor will be computed on the fly and will be factored.
UB_release_resource_ptr[in] Only significant if UB != NULL#. Points to a resource to be released when UB# is no longer needed.
UB_uplo[in] Determines if UB# is stores the upper or lower triangular elements.
update_factor[in] If true then the factor will be computed within this function call.
scale[in] Only significant if MB != NULL# or UB != NULL# (see intro).
size_type ConstrainedOptPack::MatrixSymPosDefBandedChol::kd ( ) const
inline
DMatrixSlice & ConstrainedOptPack::MatrixSymPosDefBandedChol::MB ( )
inline

Get view of MB.

Definition at line 238 of file ConstrainedOptPack_MatrixSymPosDefBandedChol.hpp.

const DMatrixSlice & ConstrainedOptPack::MatrixSymPosDefBandedChol::MB ( ) const
inline
BLAS_Cpp::Uplo ConstrainedOptPack::MatrixSymPosDefBandedChol::MB_uplo ( ) const
inline
DMatrixSlice & ConstrainedOptPack::MatrixSymPosDefBandedChol::UB ( )
inline

Get view of UB.

Definition at line 256 of file ConstrainedOptPack_MatrixSymPosDefBandedChol.hpp.

const DMatrixSlice & ConstrainedOptPack::MatrixSymPosDefBandedChol::UB ( ) const
inline
BLAS_Cpp::Uplo ConstrainedOptPack::MatrixSymPosDefBandedChol::UB_uplo ( ) const
inline
size_type ConstrainedOptPack::MatrixSymPosDefBandedChol::rows ( ) const

size_type ConstrainedOptPack::MatrixSymPosDefBandedChol::nz ( ) const

std::ostream& ConstrainedOptPack::MatrixSymPosDefBandedChol::output ( std::ostream &  out) const

void ConstrainedOptPack::MatrixSymPosDefBandedChol::Vp_StMtV ( DVectorSlice *  vs_lhs,
value_type  alpha,
BLAS_Cpp::Transp  trans_rhs1,
const DVectorSlice &  vs_rhs2,
value_type  beta 
) const

void ConstrainedOptPack::MatrixSymPosDefBandedChol::Vp_StMtV ( DVectorSlice *  vs_lhs,
value_type  alpha,
BLAS_Cpp::Transp  trans_rhs1,
const SpVectorSlice &  sv_rhs2,
value_type  beta 
) const

void ConstrainedOptPack::MatrixSymPosDefBandedChol::Vp_StPtMtV ( DVectorSlice *  vs_lhs,
value_type  alpha,
const GenPermMatrixSlice &  P_rhs1,
BLAS_Cpp::Transp  P_rhs1_trans,
BLAS_Cpp::Transp  M_rhs2_trans,
const DVectorSlice &  vs_rhs3,
value_type  beta 
) const

void ConstrainedOptPack::MatrixSymPosDefBandedChol::Vp_StPtMtV ( DVectorSlice *  vs_lhs,
value_type  alpha,
const GenPermMatrixSlice &  P_rhs1,
BLAS_Cpp::Transp  P_rhs1_trans,
BLAS_Cpp::Transp  M_rhs2_trans,
const SpVectorSlice &  sv_rhs3,
value_type  beta 
) const

void ConstrainedOptPack::MatrixSymPosDefBandedChol::V_InvMtV ( DVectorSlice *  vs_lhs,
BLAS_Cpp::Transp  trans_rhs1,
const DVectorSlice &  vs_rhs2 
) const

With throw exception if factorization is not allowed.

void ConstrainedOptPack::MatrixSymPosDefBandedChol::assert_initialized ( ) const
private
void ConstrainedOptPack::MatrixSymPosDefBandedChol::update_factorization ( ) const
private

Member Data Documentation

size_type ConstrainedOptPack::MatrixSymPosDefBandedChol::n_
private
size_type ConstrainedOptPack::MatrixSymPosDefBandedChol::kd_
private
DMatrixSlice ConstrainedOptPack::MatrixSymPosDefBandedChol::MB_
private
release_resource_ptr_t ConstrainedOptPack::MatrixSymPosDefBandedChol::MB_release_resource_ptr_
private
BLAS_Cpp::Uplo ConstrainedOptPack::MatrixSymPosDefBandedChol::MB_uplo_
private
DMatrixSlice ConstrainedOptPack::MatrixSymPosDefBandedChol::UB_
mutableprivate
release_resource_ptr_t ConstrainedOptPack::MatrixSymPosDefBandedChol::UB_release_resource_ptr_
mutableprivate
BLAS_Cpp::Uplo ConstrainedOptPack::MatrixSymPosDefBandedChol::UB_uplo_
mutableprivate
bool ConstrainedOptPack::MatrixSymPosDefBandedChol::factor_updated_
mutableprivate
value_type ConstrainedOptPack::MatrixSymPosDefBandedChol::scale_
private

The documentation for this class was generated from the following file: