CVM Class Library  8.1
This C++ class library encapsulates concepts of vector and different matrices including square, band, symmetric and hermitian ones in Euclidean space of real and complex numbers.
 All Classes Files Functions Variables Typedefs Friends Macros Pages
Public Member Functions | List of all members
basic_schmatrix< TR, TC > Class Template Reference

End-user class encapsulating hermitian matrix of complex numbers. More...

#include <cvm.h>

Inheritance diagram for basic_schmatrix< TR, TC >:
Inheritance graph
[legend]
Collaboration diagram for basic_schmatrix< TR, TC >:
Collaboration graph
[legend]

Public Member Functions

 basic_schmatrix ()
 Default constructor.
 basic_schmatrix (tint nDim)
 Constructor.
 basic_schmatrix (TC *pd, tint nDim, TR tol=basic_cvmMachSp< TR >())
 Constructor.
 basic_schmatrix (const TC *pd, tint nDim, TR tol=basic_cvmMachSp< TR >())
 Constructor.
 basic_schmatrix (const basic_schmatrix &m)
 Copy constructor.
 basic_schmatrix (basic_schmatrix &&m) noexcept
 Move constructor.
 basic_schmatrix (const BaseCMatrix &m, TR tol=basic_cvmMachSp< TR >())
 Constructor.
 basic_schmatrix (const RVector &v)
 Constructor.
 basic_schmatrix (const basic_srsmatrix< TR > &m)
 Constructor.
 basic_schmatrix (const TR *pRe, const TR *pIm, tint nDim, TR tol=basic_cvmMachSp< TR >())
 Constructor.
 basic_schmatrix (const basic_srmatrix< TR > &mRe, const basic_srmatrix< TR > &mIm, TR tol=basic_cvmMachSp< TR >())
 Constructor.
 basic_schmatrix (basic_schmatrix &m, tint nRowCol, tint nDim)
 Submatrix constructor.
TC operator() (tint nRow, tint nCol) const throw (cvmexception)
 Value of element (not l-value)
const CVector operator() (tint nCol) const throw (cvmexception)
 Column as not l-value.
const CVector operator[] (tint nRow) const throw (cvmexception)
 Row as not l-value.
const CVector diag (tint nDiag) const throw (cvmexception)
 Diagonal (not l-value)
const basic_srsmatrix< TR > real () const
 Real part (not l-value)
const basic_srmatrix< TR > imag () const
 Imaginary part (not l-value)
basic_schmatrixoperator= (const basic_schmatrix &m) throw (cvmexception)
 Assignment operator.
basic_schmatrixoperator= (basic_schmatrix &&m) throw (cvmexception)
 Move assignment operator.
basic_schmatrixassign (const CVector &v, TR tol=basic_cvmMachSp< TR >()) throw (cvmexception)
 Vector (as array) assignment.
basic_schmatrixassign (const TC *pd, TR tol=basic_cvmMachSp< TR >()) throw (cvmexception)
 External array assignment.
basic_schmatrixassign (tint nRowCol, const basic_schmatrix &m) throw (cvmexception)
 Assignment to submatrix.
basic_schmatrixset (tint nRow, tint nCol, TC c) throw (cvmexception)
 Sets both elements to keep matrix hermitian.
basic_schmatrixset_diag (tint nDiag, const CVector &vDiag) throw (cvmexception)
 Sets lower and upper diagonals.
basic_schmatrixset_main_diag (const RVector &vDiag) throw (cvmexception)
 Sets main diagonal.
basic_schmatrixassign_real (const basic_srsmatrix< TR > &mRe) throw (cvmexception)
 Assigns real symmetric matrix to real part.
basic_schmatrixset_real (TR d)
 Sets all real parts to one value.
basic_schmatrixresize (tint nNewDim) throw (cvmexception)
 Changes dimension.
bool operator== (const basic_schmatrix &m) const
 Matrix comparison.
bool operator!= (const basic_schmatrix &m) const
 Matrix comparison.
basic_schmatrixoperator<< (const basic_schmatrix &m) throw (cvmexception)
 Matrix replacement.
basic_schmatrix operator+ (const basic_schmatrix &m) const throw (cvmexception)
 Addition operator.
basic_schmatrix operator- (const basic_schmatrix &m) const throw (cvmexception)
 Subtraction operator.
basic_schmatrixsum (const basic_schmatrix &m1, const basic_schmatrix &m2) throw (cvmexception)
 Sum of matrices.
basic_schmatrixdiff (const basic_schmatrix &m1, const basic_schmatrix &m2) throw (cvmexception)
 Difference of matrices.
basic_schmatrixoperator+= (const basic_schmatrix &m) throw (cvmexception)
 Increment operator.
basic_schmatrixoperator-= (const basic_schmatrix &m) throw (cvmexception)
 Decrement operator.
basic_schmatrix operator- () const
 Unary minus operator.
basic_schmatrixoperator++ ()
 Plus identity, prefix.
basic_schmatrix operator++ (int)
 Plus identity, postfix.
basic_schmatrixoperator-- ()
 Minus identity, prefix.
basic_schmatrix operator-- (int)
 Minus identity, postfix.
basic_schmatrix operator* (TR dMult) const
 Multiply by real number operator.
basic_schmatrix operator/ (TR dDiv) const throw (cvmexception)
 Divide by real number operator.
BaseSCMatrix operator* (TC cMult) const
 Multiply by complex number operator.
BaseSCMatrix operator/ (TC cDiv) const throw (cvmexception)
 Divide by complex number operator.
basic_schmatrixoperator*= (TR dMult)
 Multiply by real number and assign.
basic_schmatrixoperator/= (TR dDiv) throw (cvmexception)
 Divide by real number and assign.
basic_schmatrixnormalize ()
 Matrix normalizer.
basic_schmatrix operator! () const throw (cvmexception)
 Matrix transposition.
basic_schmatrix operator~ () const
 Does nothing and returns copy of calling hermitian matrix.
basic_schmatrixtranspose (const basic_schmatrix &m) throw (cvmexception)
 Matrix transposition.
basic_schmatrixconj (const basic_schmatrix &m) throw (cvmexception)
 Assigns hermitian matrix m to calling one and returns reference to the matrix changed.
basic_schmatrixtranspose () throw (cvmexception)
 Matrix transposition (in-place)
basic_schmatrixconj ()
 Does nothing and returns reference to calling symmetric matrix.
CVector operator* (const CVector &v) const throw (cvmexception)
 Matrix-vector product.
BaseCMatrix operator* (const BaseCMatrix &m) const throw (cvmexception)
 Matrix-matrix product.
BaseSCMatrix operator* (const BaseSCMatrix &m) const throw (cvmexception)
 Matrix-matrix product.
CVector operator/ (const CVector &v) const throw (cvmexception)
 Linear solver operator.
basic_schmatrixherk (TR alpha, const CVector &v, TR beta) throw (cvmexception)
 Rank-1 update matrix-vector operation.
basic_schmatrixherk (bool bTransp, TR alpha, const BaseCMatrix &m, TR beta) throw (cvmexception)
 Rank-1 update matrix-matrix operation.
basic_schmatrixher2k (TC alpha, const CVector &v1, const CVector &v2, TR beta) throw (cvmexception)
 Rank-1 update matrix-vector operation.
basic_schmatrixher2k (bool bTransp, TC alpha, const BaseCMatrix &m1, const BaseCMatrix &m2, TR beta) throw (cvmexception)
 Rank-1 update matrix-matrix operation.
basic_schmatrixinv (const basic_schmatrix &m) throw (cvmexception)
 Matrix inversion.
basic_schmatrix inv () const throw (cvmexception)
 Matrix inversion.
basic_schmatrixexp (const basic_schmatrix &m, TR tol=basic_cvmMachSp< TR >()) throw (cvmexception)
 Matrix exponent.
basic_schmatrix exp (TR tol=basic_cvmMachSp< TR >()) const throw (cvmexception)
 Matrix exponent.
basic_schmatrixpolynom (const basic_schmatrix &m, const RVector &v) throw (cvmexception)
 Matrix polynomial.
basic_schmatrix polynom (const RVector &v) const
 Matrix polynomial.
RVector eig (BaseSCMatrix &mEigVect) const throw (cvmexception)
 Eigenvalues and eigenvectors.
RVector eig () const throw (cvmexception)
 Eigenvalues.
BaseSCMatrix cholesky () const throw (cvmexception)
 Cholesky factorization.
BaseSCMatrix bunch_kaufman (tint *nPivots) const throw (cvmexception)
 Bunch-Kaufman factorization.
basic_schmatrixidentity ()
 Identity matrix.
basic_schmatrixvanish ()
 Set matrix to zero.
basic_schmatrixrandomize_real (TR dFrom, TR dTo)
 Randomizer (real part)
basic_schmatrixrandomize_imag (TR dFrom, TR dTo)
 Randomizer (imaginary part)
TR norminf () const override
 Infinity norm.
bool is_positive_definite () const
 Is calling hermitian matrix positive definite?
bool equilibrate (basic_rvector< TR > &vScalings, CVector &vB) throw (cvmexception)
 Matrix equilibration.
bool equilibrate (basic_rvector< TR > &vScalings, BaseCMatrix &mB) throw (cvmexception)
 Matrix equilibration.
- Public Member Functions inherited from basic_scmatrix< TR, TC >
 basic_scmatrix ()
 Default constructor.
 basic_scmatrix (tint nDim)
 Constructor.
 basic_scmatrix (TC *pd, tint nDim)
 Constructor.
 basic_scmatrix (const TC *pd, tint nDim)
 Constructor.
 basic_scmatrix (const basic_scmatrix &m)
 Copy constructor.
 basic_scmatrix (basic_scmatrix &&m) noexcept
 Move constructor.
 basic_scmatrix (const BaseCMatrix &m)
 Constructor.
 basic_scmatrix (const CVector &v)
 Constructor.
 basic_scmatrix (const basic_srmatrix< TR > &m, bool bRealPart=true)
 Constructor.
 basic_scmatrix (const TR *pRe, const TR *pIm, tint nDim)
 Constructor.
 basic_scmatrix (const basic_srmatrix< TR > &mRe, const basic_srmatrix< TR > &mIm)
 Constructor.
 basic_scmatrix (BaseCMatrix &m, tint nRow, tint nCol, tint nDim)
 Submatrix constructor.
type_proxy< TC, TR > operator() (tint nRow, tint nCol) throw (cvmexception)
 Reference to element (l-value)
CVector operator() (tint nCol) throw (cvmexception)
 Column as l-value.
CVector operator[] (tint nRow) throw (cvmexception)
 Row as l-value.
basic_scmatrixoperator= (const basic_scmatrix &m) throw (cvmexception)
 Assignment operator.
basic_scmatrixoperator= (basic_scmatrix &&m) throw (cvmexception)
 Move assignment operator.
basic_scmatrixassign (const CVector &v) throw (cvmexception)
 Vector (as array) assignment.
basic_scmatrixassign (const TC *pd)
 External array assignment.
basic_scmatrixassign (tint nRow, tint nCol, const BaseCMatrix &m) throw (cvmexception)
 Assignment to submatrix.
basic_scmatrixset (TC c)
 Sets all elements to one value.
basic_scmatrixassign_real (const basic_srmatrix< TR > &mRe) throw (cvmexception)
 Assignment to real parts.
basic_scmatrixassign_imag (const basic_srmatrix< TR > &mIm) throw (cvmexception)
 Assignment to imaginary parts.
basic_scmatrixset_imag (TR d)
 Sets all imaginary parts to one value.
basic_scmatrixoperator<< (const basic_scmatrix &m) throw (cvmexception)
 Matrix replacement.
basic_scmatrix operator+ (const basic_scmatrix &m) const throw (cvmexception)
 Addition operator.
basic_scmatrix operator- (const basic_scmatrix &m) const throw (cvmexception)
 Subtraction operator.
basic_scmatrixsum (const basic_scmatrix &m1, const basic_scmatrix &m2) throw (cvmexception)
 Sum of matrices.
basic_scmatrixdiff (const basic_scmatrix &m1, const basic_scmatrix &m2) throw (cvmexception)
 Difference of matrices.
basic_scmatrixoperator+= (const basic_scmatrix &m) throw (cvmexception)
 Increment operator.
basic_scmatrixoperator-= (const basic_scmatrix &m) throw (cvmexception)
 Decrement operator.
basic_scmatrixoperator*= (TC cMult)
 Multiply by complex number and assign.
basic_scmatrixoperator/= (TC cDiv)
 Divide by complex number and assign.
basic_scmatrixtranspose (const basic_scmatrix &m) throw (cvmexception)
 Matrix transposition.
basic_scmatrixconj (const basic_scmatrix &m) throw (cvmexception)
 Matrix conjugation.
basic_scmatrix operator* (const basic_scmatrix &m) const throw (cvmexception)
 Matrix-matrix product.
basic_scmatrixoperator*= (const basic_scmatrix &m) throw (cvmexception)
 Matrix-matrix product with assignment.
basic_scmatrixswap_rows (tint n1, tint n2) throw (cvmexception)
 Rows swap.
basic_scmatrixswap_cols (tint n1, tint n2) throw (cvmexception)
 Columns swap.
CVector solve (const CVector &vB, TR &dErr) const throw (cvmexception)
 Linear solver.
CVector solve_tran (const CVector &vB, TR &dErr) const throw (cvmexception)
 Linear solver (transposed)
CVector solve_conj (const CVector &vB, TR &dErr) const throw (cvmexception)
 Linear solver (hermitian conjugated)
CVector solve (const CVector &vB) const throw (cvmexception)
 Linear solver.
CVector solve_tran (const CVector &vB) const throw (cvmexception)
 Linear solver (transposed)
CVector solve_conj (const CVector &vB) const throw (cvmexception)
 Linear solver (hermitian conjugated)
BaseCMatrix solve (const BaseCMatrix &mB, TR &dErr) const throw (cvmexception)
 Linear solver.
BaseCMatrix solve_tran (const BaseCMatrix &mB, TR &dErr) const throw (cvmexception)
 Linear solver (transposed)
BaseCMatrix solve_conj (const BaseCMatrix &mB, TR &dErr) const throw (cvmexception)
 Linear solver (hermitian conjugated)
BaseCMatrix solve (const BaseCMatrix &mB) const throw (cvmexception)
 Linear solver.
