43 #ifndef _TEUCHOS_SERIALSYMDENSEMATRIX_HPP_ 44 #define _TEUCHOS_SERIALSYMDENSEMATRIX_HPP_ 120 template<
typename OrdinalType,
typename ScalarType>
202 int shape(OrdinalType numRowsCols);
227 int reshape(OrdinalType numRowsCols);
292 ScalarType&
operator () (OrdinalType rowIndex, OrdinalType colIndex);
302 const ScalarType&
operator () (OrdinalType rowIndex, OrdinalType colIndex)
const;
396 virtual void print(std::ostream& os)
const;
404 void scale(
const ScalarType alpha );
407 void copyMat(
bool inputUpper, ScalarType* inputMatrix, OrdinalType inputStride,
408 OrdinalType numRowCols,
bool outputUpper, ScalarType* outputMatrix,
409 OrdinalType outputStride, OrdinalType startRowCol,
413 void copyUPLOMat(
bool inputUpper, ScalarType* inputMatrix,
414 OrdinalType inputStride, OrdinalType inputRows);
417 void checkIndex( OrdinalType rowIndex, OrdinalType colIndex = 0 )
const;
432 template<
typename OrdinalType,
typename ScalarType>
435 numRowCols_(0), stride_(0), valuesCopied_(false), values_(0), upper_(false), UPLO_(
'L')
438 template<
typename OrdinalType,
typename ScalarType>
441 numRowCols_(numRowCols_in), stride_(numRowCols_in), valuesCopied_(false), upper_(false), UPLO_(
'L')
449 template<
typename OrdinalType,
typename ScalarType>
451 DataAccess CV,
bool upper_in, ScalarType* values_in, OrdinalType stride_in, OrdinalType numRowCols_in
454 numRowCols_(numRowCols_in), stride_(stride_in), valuesCopied_(false),
455 values_(values_in), upper_(upper_in)
471 template<
typename OrdinalType,
typename ScalarType>
474 numRowCols_(Source.numRowCols_), stride_(0), valuesCopied_(true),
475 values_(0), upper_(Source.upper_), UPLO_(Source.UPLO_)
488 values_ =
new ScalarType[newsize];
498 template<
typename OrdinalType,
typename ScalarType>
501 ScalarType> &Source, OrdinalType numRowCols_in, OrdinalType startRowCol )
503 numRowCols_(numRowCols_in), stride_(Source.stride_), valuesCopied_(false), upper_(Source.upper_), UPLO_(Source.UPLO_)
518 template<
typename OrdinalType,
typename ScalarType>
528 template<
typename OrdinalType,
typename ScalarType>
532 numRowCols_ = numRowCols_in;
533 stride_ = numRowCols_;
534 values_ =
new ScalarType[stride_*numRowCols_];
536 valuesCopied_ =
true;
540 template<
typename OrdinalType,
typename ScalarType>
544 numRowCols_ = numRowCols_in;
545 stride_ = numRowCols_;
546 values_ =
new ScalarType[stride_*numRowCols_];
547 valuesCopied_ =
true;
551 template<
typename OrdinalType,
typename ScalarType>
555 ScalarType* values_tmp =
new ScalarType[numRowCols_in * numRowCols_in];
557 for(OrdinalType k = 0; k < numRowCols_in * numRowCols_in; k++)
559 values_tmp[k] = zero;
561 OrdinalType numRowCols_tmp =
TEUCHOS_MIN(numRowCols_, numRowCols_in);
564 copyMat(upper_, values_, stride_, numRowCols_tmp, upper_, values_tmp, numRowCols_in, 0);
567 numRowCols_ = numRowCols_in;
568 stride_ = numRowCols_;
569 values_ = values_tmp;
570 valuesCopied_ =
true;
578 template<
typename OrdinalType,
typename ScalarType>
582 if (upper_ !=
false) {
583 copyUPLOMat(
true, values_, stride_, numRowCols_ );
589 template<
typename OrdinalType,
typename ScalarType>
593 if (upper_ ==
false) {
594 copyUPLOMat(
false, values_, stride_, numRowCols_ );
600 template<
typename OrdinalType,
typename ScalarType>
604 if (fullMatrix ==
true) {
605 for(OrdinalType j = 0; j < numRowCols_; j++)
607 for(OrdinalType i = 0; i < numRowCols_; i++)
609 values_[i + j*stride_] = value_in;
616 for(OrdinalType j = 0; j < numRowCols_; j++) {
617 for(OrdinalType i = 0; i <= j; i++) {
618 values_[i + j*stride_] = value_in;
623 for(OrdinalType j = 0; j < numRowCols_; j++) {
624 for(OrdinalType i = j; i < numRowCols_; i++) {
625 values_[i + j*stride_] = value_in;
633 template<
typename OrdinalType,
typename ScalarType>
641 for(OrdinalType j = 0; j < numRowCols_; j++) {
642 for(OrdinalType i = 0; i < j; i++) {
650 for(OrdinalType j = 0; j < numRowCols_; j++) {
651 for(OrdinalType i = j+1; i < numRowCols_; i++) {
660 for(OrdinalType i = 0; i < numRowCols_; i++) {
661 values_[i + i*stride_] = diagSum[i] + bias;
666 template<
typename OrdinalType,
typename ScalarType>
687 UPLO_ = Source.
UPLO_;
695 UPLO_ = Source.
UPLO_;
696 const OrdinalType newsize = stride_ * numRowCols_;
698 values_ =
new ScalarType[newsize];
699 valuesCopied_ =
true;
710 UPLO_ = Source.
UPLO_;
717 UPLO_ = Source.
UPLO_;
718 const OrdinalType newsize = stride_ * numRowCols_;
720 values_ =
new ScalarType[newsize];
721 valuesCopied_ =
true;
730 template<
typename OrdinalType,
typename ScalarType>
737 template<
typename OrdinalType,
typename ScalarType>
749 template<
typename OrdinalType,
typename ScalarType>
761 template<
typename OrdinalType,
typename ScalarType>
775 copyMat(Source.
upper_, Source.
values_, Source.
stride_, numRowCols_, upper_, values_, stride_, 0 );
783 template<
typename OrdinalType,
typename ScalarType>
786 #ifdef HAVE_TEUCHOS_ARRAY_BOUNDSCHECK 787 checkIndex( rowIndex, colIndex );
789 if ( rowIndex <= colIndex ) {
792 return(values_[colIndex * stride_ + rowIndex]);
794 return(values_[rowIndex * stride_ + colIndex]);
799 return(values_[rowIndex * stride_ + colIndex]);
801 return(values_[colIndex * stride_ + rowIndex]);
805 template<
typename OrdinalType,
typename ScalarType>
808 #ifdef HAVE_TEUCHOS_ARRAY_BOUNDSCHECK 809 checkIndex( rowIndex, colIndex );
811 if ( rowIndex <= colIndex ) {
814 return(values_[colIndex * stride_ + rowIndex]);
816 return(values_[rowIndex * stride_ + colIndex]);
821 return(values_[rowIndex * stride_ + colIndex]);
823 return(values_[colIndex * stride_ + rowIndex]);
831 template<
typename OrdinalType,
typename ScalarType>
837 template<
typename OrdinalType,
typename ScalarType>
848 for (j=0; j<numRowCols_; j++) {
850 ptr = values_ + j*stride_;
851 for (i=0; i<j; i++) {
854 ptr = values_ + j + j*stride_;
855 for (i=j; i<numRowCols_; i++) {
863 for (j=0; j<numRowCols_; j++) {
865 ptr = values_ + j + j*stride_;
866 for (i=j; i<numRowCols_; i++) {
870 for (i=0; i<j; i++) {
880 template<
typename OrdinalType,
typename ScalarType>
890 for (j = 0; j < numRowCols_; j++) {
891 for (i = 0; i < j; i++) {
898 for (j = 0; j < numRowCols_; j++) {
900 for (i = j+1; i < numRowCols_; i++) {
913 template<
typename OrdinalType,
typename ScalarType>
922 for(i = 0; i < numRowCols_; i++) {
923 for(j = 0; j < numRowCols_; j++) {
924 if((*
this)(i, j) != Operand(i, j)) {
933 template<
typename OrdinalType,
typename ScalarType>
936 return !((*this) == Operand);
943 template<
typename OrdinalType,
typename ScalarType>
950 for (j=0; j<numRowCols_; j++) {
951 ptr = values_ + j*stride_;
952 for (i=0; i<= j; i++) { *ptr = alpha * (*ptr); ptr++; }
956 for (j=0; j<numRowCols_; j++) {
957 ptr = values_ + j*stride_ + j;
958 for (i=j; i<numRowCols_; i++) { *ptr = alpha * (*ptr); ptr++; }
992 template<
typename OrdinalType,
typename ScalarType>
997 os <<
"Values_copied : yes" << std::endl;
999 os <<
"Values_copied : no" << std::endl;
1000 os <<
"Rows / Columns : " << numRowCols_ << std::endl;
1001 os <<
"LDA : " << stride_ << std::endl;
1003 os <<
"Storage: Upper" << std::endl;
1005 os <<
"Storage: Lower" << std::endl;
1006 if(numRowCols_ == 0) {
1007 os <<
"(matrix is empty, no values to display)" << std::endl;
1009 for(OrdinalType i = 0; i < numRowCols_; i++) {
1010 for(OrdinalType j = 0; j < numRowCols_; j++){
1011 os << (*this)(i,j) <<
" ";
1022 template<
typename OrdinalType,
typename ScalarType>
1025 "SerialSymDenseMatrix<T>::checkIndex: " 1026 "Row index " << rowIndex <<
" out of range [0, "<< numRowCols_ <<
")");
1028 "SerialSymDenseMatrix<T>::checkIndex: " 1029 "Col index " << colIndex <<
" out of range [0, "<< numRowCols_ <<
")");
1032 template<
typename OrdinalType,
typename ScalarType>
1039 valuesCopied_ =
false;
1043 template<
typename OrdinalType,
typename ScalarType>
1045 bool inputUpper, ScalarType* inputMatrix,
1046 OrdinalType inputStride, OrdinalType numRowCols_in,
1047 bool outputUpper, ScalarType* outputMatrix,
1048 OrdinalType outputStride, OrdinalType startRowCol,
1053 ScalarType* ptr1 = 0;
1054 ScalarType* ptr2 = 0;
1056 for (j = 0; j < numRowCols_in; j++) {
1057 if (inputUpper ==
true) {
1059 ptr2 = inputMatrix + (j + startRowCol) * inputStride + startRowCol;
1060 if (outputUpper ==
true) {
1062 ptr1 = outputMatrix + j*outputStride;
1064 for(i = 0; i <= j; i++) {
1065 *ptr1++ += alpha*(*ptr2++);
1068 for(i = 0; i <= j; i++) {
1076 ptr1 = outputMatrix + j;
1078 for(i = 0; i <= j; i++) {
1079 *ptr1 += alpha*(*ptr2++);
1080 ptr1 += outputStride;
1083 for(i = 0; i <= j; i++) {
1085 ptr1 += outputStride;
1092 ptr2 = inputMatrix + (startRowCol+j) * inputStride + startRowCol + j;
1093 if (outputUpper ==
true) {
1096 ptr1 = outputMatrix + j*outputStride + j;
1098 for(i = j; i < numRowCols_in; i++) {
1099 *ptr1 += alpha*(*ptr2++);
1100 ptr1 += outputStride;
1103 for(i = j; i < numRowCols_in; i++) {
1105 ptr1 += outputStride;
1111 ptr1 = outputMatrix + j*outputStride + j;
1113 for(i = j; i < numRowCols_in; i++) {
1114 *ptr1++ += alpha*(*ptr2++);
1117 for(i = j; i < numRowCols_in; i++) {
1126 template<
typename OrdinalType,
typename ScalarType>
1128 bool inputUpper, ScalarType* inputMatrix,
1129 OrdinalType inputStride, OrdinalType inputRows
1133 ScalarType * ptr1 = 0;
1134 ScalarType * ptr2 = 0;
1137 for (j=1; j<inputRows; j++) {
1138 ptr1 = inputMatrix + j;
1139 ptr2 = inputMatrix + j*inputStride;
1140 for (i=0; i<j; i++) {
1147 for (i=1; i<inputRows; i++) {
1148 ptr1 = inputMatrix + i;
1149 ptr2 = inputMatrix + i*inputStride;
1150 for (j=0; j<i; j++) {
ScalarTraits< ScalarType >::magnitudeType normInf() const
Returns the Infinity-norm of the matrix.
virtual void print(std::ostream &os) const
Print method. Defines the behavior of the std::ostream << operator inherited from the Object class...
OrdinalType numCols() const
Returns the column dimension of this matrix.
void copyUPLOMat(bool inputUpper, ScalarType *inputMatrix, OrdinalType inputStride, OrdinalType inputRows)
T magnitudeType
Mandatory typedef for result of magnitude.
SerialSymDenseMatrix()
Default constructor; defines a zero size object.
void checkIndex(OrdinalType rowIndex, OrdinalType colIndex=0) const
Templated interface class to BLAS routines.
bool empty() const
Returns whether this matrix is empty.
ScalarType & operator()(OrdinalType rowIndex, OrdinalType colIndex)
Element access method (non-const).
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Macro for throwing an exception with breakpointing to ease debugging.
int shape(OrdinalType numRowsCols)
Set dimensions of a Teuchos::SerialSymDenseMatrix object; init values to zero.
Teuchos header file which uses auto-configuration information to include necessary C++ headers...
This class creates and provides basic support for symmetric, positive-definite dense matrices of temp...
bool operator==(const SerialSymDenseMatrix< OrdinalType, ScalarType > &Operand) const
Equality of two matrices.
Object for storing data and providing functionality that is common to all computational classes...
int reshape(OrdinalType numRowsCols)
Reshape a Teuchos::SerialSymDenseMatrix object.
This structure defines some basic traits for a scalar field type.
void setLower()
Specify that the lower triangle of the this matrix should be used.
void setUpper()
Specify that the upper triangle of the this matrix should be used.
Functionality and data that is common to all computational classes.
OrdinalType numRows() const
Returns the row dimension of this matrix.
bool upper() const
Returns true if upper triangular part of this matrix has and will be used.
SerialSymDenseMatrix< OrdinalType, ScalarType > & operator*=(const ScalarType alpha)
Inplace scalar-matrix product A = alpha*A.
int putScalar(const ScalarType value=Teuchos::ScalarTraits< ScalarType >::zero(), bool fullMatrix=false)
Set all values in the matrix to a constant value.
SerialSymDenseMatrix< OrdinalType, ScalarType > & operator+=(const SerialSymDenseMatrix< OrdinalType, ScalarType > &Source)
Add another matrix to this matrix.
#define TEUCHOS_MAX(x, y)
int random(const ScalarType bias=0.1 *Teuchos::ScalarTraits< ScalarType >::one())
Set all values in the active area (upper/lower triangle) of this matrix to be random numbers...
virtual ~SerialSymDenseMatrix()
Teuchos::SerialSymDenseMatrix destructor.
void scale(const ScalarType alpha)
static magnitudeType magnitude(T a)
Returns the magnitudeType of the scalar type a.
SerialSymDenseMatrix< OrdinalType, ScalarType > & operator-=(const SerialSymDenseMatrix< OrdinalType, ScalarType > &Source)
Subtract another matrix from this matrix.
OrdinalType ordinalType
Typedef for ordinal type.
SerialSymDenseMatrix< OrdinalType, ScalarType > & operator=(const SerialSymDenseMatrix< OrdinalType, ScalarType > &Source)
Copies values from one matrix to another.
int shapeUninitialized(OrdinalType numRowsCols)
Set dimensions of a Teuchos::SerialSymDenseMatrix object; don't initialize values.
ScalarTraits< ScalarType >::magnitudeType normOne() const
Returns the 1-norm of the matrix.
Defines basic traits for the scalar field type.
static T zero()
Returns representation of zero for this scalar type.
static T random()
Returns a random number (between -one() and +one()) of this scalar type.
char UPLO() const
Returns character value of UPLO used by LAPACK routines.
#define TEUCHOS_CHK_REF(a)
SerialSymDenseMatrix< OrdinalType, ScalarType > & assign(const SerialSymDenseMatrix< OrdinalType, ScalarType > &Source)
Copies values from one matrix to another.
#define TEUCHOS_MIN(x, y)
Teuchos::DataAccess Mode enumerable type.
ScalarType * values() const
Returns the pointer to the ScalarType data array contained in the object.
OrdinalType stride() const
Returns the stride between the columns of this matrix in memory.
void copyMat(bool inputUpper, ScalarType *inputMatrix, OrdinalType inputStride, OrdinalType numRowCols, bool outputUpper, ScalarType *outputMatrix, OrdinalType outputStride, OrdinalType startRowCol, ScalarType alpha=ScalarTraits< ScalarType >::zero())
ScalarType scalarType
Typedef for scalar type.
ScalarTraits< ScalarType >::magnitudeType normFrobenius() const
Returns the Frobenius-norm of the matrix.
bool operator!=(const SerialSymDenseMatrix< OrdinalType, ScalarType > &Operand) const
Inequality of two matrices.