BaseCMatrix solve_tran (const BaseCMatrix &mB) const throw (cvmexception)
 Linear solver (transposed)
BaseCMatrix solve_conj (const BaseCMatrix &mB) const throw (cvmexception)
 Linear solver (hermitian conjugated)
CVector operator% (const CVector &vB) const throw (cvmexception)
 Linear solver operator (transposed)
CVector solve_lu (const basic_scmatrix &mLU, const tint *pPivots, const CVector &vB, TR &dErr) const throw (cvmexception)
 LU factorization based linear solver.
CVector solve_lu (const basic_scmatrix &mLU, const tint *pPivots, const CVector &vB) const throw (cvmexception)
 LU factorization based linear solver.
BaseCMatrix solve_lu (const basic_scmatrix &mLU, const tint *pPivots, const BaseCMatrix &mB, TR &dErr) const throw (cvmexception)
 LU factorization based linear solver.
BaseCMatrix solve_lu (const basic_scmatrix &mLU, const tint *pPivots, const BaseCMatrix &mB) const throw (cvmexception)
 LU factorization based linear solver.
TC det () const throw (cvmexception)
 Matrix determinant.
basic_scmatrixlow_up (const basic_scmatrix &m, tint *nPivots) throw (cvmexception)
 Low-up (LU) factorization.
basic_scmatrix low_up (tint *nPivots) const throw (cvmexception)
 Low-up (LU) factorization.
TR cond () const throw (cvmexception)
 Condition number reciprocal.
basic_scmatrixinv (const basic_scmatrix &m) throw (cvmexception)
 Matrix inversion.
basic_scmatrixexp (const basic_scmatrix &m, TR tol=basic_cvmMachSp< TR >()) throw (cvmexception)
 Matrix exponent.
basic_scmatrixpolynom (const basic_scmatrix &m, const CVector &v) throw (cvmexception)
 Matrix polynomial.
basic_scmatrix polynom (const CVector &v) const
 Matrix polynomial.
CVector eig (basic_scmatrix &mEigVect, bool bRightVect=true) const throw (cvmexception)
 Eigenvalues and eigenvectors.
basic_scmatrixcholesky (const basic_schmatrix< TR, TC > &m) throw (cvmexception)
 Cholesky factorization.
basic_scmatrixbunch_kaufman (const basic_schmatrix< TR, TC > &m, tint *nPivots) throw (cvmexception)
 Bunch-Kaufman factorization.
void qr (basic_scmatrix< TR, TC > &mQ, basic_scmatrix< TR, TC > &mR) const throw (cvmexception)
 QR factorization.
void lq (basic_scmatrix< TR, TC > &mL, basic_scmatrix< TR, TC > &mQ) const throw (cvmexception)
 LQ factorization.
void ql (basic_scmatrix< TR, TC > &mQ, basic_scmatrix< TR, TC > &mL) const throw (cvmexception)
 QL factorization.
void rq (basic_scmatrix< TR, TC > &mR, basic_scmatrix< TR, TC > &mQ) const throw (cvmexception)
 RQ factorization.
- Public Member Functions inherited from basic_cmatrix< TR, TC >
 basic_cmatrix ()
 Default constructor.
 basic_cmatrix (tint nM, tint nN)
 Constructor.
 basic_cmatrix (TC *pd, tint nM, tint nN)
 Constructor.
 basic_cmatrix (const TC *pd, tint nM, tint nN)
 Constructor.
 basic_cmatrix (const basic_cmatrix &m)
 Copy constructor.
 basic_cmatrix (basic_cmatrix &&m) noexcept
 Move constructor.
 basic_cmatrix (const CVector &v, bool bBeColumn=true)
 Constructor.
 basic_cmatrix (const basic_rmatrix< TR > &m, bool bRealPart=true)
 Constructor.
 basic_cmatrix (const TR *pRe, const TR *pIm, tint nM, tint nN)
 Constructor.
 basic_cmatrix (const basic_rmatrix< TR > &mRe, const basic_rmatrix< TR > &mIm)
 Constructor.
 basic_cmatrix (basic_cmatrix &m, tint nRow, tint nCol, tint nHeight, tint nWidth)
 Submatrix constructor.
CVector diag (tint nDiag) throw (cvmexception)
 Diagonal as l-value.
basic_cmatrixoperator= (const basic_cmatrix &m) throw (cvmexception)
 Assignment operator.
basic_cmatrixoperator= (basic_cmatrix &&m) throw (cvmexception)
 Move assignment operator.
basic_cmatrixassign (tint nRow, tint nCol, const basic_cmatrix &m) throw (cvmexception)
 Assignment to submatrix.
basic_cmatrixassign_real (const basic_rmatrix< TR > &mRe) throw (cvmexception)
 Assignment to real parts.
basic_cmatrixassign_imag (const basic_rmatrix< TR > &mIm) throw (cvmexception)
 Assignment to imaginary parts.
basic_cmatrixresize (tint nNewM, tint nNewN) throw (cvmexception)
 Changes dimensions.
bool operator== (const basic_cmatrix &m) const
 Matrix comparison.
bool operator!= (const basic_cmatrix &m) const
 Matrix comparison.
basic_cmatrixoperator<< (const basic_cmatrix &m) throw (cvmexception)
 Matrix replacement.
basic_cmatrix operator+ (const basic_cmatrix &m) const throw (cvmexception)
 Addition operator.
basic_cmatrix operator- (const basic_cmatrix &m) const throw (cvmexception)
 Subtraction operator.
basic_cmatrixsum (const basic_cmatrix &m1, const basic_cmatrix &m2) throw (cvmexception)
 Sum of matrices.
basic_cmatrixdiff (const basic_cmatrix &m1, const basic_cmatrix &m2) throw (cvmexception)
 Difference of matrices.
basic_cmatrixoperator+= (const basic_cmatrix &m) throw (cvmexception)
 Increment operator.
basic_cmatrixoperator-= (const basic_cmatrix &m) throw (cvmexception)
 Decrement operator.
basic_cmatrixtranspose (const basic_cmatrix &m) throw (cvmexception)
 Matrix transposition.
basic_cmatrixconj (const basic_cmatrix &m) throw (cvmexception)
 Matrix conjugation.
basic_cmatrix operator* (const basic_cmatrix &m) const throw (cvmexception)
 Matrix-matrix product.
basic_cmatrixmult (const basic_cmatrix &m1, const basic_cmatrix &m2) throw (cvmexception)
 Matrix-matrix product.
basic_cmatrixrank1update_u (const CVector &vCol, const CVector &vRow) throw (cvmexception)
 Rank-1 update (unconjugated)
basic_cmatrixrank1update_c (const CVector &vCol, const CVector &vRow) throw (cvmexception)
 Rank-1 update (conjugated)
basic_cmatrixsolve (const basic_scmatrix< TR, TC > &mA, const basic_cmatrix &mB, TR &dErr) throw (cvmexception)
 Linear solver.
basic_cmatrixsolve_tran (const basic_scmatrix< TR, TC > &mA, const basic_cmatrix &mB, TR &dErr) throw (cvmexception)
 Linear solver (transposed)
basic_cmatrixsolve_conj (const basic_scmatrix< TR, TC > &mA, const basic_cmatrix &mB, TR &dErr) throw (cvmexception)
 Linear solver (conjugated)
basic_cmatrixsolve (const basic_scmatrix< TR, TC > &mA, const basic_cmatrix &mB) throw (cvmexception)
 Linear solver.
basic_cmatrixsolve_tran (const basic_scmatrix< TR, TC > &mA, const basic_cmatrix &mB) throw (cvmexception)
 Linear solver (transposed)
basic_cmatrixsolve_conj (const basic_scmatrix< TR, TC > &mA, const basic_cmatrix &mB) throw (cvmexception)
 Linear solver (conjugated)
basic_cmatrixsolve_lu (const basic_scmatrix< TR, TC > &mA, const basic_scmatrix< TR, TC > &mLU, const tint *pPivots, const basic_cmatrix &mB, TR &dErr) throw (cvmexception)
 LU factorization based linear solver.
basic_cmatrixsolve_lu (const basic_scmatrix< TR, TC > &mA, const basic_scmatrix< TR, TC > &mLU, const tint *pPivots, const basic_cmatrix &mB) throw (cvmexception)
 LU factorization based linear solver.
RVector svd () const throw (cvmexception)
 Singular value decomposition.
RVector svd (basic_scmatrix< TR, TC > &mU, basic_scmatrix< TR, TC > &mVH) const throw (cvmexception)
 Singular value decomposition.
basic_cmatrix pinv (TR threshold=basic_cvmMachSp< TR >()) const throw (cvmexception)
 Pseudo (generalized) inversion.
basic_cmatrixpinv (const basic_cmatrix &mA, TR threshold=basic_cvmMachSp< TR >()) throw (cvmexception)
 Pseudo (generalized) inversion.
basic_cmatrix gels (bool conjugate, const basic_cmatrix &mB, basic_cvector< TR, TC > &vErr) const throw (cvmexception)
 Overdetermined or underdetermined linear solver.
basic_cmatrixgels (bool conjugate, const basic_cmatrix &mA, const basic_cmatrix &mB, basic_cvector< TR, TC > &vErr) throw (cvmexception)
 Overdetermined or underdetermined linear solver.
basic_cvector< TR, TC > gels (bool conjugate, const basic_cvector< TR, TC > &vB, TC &dErr) const throw (cvmexception)
 Overdetermined or underdetermined linear solver.
basic_cmatrix gelsy (const basic_cmatrix &mB, tint &rank, TR tol=basic_cvmMachSp< TR >()) const throw (cvmexception)
 Linear least squares problem.
basic_cmatrixgelsy (const basic_cmatrix &mA, const basic_cmatrix &mB, tint &rank, TR tol=basic_cvmMachSp< TR >()) throw (cvmexception)
 Linear least squares problem.
basic_cvector< TR, TC > gelsy (const basic_cvector< TR, TC > &vB, tint &rank, TR tol=basic_cvmMachSp< TR >()) const throw (cvmexception)
 Linear least squares problem.
basic_cmatrix gelss (const basic_cmatrix &mB, basic_rvector< TR > &sv, tint &rank, TR tol=basic_cvmMachSp< TR >()) const throw (cvmexception)
 Linear least squares problem.
basic_cmatrixgelss (const basic_cmatrix &mA, const basic_cmatrix &mB, basic_rvector< TR > &sv, tint &rank, TR tol=basic_cvmMachSp< TR >()) throw (cvmexception)
 Linear least squares problem.
basic_cvector< TR, TC > gelss (const basic_cvector< TR, TC > &vB, basic_rvector< TR > &sv, tint &rank, TR tol=basic_cvmMachSp< TR >()) const throw (cvmexception)
 Linear least squares problem.
basic_cmatrix gelsd (const basic_cmatrix &mB, basic_rvector< TR > &sv, tint &rank, TR tol=basic_cvmMachSp< TR >()) const throw (cvmexception)
 Linear least squares problem.
basic_cmatrixgelsd (const basic_cmatrix &mA, const basic_cmatrix &mB, basic_rvector< TR > &sv, tint &rank, TR tol=basic_cvmMachSp< TR >()) throw (cvmexception)
 Linear least squares problem.
basic_cvector< TR, TC > gelsd (const basic_cvector< TR, TC > &vB, basic_rvector< TR > &sv, tint &rank, TR tol=basic_cvmMachSp< TR >()) const throw (cvmexception)
 Linear least squares problem.
tint rank (TR tol=basic_cvmMachSp< TR >()) const throw (cvmexception)
 Matrix rank.
void qr (basic_cmatrix< TR, TC > &mQ, basic_scmatrix< TR, TC > &mR) const throw (cvmexception)
 QR factorization ("economy" mode)
void qr (basic_scmatrix< TR, TC > &mQ, basic_cmatrix< TR, TC > &mR) const throw (cvmexception)
 QR factorization ("full" mode)
void lq (basic_scmatrix< TR, TC > &mL, basic_cmatrix< TR, TC > &mQ) const throw (cvmexception)
 LQ factorization ("economy" mode)
void lq (basic_cmatrix< TR, TC > &mL, basic_scmatrix< TR, TC > &mQ) const throw (cvmexception)
 LQ factorization ("full" mode)
void rq (basic_scmatrix< TR, TC > &mR, basic_cmatrix< TR, TC > &mQ) const throw (cvmexception)
 RQ factorization ("economy" mode)
void rq (basic_cmatrix< TR, TC > &mR, basic_scmatrix< TR, TC > &mQ) const throw (cvmexception)
 RQ factorization ("full" mode)
void ql (basic_cmatrix< TR, TC > &mQ, basic_scmatrix< TR, TC > &mL) const throw (cvmexception)
 QL factorization ("economy" mode)
void ql (basic_scmatrix< TR, TC > &mQ, basic_cmatrix< TR, TC > &mL) const throw (cvmexception)
 QL factorization ("full" mode)
basic_cmatrixgeru (TC alpha, const CVector &vCol, const CVector &vRow) throw (cvmexception)
 Rank-1 update matrix-vector operation (unconjugated)
basic_cmatrixgerc (TC alpha, const CVector &vCol, const CVector &vRow) throw (cvmexception)
 Rank-1 update matrix-vector operation (conjugated)
basic_cmatrixgemm (const basic_cmatrix &m1, bool bConj1, const basic_cmatrix &m2, bool bConj2, TC cAlpha, TC cBeta) throw (cvmexception)
 Generic matrix-matrix operation.
basic_cmatrixhemm (bool bLeft, const basic_schmatrix< TR, TC > &ms, const basic_cmatrix &m, TC cAlpha, TC cBeta) throw (cvmexception)
 Generic hermitian matrix-matrix operation.
TR norm2 () const override
 2-norm
- Public Member Functions inherited from Matrix< TR, TC >
tint msize () const
 Number of rows.
tint nsize () const
 Number of columns.
tint ld () const
 Leading dimension.
tint rowofmax () const
 Row with maximum element.
tint rowofmin () const
 Row with minimum element.
tint colofmax () const
 Column with maximum element.
tint colofmin () const
 Column with minimum element.
TR norm1 () const override
 1-norm
- Public Member Functions inherited from basic_array< TR, TC >
 basic_array ()
 Default constructor.
 basic_array (tint nSize, bool bZeroMemory=true)
 Constructor.
 basic_array (TC *pd, tint nSize, tint nIncr=1)
 Constructor.
 basic_array (const TC *pd, tint nSize, tint nIncr=1)
 Constructor.
 basic_array (const TC *begin, const TC *end)
 Constructor.
 basic_array (const basic_array &a)
 Copy constructor.
 basic_array (basic_array &&a) noexcept
 Move constructor.
basic_arrayoperator= (const basic_array &a) throw (cvmexception)
 Assignment operator.
basic_arrayoperator= (basic_array &&a) throw (cvmexception)
 Move assignment operator.
virtual ~basic_array ()
 Destructor.
tint size () const
 Size (length) of array.
TC * get ()
 Pointer to data.
const TC * get () const
 Const pointer to data.
 operator TC * ()
 Type cast to pointer to data.
 operator const TC * () const
 Type cast to constant pointer to data.
tint incr () const
 Increment between elements.
tint indofmax () const
 Index of element with maximum module.
tint indofmin () const
 Index of element with minimum module.
virtual TR norm () const
 Euclidean norm.
iterator begin ()
 (STL) iterator to begin
const_iterator begin () const
 (STL) const iterator to begin
iterator end ()
 (STL) iterator to end
const_iterator end () const
 (STL) const iterator to end
reverse_iterator rbegin ()
 (STL) iterator to begin reversed
const_reverse_iterator rbegin () const
 (STL) const iterator to begin reversed
reverse_iterator rend ()
 (STL) iterator to end reversed
const_reverse_iterator rend () const
 (STL) const iterator to end reversed
size_type max_size () const
 (STL) maximum possible size of array
size_type capacity () const
 (STL) current capacity of array, equal to size()
bool empty () const
 (STL) is array empty
reference front ()
 (STL) reference to first element
const_reference front () const
 (STL) const reference to first element
reference back ()
 (STL) reference to last element
const_reference back () const
 (STL) const reference to last element
void assign (size_type n, const TC &val) throw (cvmexception)
 (STL) assigns given value to n-th element (0-based)
void assign (const_iterator begin, const_iterator end) throw (cvmexception)
 (STL) assigns begin-end iteartor range to array
void clear ()
 (STL) clears array, dealocates memory and sets size() to zero
void swap (basic_array &v) throw (cvmexception)
 (STL) swaps array values, throws cvmexception if sizes are different
reference at (size_type n) throw (cvmexception)
 (STL) returns reference to n-th element of array (0-based), throws cvmexception if n is out of boundaries
const_reference at (size_type n) const throw (cvmexception)
 (STL) returns const reference to n-th element of array (0-based), throws cvmexception if n is out of boundaries
void push_back (const TC &x) throw (cvmexception)
 (STL) pushes new value to the end of array
void pop_back () throw (cvmexception)
 (STL) removes last element from array
iterator insert (iterator position, const TC &x) throw (cvmexception)
 (STL) inserts element to given position in array
iterator erase (iterator position) throw (cvmexception)
 (STL) removes element from given position in array

Additional Inherited Members

- Public Types inherited from basic_array< TR, TC >
typedef TC value_type
 STL-specific value type definition.
typedef value_typepointer
 STL-specific value pointer definition.
typedef value_typeiterator
 STL-specific iterator definition.
typedef const value_typeconst_iterator
 STL-specific const iterator definition.
typedef const value_typeconst_pointer
 STL-specific const pointer definition.
typedef value_typereference
 STL-specific reference definition.
typedef const value_typeconst_reference
 STL-specific const reference definition.
typedef size_t size_type
 STL-specific size type definition.
typedef ptrdiff_t difference_type
 STL-specific difference type definition.
typedef std::reverse_iterator
< const_iterator
const_reverse_iterator
 STL-specific const reverse iterator definition.
typedef std::reverse_iterator
< iterator
reverse_iterator
 STL-specific reverse iterator definition.
- Protected Member Functions inherited from Matrix< TR, TC >
 Matrix ()
 Default constructor.
 Matrix (tint nM, tint nN, tint nLD, bool bZeroMemory)
 Constructor.
 Matrix (TC *pd, tint nM, tint nN, tint nLD, tint nSize)
 Constructor.
 Matrix (const TC *pd, tint nM, tint nN, tint nLD, tint nSize)
 Constructor.
 Matrix (const BaseArray &v, bool bBeColumn)
 Constructor.
 Matrix (Matrix &&m) noexcept
 Move constructor.
Matrixoperator= (Matrix &&m) throw (cvmexception)
 Move assignment operator.
- Protected Member Functions inherited from SqMatrix< TR, TC >
 SqMatrix ()
 internal constructor
virtual ~SqMatrix ()
- Protected Attributes inherited from Matrix< TR, TC >
tint mm
 Number of rows.
tint mn
 Number of columns.
tint mld
 Leading dimension.

Detailed Description

template<typename TR, typename TC>
class basic_schmatrix< TR, TC >

End-user class encapsulating hermitian matrix of complex numbers.

TR type stands for treal, TC type stands for tcomplex. Please use predefined schmatrix class in your applications.

See Also
Matrix

Definition at line 34241 of file cvm.h.

Constructor & Destructor Documentation

template<typename TR, typename TC>
basic_schmatrix< TR, TC >::basic_schmatrix ( )
inline

Default constructor.

Creates empty hermitian matrix. No memory gets allocated.

Example:
using namespace cvm;
std::cout << m.msize() << " " << m.nsize() << " "
<< m.size() << std::endl;
m.resize (3);
std::cout << m;
prints
0 0 0
(0,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
(0,0) (0,0) (0,0)

Definition at line 34271 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix< TR, TC >::basic_schmatrix ( tint  nDim)
inlineexplicit

Constructor.

Creates $n\times n$ schmatrix object where $n$ is passed in nDim parameter. Constructor sets all elements to zero. It throws cvmexception in case of non-positive size passed or memory allocation failure. Example:

using namespace cvm;
schmatrix m (4);
std::cout << m.msize() << std::endl
<< m.nsize() << std::endl
<< m.size() << std::endl << m;

prints

4
4
16
(0,0) (0,0) (0,0) (0,0)
(0,0) (0,0) (0,0) (0,0)
(0,0) (0,0) (0,0) (0,0)
(0,0) (0,0) (0,0) (0,0)
Parameters
[in]nDimNumber of rows and columns.

Definition at line 34301 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix< TR, TC >::basic_schmatrix ( TC *  pd,
tint  nDim,
TR  tol = basic_cvmMachSp<TR>() 
)
inline

Constructor.

Creates $n\times n$ schmatrix object where $n$ is passed in nDim parameter. It throws cvmexception in case of non-positive size passed or if the matrix created doesn't appear to be hermitian. Hermiticity tolerance is set by parameter tol. Unlike others, this constructor does not allocate memory. It just shares memory with array pointed to by pd (for matrices nIncr=1 is always satisfied). If subsequent application flow would change the array passed so it becomes not hermitian matrix anymore then results are not predictable.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m ((std::complex<double>*)a, 3);
m.set(2,1,std::complex<double>(8.,8.));
std::cout << m << std::endl;
std::cout << a[0] << " " << a[1] << " " << a[2] << " "
<< a[3] << " " << a[4] << " " << a[5] << " " << std::endl;
prints
(1,0) (8,-8) (-1,-2)
(8,8) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
1 0 8 8 -1 2
See Also
http://cvmlib.com/faq.htm
basic_cmatrix::basic_cmatrix(TC*,tint,tint)
Parameters
[in]pdPointer to array to share memory with.
[in]nDimNumber of rows and columns.
[in]tolHermiticity tolerance.

Definition at line 34342 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix< TR, TC >::basic_schmatrix ( const TC *  pd,
tint  nDim,
TR  tol = basic_cvmMachSp<TR>() 
)
inline

Constructor.

Creates $n\times n$ schmatrix object where $n$ is passed in nDim parameter. Constructor throws cvmexception in case of non-positive sizes passed or if the matrix created doesn't appear to be hermitian. Hermiticity tolerance is set by parameter tol. This is const version, it allocates memory and copies every element (deep copy) from external array pointed to by pd parameter. It copies nDim*nDim elements total.

Example:
using namespace cvm;
const double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m ((const std::complex<double>*)a, 3);
m.set(2,1,std::complex<double>(8.,8.));
std::cout << m << std::endl;
std::cout << a[0] << " " << a[1] << " " << a[2] << " "
<< a[3] << " " << a[4] << " " << a[5] << " " << std::endl;
prints
(1,0) (8,-8) (-1,-2)
(8,8) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
1 0 2 1 -1 2
See Also
http://cvmlib.com/faq.htm
basic_cmatrix::basic_cmatrix(const TC*,tint,tint)
Parameters
[in]pdConst pointer to external array.
[in]nDimNumber of rows and columns.
[in]tolHermiticity tolerance.

Definition at line 34383 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix< TR, TC >::basic_schmatrix ( const basic_schmatrix< TR, TC > &  m)
inline

Copy constructor.

Creates schmatrix object as a copy of hermitian matrix m. It throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left |
std::ios::showpos);
std::cout.precision (1);
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m ((std::complex<double>*)a, 3);
scmatrix mc(m);
m.set(1,2, std::complex<double>(7.7,7.7));
std::cout << m << std::endl << mc;
prints
(+1.0e+000,+0.0e+000) (+7.7e+000,+7.7e+000) (-1.0e+000,-2.0e+000)
(+7.7e+000,-7.7e+000) (+2.0e+000,+0.0e+000) (+0.0e+000,-3.0e+000)
(-1.0e+000,+2.0e+000) (+0.0e+000,+3.0e+000) (+3.0e+000,+0.0e+000)
(+1.0e+000,+0.0e+000) (+2.0e+000,-1.0e+000) (-1.0e+000,-2.0e+000)
(+2.0e+000,+1.0e+000) (+2.0e+000,+0.0e+000) (+0.0e+000,-3.0e+000)
(-1.0e+000,+2.0e+000) (+0.0e+000,+3.0e+000) (+3.0e+000,+0.0e+000)
Parameters
[in]mschmatrix to copy from.

Definition at line 34419 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix< TR, TC >::basic_schmatrix ( basic_schmatrix< TR, TC > &&  m)
inline

Move constructor.

Implements move semantics introduced in new C++ standard. Moves data ownership from other matrix to newly created object. It's usually executed implicitly in cases like this:

rvector a(b + c);

or this

rvector a = b + c;

Here temporary result of calling b.operator+(c) will not be destroyed but rather moved to newly created object a.

Parameters
[in]mrvalue reference to other matrix.

Definition at line 34442 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix< TR, TC >::basic_schmatrix ( const BaseCMatrix m,
TR  tol = basic_cvmMachSp<TR>() 
)
inlineexplicit

Constructor.

Creates schmatrix object as a copy of complex matrix m. It's assumed that $m\times n$ matrix m must have equal sizes, i.e. $m = n$ is satisfied, and must be hermitian. Hermiticity tolerance is set by parameter tol. Constructor throws cvmexception if this is not true or in case of memory allocation failure.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
cmatrix m ((std::complex<double>*)a, 3, 3);
scmatrix mch(m);
std::cout << mch;
prints
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
Parameters
[in]mcmatrix to copy from.
[in]tolHermiticity tolerance.

Definition at line 34473 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix< TR, TC >::basic_schmatrix ( const RVector v)
inlineexplicit

Constructor.

Creates schmatrix object of size v.size() by v.size() and assigns real vector v to its main diagonal. Constructor throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
double a[] = {1., 2., 3., 4., 5.};
const rvector v (a, 5);
schmatrix m(v);
std::cout << m.msize() << " " << m.nsize() << " "
<< m.size() << std::endl << m;
prints
5 5 25
(1,0) (0,0) (0,0) (0,0) (0,0)
(0,0) (2,0) (0,0) (0,0) (0,0)
(0,0) (0,0) (3,0) (0,0) (0,0)
(0,0) (0,0) (0,0) (4,0) (0,0)
(0,0) (0,0) (0,0) (0,0) (5,0)
Parameters
[in]vrvector to copy main diagonal from.

Definition at line 34507 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix< TR, TC >::basic_schmatrix ( const basic_srsmatrix< TR > &  m)
inlineexplicit

Constructor.

Creates schmatrix object having the same dimension as real symmetric matrix m and copies matrix m to its real part.

Example:
using namespace cvm;
double a[] = {1., 2., 3., 2., 5., 6., 3., 6., 9.};
srsmatrix m(a, 3);
schmatrix mch(m);
std::cout << mch;
prints
(1,0) (2,0) (3,0)
(2,0) (5,0) (6,0)
(3,0) (6,0) (9,0)
Parameters
[in]msrsmatrix to copy from.

Definition at line 34534 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix< TR, TC >::basic_schmatrix ( const TR *  pRe,
const TR *  pIm,
tint  nDim,
TR  tol = basic_cvmMachSp<TR>() 
)
inline

Constructor.

Creates schmatrix object of dimension nDim and copies every element of arrays pointed to by pRe and pIm to real and imaginary part of the matrix created respectively. Use nullptr pointer to fill up appropriate part with zero values. Contructor checks that the matrix is hermitian and throws cvmexception otherwise. Hermiticity tolerance is set by parameter tol. Constructor also throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
double re[] = {1., 2., -1., 2., 2., 0., -1., 0., 3.};
double im[] = {0., 1., 2., -1., 0., 3., -2., -3., 0.};
schmatrix m(re, im, 3);
std::cout << m;
prints
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
Parameters
[in]pReConst pointer to external treal array to copy to real part.
[in]pImConst pointer to external treal array to copy to imaginary part.
[in]nDimMatrix dimension.
[in]tolHermiticity tolerance.

Definition at line 34569 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix< TR, TC >::basic_schmatrix ( const basic_srmatrix< TR > &  mRe,
const basic_srmatrix< TR > &  mIm,
TR  tol = basic_cvmMachSp<TR>() 
)
inline

Constructor.

Creates schmatrix object of dimension mRe.msize() (if it differs from mIm.msize() then constructor throws cvmexception) and copies matrices mRe and mIm to real and imaginary part of the matrix created respectively. Contructor checks that the matrix is hermitian and throws cvmexception otherwise. Hermiticity tolerance is set by parameter tol. Constructor also throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
double re[] = {1., 2., -1., 2., 2., 0., -1., 0., 3.};
double im[] = {0., 1., 2., -1., 0., 3., -2., -3., 0.};
srmatrix mr(re, 3);
srmatrix mi(im, 3);
schmatrix m(mr, mi);
std::cout << m;
prints
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
Parameters
[in]mResrmatrix to copy to real part.
[in]mImsrmatrix to copy to imaginary part.
[in]tolHermiticity tolerance.

Definition at line 34604 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix< TR, TC >::basic_schmatrix ( basic_schmatrix< TR, TC > &  m,
tint  nRowCol,
tint  nDim 
)
inline

Submatrix constructor.

Creates schmatrix object as submatrix of hermitian matrix m. It means that the object created shares memory with some part of m. This part is defined by its upper left corner (parameter nRowCol, CVM0 based) and its dimension (parameter nDim).

Example:
using namespace cvm;
srsmatrix m(5);
srsmatrix subm(m, 2, 2);
subm.set(1.);
std::cout << m;
prints
0 0 0 0 0
0 1 1 0 0
0 1 1 0 0
0 0 0 0 0
Parameters
[in]mParent schmatrix to attach to.
[in]nRowColRow and column to start from (CVM0 based).
[in]nDimDimension of square submatrix.

Definition at line 34636 of file cvm.h.

Member Function Documentation

template<typename TR, typename TC>
TC basic_schmatrix< TR, TC >::operator() ( tint  nRow,
tint  nCol 
) const throw (cvmexception)
inline

Value of element (not l-value)

Operator returns value of a particular element of calling matrix by its row and column index. Indexes passed are CVM0 based. Operator throws cvmexception if nRow or nCol is outside of boundaries.

Example:
using namespace cvm;
try {
double a[] = {1., 2., 3., 4., 5., 6.,
7., 8., 9., 10., 11., 12.};
const cmatrix m ((std::complex<double>*) a, 2, 3);
std::cout << m(1,1) << " " << m(2,3) << std::endl;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(1,2) (11,12)
See Also
http://cvmlib.com/faq.htm
Parameters
[in]nRowRow index (CVM0 based).
[in]nColColumn index (CVM0 based).
Returns
TC tcomplex value.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 34641 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
const CVector basic_schmatrix< TR, TC >::operator() ( tint  nCol) const throw (cvmexception)
inline

Column as not l-value.

Operator creates cvector object as a copy of nCol-th column (CVM0 based) of calling matrix. Operator throws cvmexception if nCol is outside of boundaries.

Example:
using namespace cvm;
try {
double a[] = {1., 2., 3., 4., 5., 6.,
7., 8., 9., 10., 11., 12.};
const cmatrix m ((std::complex<double>*) a, 2, 3);
std::cout << m(2) << std::endl;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(5,6) (7,8)
Parameters
[in]nColIndex of column (CVM0 based).
Returns
cvector Column value.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 34649 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
const CVector basic_schmatrix< TR, TC >::operator[] ( tint  nRow) const throw (cvmexception)
inline

Row as not l-value.

Operator creates cvector object as a copy of nRow-th row (CVM0 based) of calling matrix. Operator throws cvmexception if nRow is outside of boundaries.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
try {
double a[] = {1., 2., 3., 4., 5., 6., 7., 8., 9.};
const srmatrix m (a, 3);
std::cout << m[2] << std::endl;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
2.00e+00 5.00e+00 8.00e+00
Parameters
[in]nRowIndex of row (CVM0 based).
Returns
cvector Row value.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 34656 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
const CVector basic_schmatrix< TR, TC >::diag ( tint  nDiag) const throw (cvmexception)
inline

Diagonal (not l-value)

Operator creates cvector object as a copy of nDiag-th diagonal of calling matrix where nDiag=0 for main diagonal, nDiag<0 for lower diagonals and nDiag>0 for upper ones. Operator throws cvmexception if nDiag is outside of boundaries.

Example:
using namespace cvm;
try {
double a[] = {1., 2., 3., 4., 5., 6., 7., 8., 9.,
10., 11., 12., 13., 14., 15., 16., 17., 18.};
const scmatrix ms((std::complex<double>*)a, 3);
std::cout << ms << std::endl;
std::cout << ms.diag(0) << ms.diag(1);
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(1,2) (7,8) (13,14)
(3,4) (9,10) (15,16)
(5,6) (11,12) (17,18)
(1,2) (9,10) (17,18)
(7,8) (15,16)
Parameters
[in]nDiagIndex of diagonal (0 for main diagonal, negative for lower, positive for upper one).
Returns
cvector Diagonal value.

Reimplemented from basic_cmatrix< TR, TC >.

Definition at line 34664 of file cvm.h.

template<typename TR, typename TC>
const basic_srsmatrix<TR> basic_schmatrix< TR, TC >::real ( ) const
inline

Real part (not l-value)

Creates object of type const rmatrix as real part of calling matrix. Please note that, unlike cvector::real(), this function creates new object not sharing memory with real part of calling matrix, i.e. the matrix returned is not l-value.

Example:
using namespace cvm;
double a[] = {1., 2., 3., 4., 5., 6.,
7., 8., 9., 10., 11., 12.};
cmatrix m((std::complex<double>*) a, 2, 3);
std::cout << m << std::endl << m.real();
prints
(1,2) (5,6) (9,10)
(3,4) (7,8) (11,12)
1 5 9
3 7 11
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 34669 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
const basic_srmatrix<TR> basic_schmatrix< TR, TC >::imag ( ) const
inline

Imaginary part (not l-value)

Creates object of type const rmatrix as imaginary part of calling matrix. Please note that, unlike cvector::imag(), this function creates new object not sharing memory with imaginary part of calling matrix, i.e. the matrix returned is not l-value.

Example:
using namespace cvm;
double a[] = {1., 2., 3., 4., 5., 6.,
7., 8., 9., 10., 11., 12.};
cmatrix m((std::complex<double>*) a, 2, 3);
std::cout << m << std::endl << m.imag();
prints
(1,2) (5,6) (9,10)
(3,4) (7,8) (11,12)
2 6 10
4 8 12
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 34677 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::operator= ( const basic_schmatrix< TR, TC > &  m) throw (cvmexception)
inline

Assignment operator.

Sets every element of calling schmatrix to be equal to appropriate element of hermitian matrix m and returns reference to the matrix changed. Operator throws cvmexception in case of different matrix dimensions.

Example:
using namespace cvm;
try {
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m1((std::complex<double>*)a, 3);
schmatrix m2(3);
m2 = m1;
std::cout << m2;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
Parameters
[in]mschmatrix to assign from.
Returns
Reference to changed calling matrix.

Definition at line 34716 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::operator= ( basic_schmatrix< TR, TC > &&  m) throw (cvmexception)
inline

Move assignment operator.

Implements move semantics introduced in new C++ standard. Moves data ownership from other matrix to calling object. It's usually executed implicitly in cases like this:

a = b + c;

Here temporary result of calling b.operator+(c) will not be destroyed but rather moved to calling object a.

Parameters
[in]mrvalue reference to other matrix.

Definition at line 34736 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::assign ( const CVector v,
TR  tol = basic_cvmMachSp<TR>() 
) throw (cvmexception)
inline

Vector (as array) assignment.

Sets every element of calling hermitian matrix to be equal to appropriate element of rvector v as an array. Assignment is performed according to matrix storage (by columns). It's assumed that vector passed is long enough to fill calling matrix. Function throws cvmexception otherwise. Function also throws cvmexception if the matrix changed doesn't appear to be hermitian. Hermiticity tolerance is set by parameter tol.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
cvector v((std::complex<double>*)a, 9);
schmatrix m(3);
m.assign(v);
std::cout << m;
prints
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
Parameters
[in]vcvector to assign.
[in]tolHermiticity tolerance.
Returns
Reference to changed calling matrix.

Definition at line 34771 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::assign ( const TC *  pd,
TR  tol = basic_cvmMachSp<TR>() 
) throw (cvmexception)
inline

External array assignment.

Sets every element of calling hermitian matrix to be equal to appropriate element of array pointed to by parameter pd and returns reference to the matrix changed. Assignment is performed according to matrix storage (by columns). It's assumed that array passed is long enough to fill calling matrix. Function throws cvmexception if the matrix changed doesn't appear to be hermitian. Hermiticity tolerance is set by parameter tol.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m(3);
m.assign((std::complex<double>*)a);
std::cout << m;
prints
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
Parameters
[in]pdConst pointer to external array.
[in]tolHermiticity tolerance.
Returns
Reference to changed calling matrix.

Definition at line 34809 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::assign ( tint  nRowCol,
const basic_schmatrix< TR, TC > &  m 
) throw (cvmexception)
inline

Assignment to submatrix.

Sets main submatrix of calling hermitian matrix beginning with row and column nRowCol to hermitian matrix m and returns reference to the matrix changed. Function throws cvmexception if nRowCol is not positive or matrix m doesn't fit. Indexes are CVM0 based.

Example:
using namespace cvm;
schmatrix m1(5);
schmatrix m2(2);
m2.set_main_diag(rvector(2,2.));
m2.set(1,2,std::complex<double>(2.,2.));
m1.assign(2,m2);
std::cout << m1;
prints
(0,0) (0,0) (0,0) (0,0) (0,0)
(0,0) (2,0) (2,2) (0,0) (0,0)
(0,0) (2,-2) (2,0) (0,0) (0,0)
(0,0) (0,0) (0,0) (0,0) (0,0)
(0,0) (0,0) (0,0) (0,0) (0,0)
Parameters
[in]nRowColRow and column index (CVM0 based).
[in]mReference to hermitian matrix to assign.
Returns
Reference to changed calling matrix.

Definition at line 34846 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::set ( tint  nRow,
tint  nCol,
TC  c 
) throw (cvmexception)
inline

Sets both elements to keep matrix hermitian.

Sets main submatrix of calling hermitian matrix beginning with row and column nRowCol to hermitian matrix m and returns reference to the matrix changed. Function throws cvmexception if nRowCol is out of boundaries (it's CVM0 based).

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
m.set(1,2,std::complex<double>(7.7,7.7));
m.set(3,3,11.11);
std::cout << m;
prints
(1,0) (7.7,7.7) (-1,-2)
(7.7,-7.7) (2,0) (0,-3)
(-1,2) (0,3) (11.11,0)
Parameters
[in]nRowRow (and column) index to set (CVM0 based).
[in]nColColumn (and row) index to set (CVM0 based).
[in]cElement value to set to.
Returns
Reference to changed calling matrix.

Definition at line 34883 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::set_diag ( tint  nDiag,
const CVector vDiag 
) throw (cvmexception)
inline

Sets lower and upper diagonals.

Assigns complex vector vDiag to $i$-th diagonal of calling hermitian matrix ( $i$ is passed in nDiag parameter), where $i<0$ for lower diagonals and $i>0$ for upper ones. Function also assigns conjugated vector to $-i$-th diagonal (thus calling matrix remains hermitian). Function returns reference to the matrix changed. Function throws cvmexception if parameter nDiag is equal to zero (use set_main_diag() instead), is outside of boundaries or if vector vDiag passed has size not equal to msize()-abs(nDiag).

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
cvector v(2);
v.set(std::complex<double>(7.7,7.7));
m.set_diag(1, v);
std::cout << m;
prints
(1,0) (7.7,7.7) (-1,-2)
(7.7,-7.7) (2,0) (7.7,7.7)
(-1,2) (7.7,-7.7) (3,0)
Parameters
[in]nDiagDiagonal index to set.
[in]vDiagcvector to set diagonal(s) to.
Returns
Reference to changed calling matrix.

Definition at line 34922 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::set_main_diag ( const RVector vDiag) throw (cvmexception)
inline

Sets main diagonal.

Assigns real vector vDiag to main diagonal of calling hermitian matrix. Function returns reference to the matrix changed. Function throws cvmexception if vector vDiag passed has size not equal to msize().

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
rvector v(3);
v.set(7.7);
m.set_main_diag(v);
std::cout << m;
prints
(7.7,0) (2,-1) (-1,-2)
(2,1) (7.7,0) (0,-3)
(-1,2) (0,3) (7.7,0)
Parameters
[in]vDiagrvector to set main diagonal to.
Returns
Reference to changed calling matrix.

Definition at line 34965 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::assign_real ( const basic_srsmatrix< TR > &  mRe) throw (cvmexception)
inline

Assigns real symmetric matrix to real part.

Sets real part of calling hermitian matrix to be equal to real symmetric matrix mRe. Function returns reference to the matrix changed. It throws cvmexception in case of different sizes of the operands.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
srsmatrix ms(3);
ms.set(7.);
m.assign_real(ms);
std::cout << m;
prints
(7,0) (7,-1) (7,-2)
(7,1) (7,0) (7,-3)
(7,2) (7,3) (7,0)
Parameters
[in]mResrsmatrix to assign to real part.
Returns
Reference to changed calling matrix.

Definition at line 34998 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::set_real ( TR  d)
inline

Sets all real parts to one value.

Sets real part of every element of calling matrix to be equal to parameter d and returns reference to the matrix changed.

Example:
using namespace cvm;
cmatrix m(2,3);
m.set_real(1.);
std::cout << m;
prints
(1,0) (1,0) (1,0)
(1,0) (1,0) (1,0)
Parameters
[in]dValue to set to.
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35007 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::resize ( tint  nNewDim) throw (cvmexception)
inline

Changes dimension.

Changes dimension of calling square matrix to nNewDim and returns reference to the matrix changed. In case of increasing of its size, the matrix is filled up with zeroes. Function throws cvmexception in case of negative dimension passed or memory allocation failure.

Example:
using namespace cvm;
try {
double a[] = {1., 2., 3., 4., 5., 6., 7., 8.};
scmatrix m((std::complex<double>*) a, 2);
std::cout << m << std::endl;
m.resize (3);
std::cout << m;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(1,2) (5,6)
(3,4) (7,8)
(1,2) (5,6) (0,0)
(3,4) (7,8) (0,0)
(0,0) (0,0) (0,0)
Parameters
[in]nNewDimNew dimension.
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35013 of file cvm.h.

template<typename TR, typename TC>
bool basic_schmatrix< TR, TC >::operator== ( const basic_schmatrix< TR, TC > &  m) const
inline

Matrix comparison.

Operator compares calling hermitian matrix with hermitian matrix m and returns true if they have same dimensions and their appropriate elements differ by not more than cvmMachMin() (the smallest normalized positive number). Returns false otherwise.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m1((std::complex<double>*)a,3);
schmatrix m2(3);
m2.assign((std::complex<double>*)a);
std::cout << (m1 == m2) << std::endl;
prints
1
See Also
operator !=()
Parameters
[in]mhermitian matrix to compare to.
Returns
bool Result of comparison.

Definition at line 35045 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
bool basic_schmatrix< TR, TC >::operator!= ( const basic_schmatrix< TR, TC > &  m) const
inline

Matrix comparison.

Operator compares calling hermitian matrix with hermitian matrix m and returns true if they have different sizes or at least one of their appropriate elements differs by more than cvmMachMin() (the smallest normalized positive number). Returns false otherwise.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m1((std::complex<double>*)a,3);
schmatrix m2(3);
m2.assign((std::complex<double>*)a);
m2.set(2,1,std::complex<double>(2.,1.000001));
std::cout << (m1 != m2) << std::endl;
prints
1
See Also
operator ==()
Parameters
[in]msrsmatrix to compare to.
Returns
bool Result of comparison.

Definition at line 35076 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::operator<< ( const basic_schmatrix< TR, TC > &  m) throw (cvmexception)
inline

Matrix replacement.

Destroys calling hermitian matrix, creates a new one as a copy of hermitian matrix m and returns reference to the matrix changed. Operator throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
schmatrix mc(1);
std::cout << m << std::endl;
std::cout << mc << std::endl;
mc << m;
std::cout << mc;
prints
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
(0,0)
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
See Also
operator =()
Parameters
[in]mschmatrix to replace by.
Returns
Reference to changed calling matrix.

Definition at line 35114 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix basic_schmatrix< TR, TC >::operator+ ( const basic_schmatrix< TR, TC > &  m) const throw (cvmexception)
inline

Addition operator.

Creates object of type schmatrix as sum of calling hermitian matrix and hermitian matrix m. Operator throws cvmexception in case of different sizes of the operands or memory allocation failure.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
double b[] = {1., 0., 1., 1., 1., 1., 1., -1., 1., 0.,
1., 1., 1., -1., 1., -1., 1., 0.};
schmatrix m1((std::complex<double>*)a,3);
schmatrix m2((std::complex<double>*)b,3);
std::cout << m1 + m2;
prints
(2,0) (3,-2) (0,-3)
(3,2) (3,0) (1,-4)
(0,3) (1,4) (4,0)
See Also
sum()
Parameters
[in]mschmatrix to add to calling one.
Returns
Sum of matrices.

Definition at line 35150 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix basic_schmatrix< TR, TC >::operator- ( const basic_schmatrix< TR, TC > &  m) const throw (cvmexception)
inline

Subtraction operator.

Creates object of type schmatrix as difference of calling hermitian matrix and hermitian matrix m. It throws cvmexception in case of different sizes of the operands or memory allocation failure.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
double b[] = {1., 0., 1., 1., 1., 1., 1., -1., 1., 0.,
1., 1., 1., -1., 1., -1., 1., 0.};
schmatrix m1((std::complex<double>*)a,3);
schmatrix m2((std::complex<double>*)b,3);
std::cout << m1 - m2;
prints
(0,0) (1,0) (-2,-1)
(1,0) (1,0) (-1,-2)
(-2,1) (-1,2) (2,0)
See Also
diff()
Parameters
[in]mschmatrix to subtract from calling one.
Returns
Difference of matrices.

Definition at line 35186 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::sum ( const basic_schmatrix< TR, TC > &  m1,
const basic_schmatrix< TR, TC > &  m2 
) throw (cvmexception)
inline

Sum of matrices.

Assigns sum of hermitian matrices m1 and m2 to calling hermitian matrix and returns reference to the matrices changed. It throws cvmexception in case of different sizes of the operands.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
double b[] = {1., 0., 1., 1., 1., 1., 1., -1., 1., 0.,
1., 1., 1., -1., 1., -1., 1., 0.};
schmatrix m1((std::complex<double>*)a,3);
schmatrix m2((std::complex<double>*)b,3);
schmatrix m(3);
std::cout << m.sum(m1, m2);
prints
(2,0) (3,-2) (0,-3)
(3,2) (3,0) (1,-4)
(0,3) (1,4) (4,0)
See Also
operator +()
Parameters
[in]m1First schmatrix summand.
[in]m2Second schmatrix summand.
Returns
Reference to changed calling matrix.

Definition at line 35223 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::diff ( const basic_schmatrix< TR, TC > &  m1,
const basic_schmatrix< TR, TC > &  m2 
) throw (cvmexception)
inline

Difference of matrices.

Assigns difference of hermitian matrices m1 and m2 to calling hermitian matrix and returns reference to the matrix changed. It throws cvmexception in case of different sizes of the operands.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
double b[] = {1., 0., 1., 1., 1., 1., 1., -1., 1., 0.,
1., 1., 1., -1., 1., -1., 1., 0.};
schmatrix m1((std::complex<double>*)a,3);
schmatrix m2((std::complex<double>*)b,3);
schmatrix m(3);
std::cout << m.diff(m1, m2);
prints
(0,0) (1,0) (-2,-1)
(1,0) (1,0) (-1,-2)
(-2,1) (-1,2) (2,0)
See Also
operator -()
Parameters
[in]m1First schmatrix subtrahend.
[in]m2Second schmatrix subtrahend.
Returns
Reference to changed calling matrix.

Definition at line 35260 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::operator+= ( const basic_schmatrix< TR, TC > &  m) throw (cvmexception)
inline

Increment operator.

Adds schmatrix m to calling hermitian matrix and returns reference to the matrix changed. It throws cvmexception in case of different sizes of the operands.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
double b[] = {1., 0., 1., 1., 1., 1., 1., -1., 1., 0.,
1., 1., 1., -1., 1., -1., 1., 0.};
schmatrix m1((std::complex<double>*)a,3);
schmatrix m2((std::complex<double>*)b,3);
m1 += m2;
m2 += m2;
std::cout << m1 << std::endl;
std::cout << m2;
prints
(2,0) (3,-2) (0,-3)
(3,2) (3,0) (1,-4)
(0,3) (1,4) (4,0)
(2,0) (2,-2) (2,-2)
(2,2) (2,0) (2,-2)
(2,2) (2,2) (2,0)
See Also
operator +()
sum()
Parameters
[in]mschmatrix to increment by.
Returns
Reference to changed calling matrix.

Definition at line 35303 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::operator-= ( const basic_schmatrix< TR, TC > &  m) throw (cvmexception)
inline

Decrement operator.

Subtracts schmatrix m from calling hermitian matrix and returns reference to the matrix changed. It throws cvmexception in case of different sizes of the operands.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
double b[] = {1., 0., 1., 1., 1., 1., 1., -1., 1., 0.,
1., 1., 1., -1., 1., -1., 1., 0.};
schmatrix m1((std::complex<double>*)a,3);
schmatrix m2((std::complex<double>*)b,3);
m1 -= m2;
m2 -= m2;
std::cout << m1 << std::endl;
std::cout << m2;
prints
(0,0) (1,0) (-2,-1)
(1,0) (1,0) (-1,-2)
(-2,1) (-1,2) (2,0)
(0,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
See Also
operator -()
diff()
Parameters
[in]mschmatrix to decrement by.
Returns
Reference to changed calling matrix.

Definition at line 35345 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix basic_schmatrix< TR, TC >::operator- ( ) const
inline

Unary minus operator.

Creates object of type scmatrix as calling matrix multiplied by -1. It throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
double a[] = {1., 2., 3., 4., 5., 6., 7., 8.};
scmatrix m((std::complex<double>*) a, 2);
std::cout << -m;
prints
(-1,-2) (-5,-6)
(-3,-4) (-7,-8)
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35352 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::operator++ ( )
inline

Plus identity, prefix.

Adds identity matrix to calling square complex matrix and returns reference to the matrix changed.

Example:
using namespace cvm;
scmatrix m(3);
m.set(std::complex<double>(1.,1.));
m++;
std::cout << m << std::endl;
std::cout << ++m;
prints
(2,1) (1,1) (1,1)
(1,1) (2,1) (1,1)
(1,1) (1,1) (2,1)
(3,1) (1,1) (1,1)
(1,1) (3,1) (1,1)
(1,1) (1,1) (3,1)
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35361 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix basic_schmatrix< TR, TC >::operator++ ( int  )
inline

Plus identity, postfix.

Adds identity matrix to calling square complex matrix and returns reference to the matrix changed.

Example:
using namespace cvm;
scmatrix m(3);
m.set(std::complex<double>(1.,1.));
m++;
std::cout << m << std::endl;
std::cout << ++m;
prints
(2,1) (1,1) (1,1)
(1,1) (2,1) (1,1)
(1,1) (1,1) (2,1)
(3,1) (1,1) (1,1)
(1,1) (3,1) (1,1)
(1,1) (1,1) (3,1)
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35368 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::operator-- ( )
inline

Minus identity, prefix.

Subtracts identity matrix from calling square complex matrix and returns reference to the matrix changed.

Example:
using namespace cvm;
scmatrix m(3);
m.set(std::complex<double>(1.,1.));
m--;
std::cout << m << std::endl;
std::cout << --m;
prints
(0,1) (1,1) (1,1)
(1,1) (0,1) (1,1)
(1,1) (1,1) (0,1)
(-1,1) (1,1) (1,1)
(1,1) (-1,1) (1,1)
(1,1) (1,1) (-1,1)
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35375 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix basic_schmatrix< TR, TC >::operator-- ( int  )
inline

Minus identity, postfix.

Subtracts identity matrix from calling square complex matrix and returns reference to the matrix changed.

Example:
using namespace cvm;
scmatrix m(3);
m.set(std::complex<double>(1.,1.));
m--;
std::cout << m << std::endl;
std::cout << --m;
prints
(0,1) (1,1) (1,1)
(1,1) (0,1) (1,1)
(1,1) (1,1) (0,1)
(-1,1) (1,1) (1,1)
(1,1) (-1,1) (1,1)
(1,1) (1,1) (-1,1)
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35382 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix basic_schmatrix< TR, TC >::operator* ( TR  dMult) const
inline

Multiply by real number operator.

Creates object of type scmatrix as product of calling matrix and real number dMult. It throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
double a[] = {1., 2., 3., 4., 5., 6., 7., 8., 9.,
10., 11., 12., 13., 14., 15., 16., 17., 18.};
scmatrix m((std::complex<double>*) a, 3);
std::cout << m * 5.;
prints
(5,10) (35,40) (65,70)
(15,20) (45,50) (75,80)
(25,30) (55,60) (85,90)
See Also
operator *=()
Parameters
[in]dMultNumber to multiply by.
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35388 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix basic_schmatrix< TR, TC >::operator/ ( TR  dDiv) const throw (cvmexception)
inline

Divide by real number operator.

Creates object of type scmatrix as quotient of calling matrix and real number dDiv. It throws cvmexception if dDiv has absolute value equal or less than cvmMachMin() (the smallest normalized positive number). It also throws exception in case of memory allocation failure.

Example:
using namespace cvm;
double a[] = {1., 2., 3., 4., 5., 6., 7., 8., 9.,
10., 11., 12., 13., 14., 15., 16., 17., 18.};
scmatrix m((std::complex<double>*) a, 3);
std::cout << m / 2.;
prints
(0.5,1) (3.5,4) (6.5,7)
(1.5,2) (4.5,5) (7.5,8)
(2.5,3) (5.5,6) (8.5,9)
See Also
operator /=()
Parameters
[in]dDivNumber to divide by.
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35395 of file cvm.h.

template<typename TR, typename TC>
BaseSCMatrix basic_schmatrix< TR, TC >::operator* ( TC  cMult) const
inline

Multiply by complex number operator.

Creates object of type scmatrix as product of calling matrix and complex number cMult. It throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
double a[] = {1., 2., 3., 4., 5., 6., 7., 8., 9.,
10., 11., 12., 13., 14., 15., 16., 17., 18.};
scmatrix m((std::complex<double>*) a, 3);
std::cout << m * std::complex<double>(1.,1.);
prints
(-1,3) (-1,15) (-1,27)
(-1,7) (-1,19) (-1,31)
(-1,11) (-1,23) (-1,35)
See Also
operator *=()
Parameters
[in]cMultNumber to multiply by.
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35402 of file cvm.h.

template<typename TR, typename TC>
BaseSCMatrix basic_schmatrix< TR, TC >::operator/ ( TC  cDiv) const throw (cvmexception)
inline

Divide by complex number operator.

Creates object of type scmatrix as quotient of calling matrix and complex number cDiv. It throws cvmexception if cDiv has absolute value equal or less than cvmMachMin() (the smallest normalized positive number). It also throws exception in case of memory allocation failure.

Example:
using namespace cvm;
double a[] = {1., 2., 3., 4., 5., 6., 7., 8., 9.,
10., 11., 12., 13., 14., 15., 16., 17., 18.};
scmatrix m((std::complex<double>*) a, 3);
std::cout << m / std::complex<double>(1.,1.);
prints
(1.5,0.5) (7.5,0.5) (13.5,0.5)
(3.5,0.5) (9.5,0.5) (15.5,0.5)
(5.5,0.5) (11.5,0.5) (17.5,0.5)
See Also
operator /=()
Parameters
[in]cDivNumber to divide by.
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35409 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::operator*= ( TR  dMult)
inline

Multiply by real number and assign.

Multiplies calling square complex matrix by real number dMult and returns reference to the matrix changed.

Example:
using namespace cvm;
double a[] = {1., 2., 3., 4., 5., 6., 7., 8., 9.,
10., 11., 12., 13., 14., 15., 16., 17., 18.};
scmatrix m((std::complex<double>*) a, 3);
m *= 5.;
std::cout << m;
prints
(5,10) (35,40) (65,70)
(15,20) (45,50) (75,80)
(25,30) (55,60) (85,90)
See Also
operator *(TR) const
Parameters
[in]dMultNumber to multiply by.
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35416 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::operator/= ( TR  dDiv) throw (cvmexception)
inline

Divide by real number and assign.

Divides calling square complex matrix by real number dDiv and returns reference to the matrix changed. It throws cvmexception if dDiv has absolute value equal or less than cvmMachMin() (the smallest normalized positive number).

Example:
using namespace cvm;
double a[] = {1., 2., 3., 4., 5., 6., 7., 8., 9.,
10., 11., 12., 13., 14., 15., 16., 17., 18.};
scmatrix m((std::complex<double>*) a, 3);
m /= 2.;
std::cout << m;
prints
(0.5,1) (3.5,4) (6.5,7)
(1.5,2) (4.5,5) (7.5,8)
(2.5,3) (5.5,6) (8.5,9)
See Also
operator /(TR) const
Parameters
[in]dDivNumber to divide by.
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35422 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::normalize ( )
inline

Matrix normalizer.

Normalizes calling matrix so its Euclidean norm() becomes equal to 1 if it was greater than cvmMachMin() (the smallest normalized positive number) before the call. Does nothing otherwise.

Example:
using namespace cvm;
double a[] = {1., 2., 3., 4., 5., 6.,
7., 8., 9., 10., 11., 12.};
cmatrix ma ((std::complex<double>*) a, 2, 3);
ma.normalize();
std::cout << ma << ma.norm() << std::endl;
prints
(0.0392232,0.0784465) (0.196116,0.235339) (0.353009,0.392232)
(0.11767,0.156893) (0.274563,0.313786) (0.431455,0.470679)
1
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35429 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix basic_schmatrix< TR, TC >::operator! ( ) const throw (cvmexception)
inline

Matrix transposition.

Creates object of type schmatrix as transposed calling hermitian matrix.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
schmatrix mt(3);
std::cout << m << std::endl << !m << std::endl ;
mt.transpose(m);
std::cout << mt << std::endl;
mt.transpose();
std::cout << mt;
prints
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
(1,0) (2,1) (-1,2)
(2,-1) (2,0) (0,3)
(-1,-2) (0,-3) (3,0)
(1,0) (2,1) (-1,2)
(2,-1) (2,0) (0,3)
(-1,-2) (0,-3) (3,0)
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35473 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix basic_schmatrix< TR, TC >::operator~ ( ) const
inline

Does nothing and returns copy of calling hermitian matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35482 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::transpose ( const basic_schmatrix< TR, TC > &  m) throw (cvmexception)
inline

Matrix transposition.

Sets calling hermitian matrix to be equal to hermitian matrix m transposed. Function throws cvmexception in case of not appropriate sizes of the operands.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
schmatrix mt(3);
std::cout << m << std::endl << !m << std::endl ;
mt.transpose(m);
std::cout << mt << std::endl;
mt.transpose();
std::cout << mt;
prints
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
(1,0) (2,1) (-1,2)
(2,-1) (2,0) (0,3)
(-1,-2) (0,-3) (3,0)
(1,0) (2,1) (-1,2)
(2,-1) (2,0) (0,3)
(-1,-2) (0,-3) (3,0)
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
Parameters
[in]mschmatrix to transpose.
Returns
Reference to changed calling matrix.

Definition at line 35525 of file cvm.h.

Here is the caller graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::conj ( const basic_schmatrix< TR, TC > &  m) throw (cvmexception)
inline

Assigns hermitian matrix m to calling one and returns reference to the matrix changed.

Definition at line 35534 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::transpose ( ) throw (cvmexception)
inline

Matrix transposition (in-place)

Makes calling hermitian matrix to be equal to transposed itself. Function throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
schmatrix mt(3);
std::cout << m << std::endl << !m << std::endl ;
mt.transpose(m);
std::cout << mt << std::endl;
mt.transpose();
std::cout << mt;
prints
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
(1,0) (2,1) (-1,2)
(2,-1) (2,0) (0,3)
(-1,-2) (0,-3) (3,0)
(1,0) (2,1) (-1,2)
(2,-1) (2,0) (0,3)
(-1,-2) (0,-3) (3,0)
(1,0) (2,-1) (-1,-2)
(2,1) (2,0) (0,-3)
(-1,2) (0,3) (3,0)
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35579 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::conj ( )
inline

Does nothing and returns reference to calling symmetric matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35588 of file cvm.h.

template<typename TR, typename TC>
CVector basic_schmatrix< TR, TC >::operator* ( const CVector v) const throw (cvmexception)
inline

Matrix-vector product.

Creates object of type cvector as product of calling matrix and vector v. Function throws cvmexception if the number of columns of the calling matrix differs from size of the vector v. Use basic_cvector::mult(const basic_cmatrix<TR,TC>&,const basic_cvector<TR,TC>&) to avoid new object creation.

Example:
using namespace cvm;
cmatrix m(2, 3);
cvector v(3);
m.set(std::complex<double>(1.,1.));
v.set(std::complex<double>(1.,1.));
std::cout << m * v;
prints
(0,6) (0,6)
See Also
basic_cvector::mult(const basic_cmatrix<TR,TC>&,const basic_cvector<TR,TC>&)
Parameters
[in]vcvector to compute product with.
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35592 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
BaseCMatrix basic_schmatrix< TR, TC >::operator* ( const BaseCMatrix m) const throw (cvmexception)
inline

Matrix-matrix product.

Creates object of type cmatrix as product of calling square complex matrix and matrix m. Operator throws cvmexception if number of columns of calling matrix differs from number of rows of the matrix m. Use mult() to avoid new object creation.

Example:
using namespace cvm;
try {
scmatrix ms(3);
cmatrix m(3,2);
ms.set(std::complex<double>(1.,1.));
m.set(std::complex<double>(1.,1.));
std::cout << ms * m;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(0,6) (0,6)
(0,6) (0,6)
(0,6) (0,6)
Parameters
[in]mcmatrix to compute product with.
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35601 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
BaseSCMatrix basic_schmatrix< TR, TC >::operator* ( const BaseSCMatrix m) const throw (cvmexception)
inline

Matrix-matrix product.

Creates object of type scmatrix as product of calling hermitian matrix and square complex matrix m. Operator throws cvmexception if dimension of calling matrix differs from dimension the matrix m. Use mult() to avoid new object creation.

Example:
using namespace cvm;
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
try {
schmatrix ms((std::complex<double>*)a,3);
scmatrix m(3);
m.set(std::complex<double>(1.,1.));
std::cout << ms * m;
}
catch (std::exception& e) {
std::cout << "Exception: " << e.what () << std::endl;
}
prints
(5,-1) (5,-1) (5,-1)
(6,2) (6,2) (6,2)
(-3,7) (-3,7) (-3,7)
Parameters
[in]mscmatrix to compute product with.
Returns
Result object.

Definition at line 35643 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
CVector basic_schmatrix< TR, TC >::operator/ ( const CVector vB) const throw (cvmexception)
inline

Linear solver operator.

Creates cvector object as solution $x$ of linear equation $A*x=b$ where calling matrix is square complex matrix $A$ and parameter vB is vector $b$. Operator throws cvmexception in case of inappropriate sizes of the operands or when matrix $A$ is close to singular.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::showpos);
std::cout.precision (7);
try {
double m[] = {1., -1., 1., 2., -2., 1., 3., -3.};
double b[] = {1., 2., 5., -3.};
scmatrix ma((std::complex<double>*) m, 2);
cvector vb((std::complex<double>*) b, 2);
cvector vx(2);
vx = ma / vb;
std::cout << ma * vx - vb;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(-6.6613381e-016,+4.4408921e-016) (+0.0000000e+000,+0.0000000e+000)
Parameters
[in]vBcvector $b$.
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 35655 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::herk ( TR  alpha,
const CVector v,
TR  beta 
) throw (cvmexception)
inline

Rank-1 update matrix-vector operation.

Calls one of ZHERK routines of the BLAS library performing rank-1 update matrix-vector operation defined as

\[ A=\alpha\,\begin{pmatrix} v_1 \\ v_2 \\ \vdots \\ v_n \end{pmatrix} \begin{pmatrix} v_1^* & v_2^* & \cdots & v_n^* \end{pmatrix} + \beta A, \]

where $\alpha$ and $\beta$ are real numbers (parameters alpha and beta), $A$ is calling hermitian matrix and $v$ is complex vector (parameter v). Function returns reference to the matrix changed and throws cvmexception in case of inappropriate sizes of the operands.

Example:
using namespace cvm;
const treal a[] = {1., -1., 2., 2., 3., -3.};
const cvector v((std::complex<double>*)a, 3);
const treal alpha = 2.12;
const treal beta = -3.07;
schmatrix mh(3), mh2(3);
mh.randomize_real(-1.,2.);
mh.randomize_imag(-2.,3.);
mh2 = mh;
mh.herk (alpha, v, beta);
mh2 = alpha * schmatrix(v.rank1update_c(v)) + beta * mh2;
std::cout << mh - mh2;
prints
(0,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
See Also
http://www.netlib.org/blas
Parameters
[in]alphaMultiplier $\alpha$.
[in]vcvector $v$.
[in]betaMultiplier $\beta$.
Returns
Reference to changed calling matrix.

Definition at line 35707 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::herk ( bool  bTransp,
TR  alpha,
const BaseCMatrix m,
TR  beta 
) throw (cvmexception)
inline

Rank-1 update matrix-matrix operation.

Calls ZHERK routine of the BLAS library performing rank-1 update matrix-vector operation defined as

\[ A=\alpha MM^H + \beta A\quad\text{or}\quad A=\alpha M^HM + \beta A. \]

where $\alpha$ and $\beta$ are real numbers (parameters alpha and beta), $A$ is calling hermitian matrix and $M$ is complex matrix (parameter m). First operation is performed if bTransp passed is false and second one otherwise. Function returns reference to the matrix changed and throws cvmexception in case of inappropriate sizes of the operands.

Example 1:
using namespace cvm;
cmatrix m(3,3);
const treal alpha = 2.12;
const treal beta = -3.07;
schmatrix mh(3), mh2(3);
m.randomize_real(-1.,2.);
m.randomize_imag(-2.,3.);
mh.randomize_real(-1.,2.);
mh.randomize_imag(-2.,3.);
mh2 = mh;
mh.herk (false, alpha, m, beta);
mh2 = alpha * schmatrix (m * ~m) + beta * mh2;
std::cout << mh - mh2;
prints
(0,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
Example 2:
using namespace cvm;
cmatrix m(3,3);
const treal alpha = 2.12;
const treal beta = -3.07;
schmatrix mh(3), mh2(3);
m.randomize_real(-1.,2.);
m.randomize_imag(-2.,3.);
mh.randomize_real(-1.,2.);
mh.randomize_imag(-2.,3.);
mh2 = mh;
mh.herk (true, alpha, m, beta);
mh2 = alpha * schmatrix (~m * m) + beta * mh2;
std::cout << mh - mh2;
prints
(3.55271e-015,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
(0,0) (0,0) (1.77636e-015,0)
See Also
http://www.netlib.org/blas
Parameters
[in]bTranspUse transposed operation.
[in]alphaMultiplier $\alpha$.
[in]mcmatrix $M$.
[in]betaMultiplier $\beta$.
Returns
Reference to changed calling matrix.

Definition at line 35781 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::her2k ( TC  alpha,
const CVector v1,
const CVector v2,
TR  beta 
) throw (cvmexception)
inline

Rank-1 update matrix-vector operation.

Calls ZHER2K routine of the BLAS library performing rank-1 update matrix-vector operation defined as

\[ A=\alpha v_1 v_2^* + \alpha^* v_2 v_1^* + \beta A, \]

where $\alpha$ is complex and $\beta$ is real number (parameters alpha and beta), $A$ is calling hermitian matrix and $v_1$ and $v_2$ are complex vectors (parameters v1 and v2). Function returns reference to the matrix changed and throws cvmexception in case of inappropriate sizes of the operands.

Example:
using namespace cvm;
std::cout.setf(std::ios::scientific | std::ios::showpos);
std::cout.precision(2);
const tcomplex alpha(2.12,-0.14);
const tcomplex alphac(2.12,0.14);
const treal beta = -3.07;
cvector v1(3), v2(3);
schmatrix mh(3), mh2(3);
v1.randomize_real(-1.,2.);
v1.randomize_imag(-2.,3.);
v2.randomize_real(-1.,2.);
v2.randomize_imag(-2.,3.);
mh.randomize_real(-1.,2.);
mh.randomize_imag(-2.,3.);
mh2 = mh;
mh.her2k (alpha, v1, v2, beta);
mh2 = schmatrix(alpha * v1.rank1update_c(v2) + alphac * v2.rank1update_c(v1)) + beta * mh2;
std::cout << mh - mh2;
prints
(+1.78e-015,+0.00e+000) (+0.00e+000,+0.00e+000) (+0.00e+000,+0.00e+000)
(+0.00e+000,+0.00e+000) (-1.78e-015,+0.00e+000) (+0.00e+000,+1.78e-015)
(+0.00e+000,+0.00e+000) (+0.00e+000,-1.78e-015) (+4.44e-016,+0.00e+000)
See Also
http://www.netlib.org/blas
Parameters
[in]alphaComplex multiplier $\alpha$.
[in]v1cvector $v_1$.
[in]v2cvector $v_2$.
[in]betaReal multiplier $\beta$.
Returns
Reference to changed calling matrix.

Definition at line 35836 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::her2k ( bool  bTransp,
TC  alpha,
const BaseCMatrix m1,
const BaseCMatrix m2,
TR  beta 
) throw (cvmexception)
inline

Rank-1 update matrix-matrix operation.

Calls ZHE2K routine of the BLAS library performing rank-1 update matrix-matrix operation defined as

\[ A=\alpha M_1 M_2^H + \alpha^* M_2 M_1^H + \beta A\quad \text{or}\quad A=\alpha M_1^H M_2 + \alpha^* M_2^H M_1 + \beta A. \]

where $\alpha$ is complex and $\beta$ is real number (parameters alpha and beta), $A$ is calling hermitian matrix and $M_1$ and $M_2$ are complex matrices (parameters m1 and m2). First operation is performed if bTransp passed is false and second one otherwise. Function returns reference to the matrix changed and throws cvmexception in case of inappropriate sizes of the operands.

Example:
using namespace cvm;
std::cout.setf(std::ios::scientific | std::ios::showpos);
std::cout.precision(2);
const tcomplex alpha(2.12,-0.14);
const tcomplex alphac(2.12,0.14);
const treal beta = -3.07;
cmatrix m1(3,3), m2(3,3);
schmatrix mh(3), mh2(3);
m1.randomize_real(-1.,2.);
m1.randomize_imag(-2.,3.);
m2.randomize_real(-1.,2.);
m2.randomize_imag(-2.,3.);
mh.randomize_real(-1.,2.);
mh.randomize_imag(-2.,3.);
mh2 = mh;
mh.her2k (false, alpha, m1, m2, beta);
mh2 = schmatrix(alpha * m1 * ~m2 + alphac * m2 * ~m1, 1.e-14) + beta * mh2;
std::cout << mh - mh2;
prints
(+0.00e+000,+4.44e-016) (+0.00e+000,+0.00e+000) (+0.00e+000,+1.78e-015)
(+3.55e-015,+0.00e+000) (+4.44e-016,-8.88e-016) (+0.00e+000,+0.00e+000)
(+8.88e-016,+8.88e-016) (+0.00e+000,+0.00e+000) (+8.88e-016,+1.78e-015)
Example:
using namespace cvm;
std::cout.setf(std::ios::scientific | std::ios::showpos);
std::cout.precision(2);
const tcomplex alpha(2.12,-0.14);
const tcomplex alphac(2.12,0.14);
const treal beta = -3.07;
cmatrix m1(3,3), m2(3,3);
schmatrix mh(3), mh2(3);
m1.randomize_real(-1.,2.);
m1.randomize_imag(-2.,3.);
m2.randomize_real(-1.,2.);
m2.randomize_imag(-2.,3.);
mh.randomize_real(-1.,2.);
mh.randomize_imag(-2.,3.);
mh2 = mh;
mh.her2k (true, alpha, m1, m2, beta);
mh2 = schmatrix(alpha * ~m1 * m2 + alphac * ~m2 * m1, 1.e-14) + beta * mh2;
std::cout << mh - mh2;
prints
(+3.55e-015,-8.88e-016) (+0.00e+000,+0.00e+000) (+0.00e+000,-1.78e-015)
(-8.88e-016,+1.78e-015) (+0.00e+000,-4.44e-016) (+8.88e-016,+0.00e+000)
(+0.00e+000,+3.55e-015) (-1.33e-015,-1.78e-015) (+0.00e+000,+0.00e+000)
See Also
http://www.netlib.org/blas
Parameters
[in]bTranspUse transposed version of the operation.
[in]alphaComplex multiplier $\alpha$.
[in]m1cmatrix $M_1$.
[in]m2cmatrix $M_2$.
[in]betaReal multiplier $\beta$.
Returns
Reference to changed calling matrix.

Definition at line 35925 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::inv ( const basic_schmatrix< TR, TC > &  m) throw (cvmexception)
inline

Matrix inversion.

This version sets calling matrix to be equal to hermitian matrix m inverted. Function throws cvmexception in case of inappropriate sizes of the operands or when the matrix to be inverted is close to singular.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
schmatrix mi(3);
mi.inv(m);
std::cout << mi << std::endl;
std::cout << mi * m - eye_complex(3);
prints
(-1.50e+000,0.00e+000) (-1.67e-016,8.33e-017) (-5.00e-001,-1.00e+000)
(-1.67e-016,-8.33e-017) (-1.00e+000,0.00e+000) (0.00e+000,-1.00e+000)
(-5.00e-001,1.00e+000) (0.00e+000,1.00e+000) (-1.50e+000,0.00e+000)
(2.22e-016,0.00e+000) (1.11e-016,0.00e+000) (2.78e-017,0.00e+000)
(2.78e-016,2.22e-016) (4.44e-016,0.00e+000) (0.00e+000,-4.44e-016)
(2.22e-016,-4.44e-016) (0.00e+000,-8.88e-016) (-2.22e-016,0.00e+000)
Parameters
[in]mschmatrix to invert.
Returns
Reference to changed calling matrix.

Definition at line 35973 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix basic_schmatrix< TR, TC >::inv ( ) const throw (cvmexception)
inline

Matrix inversion.

This version creates schmatrix object equal to calling hermitian matrix inverted. Function throws cvmexception in case of memory allocation failure or when the matrix to be inverted is close to singular.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
const schmatrix mi = m.inv();
std::cout << mi << std::endl;
std::cout << mi * m - eye_complex(3);
prints
(-1.50e+000,0.00e+000) (-1.67e-016,8.33e-017) (-5.00e-001,-1.00e+000)
(-1.67e-016,-8.33e-017) (-1.00e+000,0.00e+000) (0.00e+000,-1.00e+000)
(-5.00e-001,1.00e+000) (0.00e+000,1.00e+000) (-1.50e+000,0.00e+000)
(2.22e-016,0.00e+000) (1.11e-016,0.00e+000) (2.78e-017,0.00e+000)
(2.78e-016,2.22e-016) (4.44e-016,0.00e+000) (0.00e+000,-4.44e-016)
(2.22e-016,-4.44e-016) (0.00e+000,-8.88e-016) (-2.22e-016,0.00e+000)
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 36009 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::exp ( const basic_schmatrix< TR, TC > &  m,
TR  tol = basic_cvmMachSp<TR>() 
) throw (cvmexception)
inline

Matrix exponent.

Computes exponent of hermitian matrix using Pade approximation defined as

\[ R_{pq}(z)=D_{pq}(z)^{-1}N_{pq}(z)=1+z+\dots+z^p/p!\,, \]

where

\[\begin{aligned} N_{pq}(z)&=\sum_{k=0}^p\frac{(p+q-k)!p!}{(p+q)!k!(q-k)!}z^k,\\ D_{pq}(z)&=\sum_{k=0}^q\frac{(p+q-k)!p!}{(p+q)!k!(q-k)!}(-z)^k \end{aligned}\]

along with the matrix normalizing as described in Gene H. Golub, Charles F. Van Loan. Matrix Computations, The Johns Hopkins University Press, 1996, 694~p., ISBN 0-8018-5413-X, p. 572. Function uses DMEXP (or SMEXP for float version) FORTRAN subroutine implementing the algorithm. This version sets calling matrix to be equal to exponent of hermitian matrix m and returns reference to the matrix changed. The algorithm uses parameter tol as $\varepsilon(p,q)$ in order to choose constants $p$ and $q$ so that

\[ \varepsilon(p,q)\ge 2^{3-(p+q)}\frac{p!q!}{(p+q)!(p+q+1)!}. \]

This parameter is equal to the largest relative spacing (see cvmMachSp() ) by default. Function throws cvmexception in case of inappropriate sizes of the operands or when LAPACK subroutine fails.

Example
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left |
std::ios::showpos);
std::cout.precision (15);
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
schmatrix me(3);
me.exp(m);
std::cout << "Column 1" << std::endl
<< me(1,1) << std::endl << me(2,1) << std::endl << me(3,1) << std::endl
<< "Column 2" << std::endl
<< me(1,2) << std::endl << me(2,2) << std::endl << me(3,2) << std::endl
<< "Column 3" << std::endl
<< me(1,3) << std::endl << me(2,3) << std::endl << me(3,3) << std::endl;
prints
Column 1
(+2.673228708372002e+002,+1.091141066389412e-014)
(+3.071187567026803e+002,+1.535593783513402e+002)
(-1.749365628720766e+002,+3.498731257441531e+002)
Column 2
(+3.071187567026803e+002,-1.535593783513401e+002)
(+4.422594337092766e+002,+1.919736460939478e-015)
(-9.600094996571151e-015,+5.034325040954932e+002)
Column 3
(-1.749365628720765e+002,-3.498731257441531e+002)
(+6.184072406183948e-015,-5.034325040954932e+002)
(+5.744416275398805e+002,+1.540673883337074e-014)
Matlab output:
Column 1
2.673228708371998e+002 -7.105427357601002e-015i
3.071187567026802e+002 +1.535593783513401e+002i
-1.749365628720764e+002 +3.498731257441527e+002i
Column 2
3.071187567026802e+002 -1.535593783513401e+002i
4.422594337092769e+002 -5.489286670342458e-016i
3.549798266275454e-015 +5.034325040954932e+002i
Column 3
-1.749365628720763e+002 -3.498731257441526e+002i
-1.776065298147746e-014 -5.034325040954931e+002i
5.744416275398801e+002 -2.096383162906490e-014i
Parameters
[in]mschmatrix to compute exponent for.
[in]tolComputation tolerance.
Returns
Reference to changed calling matrix.

Definition at line 36100 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix basic_schmatrix< TR, TC >::exp ( TR  tol = basic_cvmMachSp<TR>()) const throw (cvmexception)
inline

Matrix exponent.

Computes exponent of hermitian matrix using Pade approximation defined as

\[ R_{pq}(z)=D_{pq}(z)^{-1}N_{pq}(z)=1+z+\dots+z^p/p!\,, \]

where

\[\begin{aligned} N_{pq}(z)&=\sum_{k=0}^p\frac{(p+q-k)!p!}{(p+q)!k!(q-k)!}z^k,\\ D_{pq}(z)&=\sum_{k=0}^q\frac{(p+q-k)!p!}{(p+q)!k!(q-k)!}(-z)^k \end{aligned}\]

along with the matrix normalizing as described in Gene H. Golub, Charles F. Van Loan. Matrix Computations, The Johns Hopkins University Press, 1996, 694~p., ISBN 0-8018-5413-X, p. 572. Function uses DMEXP (or SMEXP for float version) FORTRAN subroutine implementing the algorithm. This version creates object of type schmatrix as exponent of calling hermitian matrix. The algorithm uses parameter tol as $\varepsilon(p,q)$ in order to choose constants $p$ and $q$ so that

\[ \varepsilon(p,q)\ge 2^{3-(p+q)}\frac{p!q!}{(p+q)!(p+q+1)!}. \]

This parameter is equal to the largest relative spacing (see cvmMachSp() ) by default. Function throws cvmexception in case of inappropriate sizes of the operands or when LAPACK subroutine fails.

Example
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left |
std::ios::showpos);
std::cout.precision (15);
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
schmatrix me = m.exp();
std::cout << "Column 1" << std::endl
<< me(1,1) << std::endl << me(2,1) << std::endl << me(3,1) << std::endl
<< "Column 2" << std::endl
<< me(1,2) << std::endl << me(2,2) << std::endl << me(3,2) << std::endl
<< "Column 3" << std::endl
<< me(1,3) << std::endl << me(2,3) << std::endl << me(3,3) << std::endl;
prints
Column 1
(+2.673228708372002e+002,+1.091141066389412e-014)
(+3.071187567026803e+002,+1.535593783513402e+002)
(-1.749365628720766e+002,+3.498731257441531e+002)
Column 2
(+3.071187567026803e+002,-1.535593783513401e+002)
(+4.422594337092766e+002,+1.919736460939478e-015)
(-9.600094996571151e-015,+5.034325040954932e+002)
Column 3
(-1.749365628720765e+002,-3.498731257441531e+002)
(+6.184072406183948e-015,-5.034325040954932e+002)
(+5.744416275398805e+002,+1.540673883337074e-014)
Matlab output:
Column 1
2.673228708371998e+002 -7.105427357601002e-015i
3.071187567026802e+002 +1.535593783513401e+002i
-1.749365628720764e+002 +3.498731257441527e+002i
Column 2
3.071187567026802e+002 -1.535593783513401e+002i
4.422594337092769e+002 -5.489286670342458e-016i
3.549798266275454e-015 +5.034325040954932e+002i
Column 3
-1.749365628720763e+002 -3.498731257441526e+002i
-1.776065298147746e-014 -5.034325040954931e+002i
5.744416275398801e+002 -2.096383162906490e-014i
Parameters
[in]tolComputation tolerance.
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 36187 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::polynom ( const basic_schmatrix< TR, TC > &  m,
const RVector v 
) throw (cvmexception)
inline

Matrix polynomial.

Computes hermitian matrix polynomial defined as

\[ p(A)=b_0I+b_1A+\dots+b_qA^q \]

using the Horner's rule:

\[ \newcommand{\floor}{\operatorname{floor}}p(A)=\sum_{k=0}^r B_k(A^s)^k,\quad s=\floor(\!\sqrt{q}\,),\ r=\floor(q/s) \]

where

\[ B_k=\begin{cases} \sum\limits_{i=0}^{s-1}b_{sk+i} A^i, & k=0,1,\dots,r-1\\ \sum\limits_{i=0}^{q-sr}b_{sr+i} A^i, & k=r. \end{cases} \]

See also Gene H. Golub, Charles F. Van Loan. Matrix Computations, The Johns Hopkins University Press, 1996, 694~p., ISBN 0-8018-5413-X, p. 568. Real coefficients $b_0,b_1,\dots,b_q$ are passed in parameter v, where $q$ is equal to v.size()-1, so the function computes matrix polynomial equal to

\[ \mathtt{v[1]*I + v[2]*m +\cdots + v[v.size()]*m^{v.size()-1}} \]

This version sets calling matrix to be equal to the polynomial of hermitian matrix m. Function uses ZPOLY (or CPOLY for float version) FORTRAN subroutine implementing the Horner's algorithm. Function throws cvmexception in case of inappropriate sizes of the operands.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left |
std::ios::showpos);
std::cout.precision (10);
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
double re[]={2.2,1.3,1.1,-0.9,0.2,-0.45,45.,-30.,10.,3.,1.13};
const rvector vr(re, 11);
mp.polynom (m, vr);
std::cout << "Column 1" << std::endl
<< mp(1,1) << std::endl << mp(2,1) << std::endl << mp(3,1) << std::endl
<< "Column 2" << std::endl
<< mp(1,2) << std::endl << mp(2,2) << std::endl << mp(3,2) << std::endl
<< "Column 3" << std::endl
<< mp(1,3) << std::endl << mp(2,3) << std::endl << mp(3,3) << std::endl;
prints
Column 1
(+1.2319548758e+008,+0.0000000000e+000)
(+1.4179323916e+008,+7.0896619580e+007)
(-8.0802738460e+007,+1.6160547692e+008)
Column 2
(+1.4179323916e+008,-7.0896619580e+007)
(+2.0399822604e+008,+0.0000000000e+000)
(+0.0000000000e+000,+2.3250209650e+008)
Column 3
(-8.0802738460e+007,-1.6160547692e+008)
(+0.0000000000e+000,-2.3250209650e+008)
(+2.6498872674e+008,+0.0000000000e+000)
Matlab output:
Column 1
1.231954875800000e+008
1.417932391600000e+008 +7.089661958000000e+007i
-8.080273845999999e+007 +1.616054769200000e+008i
Column 2
1.417932391600000e+008 -7.089661958000000e+007i
2.039982260400000e+008
0 +2.325020965000000e+008i
Column 3
-8.080273845999999e+007 -1.616054769200000e+008i
0 -2.325020965000000e+008i
2.649887267400000e+008
Parameters
[in]mschmatrix to compute polynomial for.
[in]vReal vector of coefficients.
Returns
Reference to changed calling matrix.

Definition at line 36284 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_schmatrix basic_schmatrix< TR, TC >::polynom ( const RVector v) const
inline

Matrix polynomial.

Computes hemitian matrix polynomial defined as

\[ p(A)=b_0I+b_1A+\dots+b_qA^q \]

using the Horner's rule:

\[ \newcommand{\floor}{\operatorname{floor}}p(A)=\sum_{k=0}^r B_k(A^s)^k,\quad s=\floor(\!\sqrt{q}\,),\ r=\floor(q/s) \]

where

\[ B_k=\begin{cases} \sum\limits_{i=0}^{s-1}b_{sk+i} A^i, & k=0,1,\dots,r-1\\ \sum\limits_{i=0}^{q-sr}b_{sr+i} A^i, & k=r. \end{cases} \]

See also Gene H. Golub, Charles F. Van Loan. Matrix Computations, The Johns Hopkins University Press, 1996, 694~p., ISBN 0-8018-5413-X, p. 568. The coefficients $b_0,b_1,\dots,b_q$ are passed in parameter v, where $q$ is equal to v.size()-1, so the function computes matrix polynomial equal to

\[ \mathtt{v[1]*I + v[2]*m +\cdots + v[v.size()]*m^{v.size()-1}} \]

This version creates object of type schmatrix as the polynomial of hemitian symmetric matrix. Function uses ZPOLY (or CPOLY for float version) FORTRAN subroutine implementing the Horner's algorithm. Function throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left |
std::ios::showpos);
std::cout.precision (10);
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
double re[]={2.2,1.3,1.1,-0.9,0.2,-0.45,45.,-30.,10.,3.,1.13};
const rvector vr(re, 11);
schmatrix mp = m.polynom(vr);
std::cout << "Column 1" << std::endl
<< mp(1,1) << std::endl << mp(2,1) << std::endl << mp(3,1) << std::endl
<< "Column 2" << std::endl
<< mp(1,2) << std::endl << mp(2,2) << std::endl << mp(3,2) << std::endl
<< "Column 3" << std::endl
<< mp(1,3) << std::endl << mp(2,3) << std::endl << mp(3,3) << std::endl;
prints
Column 1
(+1.2319548758e+008,+0.0000000000e+000)
(+1.4179323916e+008,+7.0896619580e+007)
(-8.0802738460e+007,+1.6160547692e+008)
Column 2
(+1.4179323916e+008,-7.0896619580e+007)
(+2.0399822604e+008,+0.0000000000e+000)
(+0.0000000000e+000,+2.3250209650e+008)
Column 3
(-8.0802738460e+007,-1.6160547692e+008)
(+0.0000000000e+000,-2.3250209650e+008)
(+2.6498872674e+008,+0.0000000000e+000)
Matlab output:
Column 1
1.231954875800000e+008
1.417932391600000e+008 +7.089661958000000e+007i
-8.080273845999999e+007 +1.616054769200000e+008i
Column 2
1.417932391600000e+008 -7.089661958000000e+007i
2.039982260400000e+008
0 +2.325020965000000e+008i
Column 3
-8.080273845999999e+007 -1.616054769200000e+008i
0 -2.325020965000000e+008i
2.649887267400000e+008
Parameters
[in]vReal vector of coefficients.
Returns
Result object.

Definition at line 36380 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
RVector basic_schmatrix< TR, TC >::eig ( BaseSCMatrix mEigVect) const throw (cvmexception)
inline

Eigenvalues and eigenvectors.

Solves symmetric eigenvalue problem and creates rvector object equal to eigenvalues of calling hermitian matrix. The symmetric eigenvalue problem is defined as follows: given symmetric (hermitian) matrix $A$, find the eigenvalues $\lambda$ and the corresponding eigenvectors $z$ that satisfy the equation

\[ Az = \lambda z. \]

Eigenvalues of hermitian complex matrix $A$ are real and eigenvectors are complex. Function sets output parameter mEigVect to be equal to square matrix containing eigenvectors as columns. Function throws cvmexception in case of inappropriate calling object sizes or in case of convergence error.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
scmatrix me(3);
rvector v(3);
v = m.eig(me);
std::cout << v << std::endl;
cvector vc(v);
std::cout << m * me(1) - me(1) * vc(1);
std::cout << m * me(2) - me(2) * vc(2);
std::cout << m * me(3) - me(3) * vc(3);
// orthogonality check:
std::cout << me(1) % me(2) << std::endl;
std::cout << me(2) % me(3) << std::endl;
std::cout << me(1) % me(3) << std::endl;
prints
-8.13e-001 -3.44e-001 7.16e+000
(1.39e-016,2.22e-016) (5.25e-017,-1.11e-016) (1.94e-016,1.67e-016)
(4.86e-016,4.44e-016) (7.63e-017,0.00e+000) (3.33e-016,2.22e-016)
(-2.22e-016,-8.88e-016) (-5.55e-017,-8.88e-016) (8.88e-016,-5.55e-017)
(-5.17e-017,-9.74e-017)
(-5.81e-017,-5.40e-017)
(2.37e-018,-3.56e-017)
See Also
rvector::eig()
Parameters
[out]mEigVectEigenvectors of calling symmetric matrix.
Returns
Result object.

Definition at line 36441 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
RVector basic_schmatrix< TR, TC >::eig ( ) const throw (cvmexception)
inline

Eigenvalues.

Solves symmetric eigenvalue problem and creates rvector object equal to eigenvalues of calling hermitian matrix. The symmetric eigenvalue problem is defined as follows: given symmetric (hermitian) matrix $A$, find the eigenvalues $\lambda$ and the corresponding eigenvectors $z$ that satisfy the equation

\[ Az = \lambda z. \]

Eigenvalues of hermitian complex matrix $A$ are real. Function throws cvmexception in case of memory allocation failure or convergence error.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 0., 2., 1., -1., 2., 2., -1., 2., 0.,
0., 3., -1., -2., 0., -3., 3., 0.};
schmatrix m((std::complex<double>*)a,3);
std::cout << m.eig();
prints
-8.13e-001 -3.44e-001 7.16e+000
See Also
rvector::eig()
Returns
Result object.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 36480 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
BaseSCMatrix basic_schmatrix< TR, TC >::cholesky ( ) const throw (cvmexception)
inline

Cholesky factorization.

Forms the Cholesky factorization of hermitian positive-definite calling matrix $A$ defined as

\[ A=U^H U, \]

where $U$ is upper triangular matrix. It utilizes one of ZPOTRF routines of the LAPACK library. Function creates scmatrix object equal to the factorization. Function throws cvmexception in case of convergence error.

Example:
using namespace cvm;
try {
double a[] = {3., 0., 2., 1., -1., 2., 2., -1., 3., 0.,
0., 3., -1., -2., 0., -3., 5., 0.};
const schmatrix m((std::complex<double>*)a,3);
scmatrix h = m.cholesky();
std::cout << h << std::endl;
std::cout << ~h * h - m;
}
catch (std::exception& e) {
std::cout << "Exception: " << e.what () << std::endl;
}
prints
(1.73205,0) (1.1547,-0.57735) (-0.57735,-1.1547)
(0,0) (1.1547,0) (0,-1.1547)
(0,0) (0,0) (1.41421,0)
(-4.44089e-016,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
See Also
http://www.netlib.org/lapack
basic_scmatrix::cholesky()
Returns
Result object.

Definition at line 36528 of file cvm.h.

template<typename TR, typename TC>
BaseSCMatrix basic_schmatrix< TR, TC >::bunch_kaufman ( tint nPivots) const throw (cvmexception)
inline

Bunch-Kaufman factorization.

Forms the Bunch-Kaufman factorization of hermitian calling matrix (cited from the MKL library documentation):

\[ A=PUDU^HP^H, \]

where $A$ is calling hermitian matrix, $P$ is a permutation matrix, $U$ and $L$ are upper and lower triangular matrices with unit diagonal, and $D$ is a hermitian block-diagonal matrix with 1-by-1 and 2-by-2 diagonal blocks. $U$ and $L$ have 2-by-2 unit diagonal blocks corresponding to the 2-by-2 blocks of $D$. It utilizes one of ZHETRF routines of the LAPACK library. Function creates scmatrix object equal to the factorization of hermitian calling matrix. Function throws cvmexception in case of convergence error. Function is mostly designed to be used for subsequent calls of ZHETRS, ZHECON and ZHETRI routines of the LAPACK library. Currently it's used internally in det() flow when argument is symmetric but not positive-definite.

See Also
http://www.netlib.org/lapack
basic_scmatrix::bunch_kaufman()
Parameters
[out]nPivotsPivot indices array.
Returns
Result object.

Definition at line 36566 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::identity ( )
inline

Identity matrix.

Sets calling square complex matrix to be equal to identity matrix and returns reference to the matrix changed.

Example:
using namespace cvm;
scmatrix m(3);
m.randomize_real(0.,3.);
m.randomize_imag(-1.,2.);
std::cout << m << std::endl;
std::cout << m.identity();
prints
(1.31162,-0.52501) (2.8612,-0.531144) (1.31849,0.547838)
(1.19929,1.48253) (0.535417,0.41316) (0.459883,1.7019)
(0.415937,-0.491134) (2.0969,-0.218024) (0.545305,1.17866)
(1,0) (0,0) (0,0)
(0,0) (1,0) (0,0)
(0,0) (0,0) (1,0)
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 36573 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::vanish ( )
inline

Set matrix to zero.

Sets every element of calling square complex matrix to be equal to zero and returns reference to the matrix changed. This function is faster than, for example, set(TR) with zero parameter passed.

Example:
using namespace cvm;
scmatrix m(3);
m.randomize_real(0.,3.);
m.randomize_imag(-1.,2.);
std::cout << m << std::endl;
std::cout << m.vanish();
prints
(1.34834,-0.758385) (0.837825,-0.225532) (0.367687,0.791833)
(2.23698,-0.183142) (2.6878,0.741111) (0.495865,0.698904)
(0.584124,0.00491348) (1.31574,0.687643) (0.482131,1.66482)
(0,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
(0,0) (0,0) (0,0)
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 36580 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::randomize_real ( TR  dFrom,
TR  dTo 
)
inline

Randomizer (real part)

Fills real part of calling complex matrix with pseudo-random numbers distributed between dFrom and dTo. It returns reference to the matrix changed.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (3);
cmatrix m(2,3);
m.randomize_real(-1., 2.);
std::cout << m;
prints
(1.090e+00,0.000e+00) (-6.375e-01,0.000e+00) (1.248e+00,0.000e+00)
(-1.272e-01,0.000e+00) (-8.557e-01,0.000e+00) (4.848e-01,0.000e+00)
Parameters
[in]dFromFirst limit.
[in]dToSecond limit.
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 36586 of file cvm.h.

template<typename TR, typename TC>
basic_schmatrix& basic_schmatrix< TR, TC >::randomize_imag ( TR  dFrom,
TR  dTo 
)
inline

Randomizer (imaginary part)

Fills imaginary part of calling complex matrix with pseudo-random numbers distributed between dFrom and dTo. It returns reference to the matrix changed.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (3);
cmatrix m(2,3);
m.randomize_imag(-1., 2.);
std::cout << m;
prints
(0.000e+00,1.113e+00) (0.000e+00,6.615e-01) (0.000e+00,1.017e+00)
(0.000e+00,-3.397e-01) (0.000e+00,1.577e+00) (0.000e+00,8.071e-01)
Parameters
[in]dFromFirst limit.
[in]dToSecond limit.
Returns
Reference to changed calling matrix.

Reimplemented from basic_scmatrix< TR, TC >.

Definition at line 36592 of file cvm.h.

template<typename TR, typename TC>
TR basic_schmatrix< TR, TC >::norminf ( ) const
inlineoverridevirtual

Infinity norm.

Infinity norm of calling array that for vectors is defined as

\[ {\|x\|}_\infty=\max_{i=1,\dots,n} |x_i| \]

and for matrices as

\[ {\|A\|}_\infty=\max_{i=1,\dots,m} \sum_{j=1}^{n} |a_{ij}|, \]

Example:
using namespace cvm;
double a[] = {1., 2., 3., -4., 5., -6.};
const rvector v (a, 3);
const rmatrix m (a, 2, 3);
std::cout << v << v.norminf () << std::endl;
std::cout << m << m.norminf () << std::endl;
prints
1 2 3
3
1 3 5
2 -4 -6
12
Returns
treal Norm value

Reimplemented from Matrix< TR, TC >.

Definition at line 36599 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
bool basic_schmatrix< TR, TC >::is_positive_definite ( ) const
inline

Is calling hermitian matrix positive definite?

Returns
true if calling hermitian matrix is positive definite.

Definition at line 36609 of file cvm.h.

template<typename TR, typename TC>
bool basic_schmatrix< TR, TC >::equilibrate ( basic_rvector< TR > &  vScalings,
CVector vB 
) throw (cvmexception)
inline

Matrix equilibration.

Useful for further solve and solve_lu calling.

Returns
true if equilibration was needed and performed.

Definition at line 36631 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
bool basic_schmatrix< TR, TC >::equilibrate ( basic_rvector< TR > &  vScalings,
BaseCMatrix mB 
) throw (cvmexception)
inline

Matrix equilibration.

Useful for further solve and solve_lu calling.

Returns
true if equilibration was needed and performed.

Definition at line 36650 of file cvm.h.

Here is the call graph for this function:


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