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_cvector< TR, TC > Class Template Reference

End-user class encapsulating vector of complex numbers. More...

#include <cvm.h>

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

Public Member Functions

 basic_cvector ()
 Default constructor.
 basic_cvector (tint nSize)
 Constructor.
 basic_cvector (std::initializer_list< TC > list)
 Constructor.
 basic_cvector (tint nSize, TC c)
 Constructor.
 basic_cvector (TC *pd, tint nSize, tint nIncr=1)
 Constructor.
 basic_cvector (const TC *pd, tint nSize, tint nIncr=1)
 Constructor.
 basic_cvector (const basic_cvector &v)
 Copy constructor.
 basic_cvector (basic_cvector &&a) noexcept
 Move constructor.
 basic_cvector (const TR *pRe, const TR *pIm, tint nSize, tint nIncrRe=1, tint nIncrIm=1)
 Constructor.
 basic_cvector (const basic_rvector< TR > &vRe, const basic_rvector< TR > &vIm)
 Constructor.
 basic_cvector (const TR *pA, tint nSize, bool bRealPart=true, tint nIncr=1)
 Constructor.
 basic_cvector (const basic_rvector< TR > &v, bool bRealPart=true)
 Constructor.
basic_rvector< TR > real ()
 Real part (l-value)
basic_rvector< TR > imag ()
 Imaginary part (l-value)
basic_cvectoroperator= (const basic_cvector &v) throw (cvmexception)
 Assignment operator.
basic_cvectoroperator= (basic_cvector &&a) throw (cvmexception)
 Move assignment operator.
basic_cvectorassign (const TC *pd, tint nIncr=1)
 External array assignment.
basic_cvectorassign (tint n, const TC *pd, tint nIncr=1)
 External array assignment (to tail)
basic_cvectorassign (tint n, const TC *pd, tint nSize, tint nIncr)
 External array assignment (range)
basic_cvectorassign (tint n, const basic_cvector &v) throw (cvmexception)
 Subvector assignment.
basic_cvectorset (TC c)
 Sets all elements to one value.
basic_cvectorassign_real (const basic_rvector< TR > &vRe) throw (cvmexception)
 Assignment to real part.
basic_cvectorassign_imag (const basic_rvector< TR > &vIm) throw (cvmexception)
 Assignment to imaginary part.
basic_cvectorset_real (TR d)
 Fills real part.
basic_cvectorset_imag (TR d)
 Fills imaginary part.
basic_cvectorresize (tint nNewSize) throw (cvmexception)
 Changes size of vector.
bool operator== (const basic_cvector &v) const
 Vector comparison.
bool operator!= (const basic_cvector &v) const
 Vector comparison.
basic_cvectoroperator<< (const basic_cvector &v) throw (cvmexception)
 Vector replacement.
basic_cvector operator+ (const basic_cvector &v) const throw (cvmexception)
 Addition operator.
basic_cvector operator- (const basic_cvector &v) const throw (cvmexception)
 Subtraction operator.
basic_cvectorsum (const basic_cvector &v1, const basic_cvector &v2) throw (cvmexception)
 Sum of vectors.
basic_cvectordiff (const basic_cvector &v1, const basic_cvector &v2) throw (cvmexception)
 Difference of vectors.
basic_cvectoroperator+= (const basic_cvector &v) throw (cvmexception)
 Increment operator.
basic_cvectoroperator-= (const basic_cvector &v) throw (cvmexception)
 Decrement operator.
basic_cvector operator- () const throw (cvmexception)
 Unary minus operator.
basic_cvector operator* (TR dMult) const
 Multiply by real number operator.
basic_cvector operator/ (TR dDiv) const throw (cvmexception)
 Divide by real number operator.
basic_cvector operator* (TC cMult) const
 Multiply by complex number operator.
basic_cvector operator/ (TC cDiv) const throw (cvmexception)
 Divide by complex number operator.
basic_cvectoroperator*= (TR dMult)
 Multiply by real number and assign.
basic_cvectoroperator/= (TR dDiv)
 Divide by real number and assign.
basic_cvectoroperator*= (TC cMult)
 Multiply by complex number and assign.
basic_cvectoroperator/= (TC cDiv)
 Divide by complex number and assign.
basic_cvectornormalize ()
 Vector normalizer.
basic_cvector operator~ () const throw (cvmexception)
 Conjugation operator.
basic_cvectorconj (const basic_cvector &v) throw (cvmexception)
 Conjugation.
basic_cvectorconj ()
 Conjugation.
TC operator* (const basic_cvector &v) const throw (cvmexception)
 Scalar product.
TC operator% (const basic_cvector &v) const throw (cvmexception)
 Scalar product.
basic_cvector operator* (const basic_cmatrix< TR, TC > &m) const throw (cvmexception)
 Vector-matrix product.
basic_cvectormult (const basic_cvector &v, const basic_cmatrix< TR, TC > &m) throw (cvmexception)
 Vector-matrix product.
basic_cvectormult (const basic_cmatrix< TR, TC > &m, const basic_cvector &v) throw (cvmexception)
 Matrix-vector product.
basic_cmatrix< TR, TC > rank1update_u (const basic_cvector &v) const
 Rank-1 update (unconjugated)
basic_cmatrix< TR, TC > rank1update_c (const basic_cvector &v) const
 Rank-1 update (conjugated)
basic_cvectorsolve (const basic_scmatrix< TR, TC > &mA, const basic_cvector &vB, TR &dErr) throw (cvmexception)
 Linear solver.
basic_cvectorsolve_tran (const basic_scmatrix< TR, TC > &mA, const basic_cvector &vB, TR &dErr) throw (cvmexception)
 Linear solver (transposed)
basic_cvectorsolve_conj (const basic_scmatrix< TR, TC > &mA, const basic_cvector &vB, TR &dErr) throw (cvmexception)
 Linear solver (conjugated)
basic_cvectorsolve (const basic_scmatrix< TR, TC > &mA, const basic_cvector &vB) throw (cvmexception)
 Linear solver.
basic_cvectorsolve_tran (const basic_scmatrix< TR, TC > &mA, const basic_cvector &vB) throw (cvmexception)
 Linear solver (transposed)
basic_cvectorsolve_conj (const basic_scmatrix< TR, TC > &mA, const basic_cvector &vB) throw (cvmexception)
 Linear solver (conjugated)
basic_cvector operator/ (const basic_scmatrix< TR, TC > &mA) const throw (cvmexception)
 Linear solver operator (transposed)
basic_cvector operator% (const basic_scmatrix< TR, TC > &mA) const throw (cvmexception)
 Linear solver operator.
basic_cvectorsolve_lu (const basic_scmatrix< TR, TC > &mA, const basic_scmatrix< TR, TC > &mLU, const tint *pPivots, const basic_cvector &vB, TR &dErr) throw (cvmexception)
 LU factorization based linear solver.
basic_cvectorsolve_lu (const basic_scmatrix< TR, TC > &mA, const basic_scmatrix< TR, TC > &mLU, const tint *pPivots, const basic_cvector &vB) throw (cvmexception)
 LU factorization based linear solver.
basic_cvectorgels (bool conjugate, const basic_cmatrix< TR, TC > &mA, const basic_cvector &vB, TC &cErr) throw (cvmexception)
 Overdetermined or underdetermined linear solver.
basic_cvectorgelsy (const basic_cmatrix< TR, TC > &mA, const basic_cvector &vB, tint &rank, TR tol=basic_cvmMachSp< TR >()) throw (cvmexception)
 Linear least squares problem.
basic_cvectorgelss (const basic_cmatrix< TR, TC > &mA, const basic_cvector &vB, basic_rvector< TR > &sv, tint &rank, TR tol=basic_cvmMachSp< TR >()) throw (cvmexception)
 Linear least squares problem.
basic_cvectorgelsd (const basic_cmatrix< TR, TC > &mA, const basic_cvector &vB, basic_rvector< TR > &sv, tint &rank, TR tol=basic_cvmMachSp< TR >()) throw (cvmexception)
 Linear least squares problem.
basic_cvectoreig (const basic_srmatrix< TR > &mA) throw (cvmexception)
 Eigenvalues.
basic_cvectoreig (const basic_srmatrix< TR > &mA, basic_scmatrix< TR, TC > &mEigVect, bool bRightVect=true) throw (cvmexception)
 Eigenvalues anf eigenvectors.
basic_cvectoreig (const basic_scmatrix< TR, TC > &mA) throw (cvmexception)
 Eigenvalues.
basic_cvectoreig (const basic_scmatrix< TR, TC > &mA, basic_scmatrix< TR, TC > &mEigVect, bool bRightVect=true) throw (cvmexception)
 Eigenvalues.
basic_cvectorgeneig (const basic_srmatrix< TR > &mA, const basic_srmatrix< TR > &mB, basic_rvector< TR > &vBeta) throw (cvmexception)
 Generalized eigenvalues of real matrices.
basic_cvectorgeneig (const basic_srmatrix< TR > &mA, const basic_srmatrix< TR > &mB, basic_rvector< TR > &vBeta, basic_scmatrix< TR, TC > &mEigVectLeft, basic_scmatrix< TR, TC > &mEigVectRight) throw (cvmexception)
 Generalized eigenvalues of real matrices.
basic_cvectorgeneig (const basic_srmatrix< TR > &mA, const basic_srmatrix< TR > &mB, basic_rvector< TR > &vBeta, basic_scmatrix< TR, TC > &mEigVect, bool bRightVect=true) throw (cvmexception)
 Generalized eigenvalues of real matrices.
basic_cvectorgeneig (const basic_scmatrix< TR, TC > &mA, const basic_scmatrix< TR, TC > &mB, basic_cvector &vBeta) throw (cvmexception)
 Generalized eigenvalues of complex matrices.
basic_cvectorgeneig (const basic_scmatrix< TR, TC > &mA, const basic_scmatrix< TR, TC > &mB, basic_cvector &vBeta, basic_scmatrix< TR, TC > &mEigVectLeft, basic_scmatrix< TR, TC > &mEigVectRight) throw (cvmexception)
 Generalized eigenvalues of complex matrices.
basic_cvectorgeneig (const basic_scmatrix< TR, TC > &mA, const basic_scmatrix< TR, TC > &mB, basic_cvector &vBeta, basic_scmatrix< TR, TC > &mEigVect, bool bRightVect=true) throw (cvmexception)
 Generalized eigenvalues of complex matrices.
basic_cvectorgemv (bool bLeft, const basic_cmatrix< TR, TC > &m, TC dAlpha, const basic_cvector &v, TC dBeta) throw (cvmexception)
 Generic matrix-vector operation.
basic_cvectorgbmv (bool bLeft, const basic_scbmatrix< TR, TC > &m, TC dAlpha, const basic_cvector &v, TC dBeta) throw (cvmexception)
 Generic band matrix-vector operation.
basic_cvectorrandomize_real (TR dFrom, TR dTo)
 Randomizer (real part)
basic_cvectorrandomize_imag (TR dFrom, TR dTo)
 Randomizer (imaginary part)
- 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.
TC & operator() (tint n) throw (cvmexception)
 Reference to element (l-value)
TC operator() (tint n) const throw (cvmexception)
 Value of element (not l-value)
TC & operator[] (tint n) throw (cvmexception)
 Reference to element (l-value)
TC operator[] (tint n) const throw (cvmexception)
 Value of element (not l-value)
basic_arrayassign (const TC *p)
 Assignment from external array.
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.
virtual TR norminf () const
 Infinity norm.
virtual TR norm1 () const
 1-norm
virtual TR norm2 () const
 2-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 Attributes inherited from basic_array< TR, TC >
tint msz
 Number of elements of type TC allocated.
tint mincr
 Increment (distance) between elements (default is 1, i.e. elements follow each other)
std::shared_ptr< TC > mp
 native data pointer
TC * mpf
 Foreign data pointer.

Detailed Description

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

End-user class encapsulating vector of complex numbers.

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

See Also
basic_array

Definition at line 6155 of file cvm.h.

Constructor & Destructor Documentation

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

Default constructor.

Creates zero size cvector. No memory gets allocated.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
std::cout << v.size() << std::endl;
v.resize (3);
v(1) = std::complex<double>(1.5, -1.);
std::cout << v;
prints
0
(1.50e+00,-1.00e+00) (0.00e+00,0.00e+00) (0.00e+00,0.00e+00)
See Also
resize()

Definition at line 6184 of file cvm.h.

Here is the caller graph for this function:

template<typename TR, typename TC>
basic_cvector< TR, TC >::basic_cvector ( tint  nSize)
inlineexplicit

Constructor.

Creates cvector object of given size. It throws cvmexception in case of non-positive size passed or memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
cvector v(3);
std::cout << v.size() << std::endl;
std::cout << v;
prints
3
(0.00e+00,0.00e+00) (0.00e+00,0.00e+00) (0.00e+00,0.00e+00)
Parameters
[in]nSizeSize of cvector (must be positive).

Definition at line 6210 of file cvm.h.

template<typename TR, typename TC>
basic_cvector< TR, TC >::basic_cvector ( std::initializer_list< TC >  list)
inline

Constructor.

Creates cvector object and fills it with values provided in the initializer list. This constructor is available only if your compiler supports initializer lists (GCC 4.4 and higher, MS Visual Studio 2013 and higher, Apple LLVM 5.0 and higher).

Example:
using namespace cvm;
std::cout.setf(std::ios::scientific | std::ios::left);
std::cout.precision(3);
cvector v = { tcomplex(1.2, 3.4), tcomplex(3.4, 5.6), 99.99 };
std::cout << v;
prints
(1.200e+000,3.400e+000) (3.400e+000,5.600e+000) (9.999e+001,0.000e+000)
Also, if your compiler supports user-defined literals (GCC 4.7 and higher, Apple LLVM 5.0 and higher) you can write
cvector v = { 1.2+3.4_i, 3.4+5.6_i, 99.99 };
Parameters
[in]listInitializer list as shown above.

Definition at line 6242 of file cvm.h.

template<typename TR, typename TC>
basic_cvector< TR, TC >::basic_cvector ( tint  nSize,
TC  c 
)
inline

Constructor.

Creates cvector object of given size and fills it with given complex number. It throws cvmexception in case of non-positive size passed or memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
cvector v (3, std::complex<double>(1.5, -1.));
std::cout << v;
prints
(1.50e+00,-1.00e+00) (1.50e+00,-1.00e+00) (1.50e+00,-1.00e+00)
Parameters
[in]nSizeSize of cvector (must be positive).
[in]cNumber to fill cvector in.

Definition at line 6274 of file cvm.h.

template<typename TR, typename TC>
basic_cvector< TR, TC >::basic_cvector ( TC *  pd,
tint  nSize,
tint  nIncr = 1 
)
inline

Constructor.

Creates cvector object of given size with given increment between elements (default is 1). It throws cvmexception in case of non-positive size passed. Unlike others, this constructor does not allocate memory. It just shares memory with array pointed to by pd using distance between elements equal to nIncr. It is intented to make possible the following syntax:

cmatrix m (10, 20);
cvector v (20);
m[1] = v; // assigns v to the 1st row of m

And for example this code...

cmatrix m (10,20);
cvector vRow = m[1];

will also call this constructor and memory will be shared.

If you need the code like this with memory allocation, use the following:

cmatrix m (10,20);
cvector vRow (m.msize());
vRow = m[1];
Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4., 5., 6., 7., 8.};
cvector v1 ((std::complex<double>*) a, 2, 2);
std::cout << v1;
v1(2) = std::complex<double> (9.99, 9.99);
std::cout << v1 << std::endl;
for (int i = 0; i < 6; i++) {
std::cout << a[i] << " ";
}
std::cout << std::endl;
cvector v2 ((std::complex<double>*) a, 3);
std::cout << v2;
prints
(1.00e+000,2.00e+000) (5.00e+000,6.00e+000)
(1.00e+000,2.00e+000) (9.99e+000,9.99e+000)
1.00e+000 2.00e+000 3.00e+000 4.00e+000 9.99e+000 9.99e+000
(1.00e+000,2.00e+000) (3.00e+000,4.00e+000) (9.99e+000,9.99e+000)
See Also
http://cvmlib.com/faq.htm
Parameters
[in]pdPointer to array to share memory with.
[in]nSizeSize of cvector (must be positive).
[in]nIncrIncrement between elements, default is 1 (one after another).

Definition at line 6343 of file cvm.h.

template<typename TR, typename TC>
basic_cvector< TR, TC >::basic_cvector ( const TC *  pd,
tint  nSize,
tint  nIncr = 1 
)
inline

Constructor.

Creates cvector of nSize elements with increment 1 in between. This is const version, it allocates memory and copies every nIncr-th element (deep copy) from external array pointed to by pd parameter. It copies nSize elements total.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
const double a[] = {1., 2., 3., 4., 5., 6., 7., 8.};
cvector v1 ((const std::complex<double>*) a, 2, 2);
std::cout << v1;
v1(2) = std::complex<double> (9.99, 9.99);
std::cout << v1 << std::endl;
for (int i = 0; i < 6; i++) {
std::cout << a[i] << " ";
}
std::cout << std::endl;
cvector v2 ((std::complex<double>*) a, 3);
std::cout << v2;
prints
(1.00e+000,2.00e+000) (5.00e+000,6.00e+000)
(1.00e+000,2.00e+000) (9.99e+000,9.99e+000)
1.00e+000 2.00e+000 3.00e+000 4.00e+000 5.00e+000 6.00e+000
(1.00e+000,2.00e+000) (3.00e+000,4.00e+000) (5.00e+000,6.00e+000)
See Also
http://cvmlib.com/faq.htm
Parameters
[in]pdConst pointer to external array.
[in]nSizeNumber of array elements.
[in]nIncrIncrement between external array elements.

Definition at line 6388 of file cvm.h.

template<typename TR, typename TC>
basic_cvector< TR, TC >::basic_cvector ( const basic_cvector< TR, TC > &  v)
inline

Copy constructor.

Creates cvector object of size equal to size of vector v and sets every element of created vector to be equal to appropriate element of v. Constructor throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4., 5., 6., 7., 8., 9., 10.};
const cvector v ((std::complex<double>*) a, 3, 2);
cvector vc (v);
vc(1) = std::complex<double>(7.77,8.88);
std::cout << vc;
std::cout << v;
prints
(7.77e+00,8.88e+00) (5.00e+00,6.00e+00) (9.00e+00,1.00e+01)
(1.00e+00,2.00e+00) (5.00e+00,6.00e+00) (9.00e+00,1.00e+01)
Parameters
[in]vVector to copy from.

Definition at line 6421 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector< TR, TC >::basic_cvector ( basic_cvector< TR, TC > &&  a)
inline

Move constructor.

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

cvector a(b + c);

or this

cvector 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]arvalue reference to other array.

Definition at line 6444 of file cvm.h.

template<typename TR, typename TC>
basic_cvector< TR, TC >::basic_cvector ( const TR *  pRe,
const TR *  pIm,
tint  nSize,
tint  nIncrRe = 1,
tint  nIncrIm = 1 
)
inline

Constructor.

Creates cvector of nSize elements with increment 1 in between. Constructor copies every nIncrRe-th element of array pointed to by pRe and every nIncrIm-th element of array pointed to by pIm to real and imaginary part of the object created. It throws cvmexception in case of non-positive size passed or memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double re[] = {1., 2., 3., 4., 5.};
double im[] = {5., 4., 3., 2., 1.};
cvector v (re, im, 3, 2);
std::cout << v;
re[0] = 7.77;
std::cout << v;
const double rec[] = {1., 2., 3.};
const double imc[] = {5., 4., 3.};
const cvector vc (rec, imc, 3);
std::cout << vc;
prints
(1.00e+00,5.00e+00) (3.00e+00,4.00e+00) (5.00e+00,3.00e+00)
(1.00e+00,5.00e+00) (3.00e+00,4.00e+00) (5.00e+00,3.00e+00)
(1.00e+00,5.00e+00) (2.00e+00,4.00e+00) (3.00e+00,3.00e+00)
Parameters
[in]pReConst pointer to external real array.
[in]pImConst pointer to external imaginary array.
[in]nSizeNumber of array elements.
[in]nIncrReIncrement between external real array elements.
[in]nIncrImIncrement between external iumaginary array elements.

Definition at line 6488 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector< TR, TC >::basic_cvector ( const basic_rvector< TR > &  vRe,
const basic_rvector< TR > &  vIm 
)
inline

Constructor.

Creates cvector object of size equal to vRe.size() and vIm.size() and copies vectors vRe and vIm to real and imaginary part of the object created. It throws cvmexception in case of non-equal sizes of the parameters passed or memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
rvector vr(3), vi(3);
vr[1] = 1.;
vr[2] = 2.;
vr[3] = 3.;
vi[1] = 5.;
vi[2] = 4.;
vi[3] = 3.;
const cvector vc(vr, vi);
std::cout << vc;
prints
(1.00e+00,5.00e+00) (2.00e+00,4.00e+00) (3.00e+00,3.00e+00)
Parameters
[in]vRervector of real values.
[in]vImrvector of imaginary values.

Definition at line 6527 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector< TR, TC >::basic_cvector ( const TR *  pA,
tint  nSize,
bool  bRealPart = true,
tint  nIncr = 1 
)
inline

Constructor.

Creates cvector object of size equal to nSize and copies every nIncr-th element of array pointed to by pA to real (if bRealPart is true) or imaginary (if bRealPart is false) part of the object created. Constructor throws cvmexception in case of non-positive size passed or memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
const double a[] = {1., 2., 3., 4., 5.};
cvector v1 (a, 3, false, 2);
cvector v2 (a, 2);
std::cout << v1 << v2;
prints
(0.00e+00,1.00e+00) (0.00e+00,3.00e+00) (0.00e+00,5.00e+00)
(1.00e+00,0.00e+00) (2.00e+00,0.00e+00)
Parameters
[in]pAbasic_array of real values.
[in]nSizeSize of cvector (must be positive).
[in]bRealPartTrue to copy pA values to real part of the object created.
[in]nIncrCopy every nIncr-th element of pA array.

Definition at line 6566 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector< TR, TC >::basic_cvector ( const basic_rvector< TR > &  v,
bool  bRealPart = true 
)
inlineexplicit

Constructor.

Creates cvector object of size equal to v.size() and copies every element of vector v to real (if bRealPart is true) or imaginary (if bRealPart is false) part of the object created. Constructor throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
rvector vr (3);
vr(1) = 1.;
vr(2) = 2.;
vr(3) = 3.;
cvector v1 (vr);
cvector v2 (vr, false);
std::cout << v1 << v2;
prints
(1.00e+00,0.00e+00) (2.00e+00,0.00e+00) (3.00e+00,0.00e+00)
(0.00e+00,1.00e+00) (0.00e+00,2.00e+00) (0.00e+00,3.00e+00)
Parameters
[in]vrvector to copy values from.
[in]bRealPartTrue to copy v values to real part of the object created.

Definition at line 6603 of file cvm.h.

Here is the call graph for this function:

Member Function Documentation

template<typename TR, typename TC>
basic_rvector<TR> basic_cvector< TR, TC >::real ( )
inline

Real part (l-value)

Creates rvector object of size equal to size of calling vector sharing memory with its real part. In other words, the vector returned is l-value.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
cvector vc(3);
vc.set(std::complex<double>(1.,1.));
std::cout << vc << vc.real();
vc.real()(1) = 7.77;
std::cout << vc;
prints
(1.00e+00,1.00e+00) (1.00e+00,1.00e+00) (1.00e+00,1.00e+00)
1.00e+00 1.00e+00 1.00e+00
(7.77e+00,1.00e+00) (1.00e+00,1.00e+00) (1.00e+00,1.00e+00)
See Also
imag()
Returns
Result object (vector of real parts).

Definition at line 6637 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_rvector<TR> basic_cvector< TR, TC >::imag ( )
inline

Imaginary part (l-value)

Creates rvector object of size equal to size of calling vector sharing memory with its imaginary part. In other words, the vector returned is l-value.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
cvector vc(3);
vc.set(std::complex<double>(1.,1.));
std::cout << vc << vc.imag();
vc.imag()(1) = 7.77;
std::cout << vc;
prints
(1.00e+00,1.00e+00) (1.00e+00,1.00e+00) (1.00e+00,1.00e+00)
1.00e+00 1.00e+00 1.00e+00
(1.00e+00,7.77e+00) (1.00e+00,1.00e+00) (1.00e+00,1.00e+00)
See Also
real()
Returns
Result object (vector of imaginary parts).

Definition at line 6668 of file cvm.h.

Here is the call graph for this function:

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

Assignment operator.

Sets every element of calling cvector to be equal to appropriate element of vector v and returns reference to the vector changed. Operator throws cvmexception in case of different vector sizes.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
cvector vc(3);
vc.set(std::complex<double>(1.,1.));
std::cout << vc << vc.imag();
vc.imag()(1) = 7.77;
std::cout << vc;
prints
(1.00e+00,1.00e+00) (1.00e+00,1.00e+00) (1.00e+00,1.00e+00)
1.00e+00 1.00e+00 1.00e+00
(1.00e+00,7.77e+00) (1.00e+00,1.00e+00) (1.00e+00,1.00e+00)
Parameters
[in]vcvector to assign from.
Returns
Reference to changed calling vector.

Definition at line 6700 of file cvm.h.

Here is the call graph for this function:

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

Move assignment operator.

Implements move semantics introduced in new C++ standard. Moves data ownership from other vector 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]arvalue reference to other array.

Definition at line 6720 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::assign ( const TC *  pd,
tint  nIncr = 1 
)
inline

External array assignment.

Sets every element of calling vector to be equal to every nIncr-th element of array pointed to by parameter pd and returns reference to the vector changed.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
const double a[] = {1., 2., 3., 4., 5., 6., 7.};
cvector v1(3);
cvector v2(2);
v1.assign((const std::complex<double>*) a);
v2.assign((const std::complex<double>*) a, 2);
std::cout << v1;
std::cout << v2;
prints
(1.00e+00,2.00e+00) (3.00e+00,4.00e+00) (5.00e+00,6.00e+00)
(1.00e+00,2.00e+00) (5.00e+00,6.00e+00)
Parameters
[in]pdConst pointer to external array.
[in]nIncrIncrement between elements in external array.
Returns
Reference to changed calling vector.

Definition at line 6757 of file cvm.h.

Here is the caller graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::assign ( tint  n,
const TC *  pd,
tint  nIncr = 1 
)
inline

External array assignment (to tail)

Sets every element of calling vector, starting from CVM0 based n-th one till the end, to be equal to every nIncr-th element of array pointed to by parameter pd and returns reference to the vector changed.

Example:
using namespace cvm;
const double a[] = {1., 2., 3., 4., 5., 6., 7., 8., 9., 10.};
cvector v1(4);
cvector v2(4);
v1.assign(3, (const std::complex<double>*) a);
v2.assign(2, (const std::complex<double>*) a, 2);
std::cout << v1;
std::cout << v2;
prints
(0,0) (0,0) (1,2) (3,4)
(0,0) (1,2) (5,6) (9,10)
Parameters
[in]nStart cvector index to assign from.
[in]pdConst pointer to external array.
[in]nIncrIncrement between elements in external array.
Returns
Reference to changed calling vector.

Definition at line 6792 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::assign ( tint  n,
const TC *  pd,
tint  nSize,
tint  nIncr 
)
inline

External array assignment (range)

Sets every element of calling vector, starting from CVM0 based n-th one, up to nSize total, to be equal to every nIncr-th element of array pointed to by parameter pd and returns reference to the vector changed. If n+nSize goes beyond rvector boundaries assignment stops at the last element.

Example:
using namespace cvm;
const double a[] = {1., 2., 3., 4., 5., 6., 7., 8., 9., 10.};
cvector v1(4);
v1.assign(2, (const std::complex<double>*) a, 2, 3);
std::cout << v1;
prints
(0,0) (1,2) (7,8) (0,0)
Parameters
[in]nStart cvector index to assign from.
[in]pdConst pointer to external array.
[in]nSizeMaximum number of elements to be assigned.
[in]nIncrIncrement between elements in external array.
Returns
Reference to changed calling vector.

Definition at line 6827 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::assign ( tint  n,
const basic_cvector< TR, TC > &  v 
) throw (cvmexception)
inline

Subvector assignment.

Sets every element of calling vector's sub-vector beginning with CVM0 based index n to vector v and returns reference to the vector changed. Function throws cvmexception if n is not positive or v.size()+n-1 is greater than calling vector's size.

Example:
using namespace cvm;
cvector v1(5);
cvector v2(2);
v1.set(std::complex<double>(1.,1.));
v2.set(std::complex<double>(2.,2.));
v1.assign(3, v2);
std::cout << v1;
prints
(1,1) (1,1) (2,2) (2,2) (1,1)
Parameters
[in]nStart cvector index to assign from.
[in]vConst reference to vector.
Returns
Reference to changed calling vector.

Definition at line 6862 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::set ( TC  c)
inline

Sets all elements to one value.

Sets every element of calling vector to be equal to parameter c and returns reference to the vector changed.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
cvector v(3);
v.set(std::complex<double>(3.,1.));
std::cout << v;
prints
(3.00e+00,1.00e+00) (3.00e+00,1.00e+00) (3.00e+00,1.00e+00)
Parameters
[in]cValue to set to.
Returns
Reference to changed calling vector.

Reimplemented from basic_array< TR, TC >.

Definition at line 6891 of file cvm.h.

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::assign_real ( const basic_rvector< TR > &  vRe) throw (cvmexception)
inline

Assignment to real part.

Sets real part of every element of calling vector to be equal to appropriate element of vector vRe and returns reference to the vector changed. Function throws cvmexception in case of different sizes of the operands.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
rvector v(3);
cvector vc(3);
v(1) = 1.; v(2) = 2.; v(3) = 3.;
std::cout << vc;
prints
(1.00e+00,0.00e+00) (2.00e+00,0.00e+00) (3.00e+00,0.00e+00)
Parameters
[in]vRervector to assign to real part.
Returns
Reference to changed calling vector.

Definition at line 6923 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::assign_imag ( const basic_rvector< TR > &  vIm) throw (cvmexception)
inline

Assignment to imaginary part.

Sets imaginary part of every element of calling vector to be equal to appropriate element of vector vIm and returns reference to the vector changed. Function throws cvmexception in case of different sizes of the operands.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
rvector v(3);
cvector vc(3);
v(1) = 1.; v(2) = 2.; v(3) = 3.;
std::cout << vc;
prints
(0.00e+00,1.00e+00) (0.00e+00,2.00e+00) (0.00e+00,3.00e+00)
Parameters
[in]vImrvector to assign to imaginary part.
Returns
Reference to changed calling vector.

Definition at line 6956 of file cvm.h.

Here is the call graph for this function:

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

Fills real part.

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

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
cvector v(3);
v.set_real(1.);
std::cout << v;
prints
(1.00e+00,0.00e+00) (1.00e+00,0.00e+00) (1.00e+00,0.00e+00)
Parameters
[in]dtreal number to fill real part with.
Returns
Reference to changed calling vector.

Definition at line 6985 of file cvm.h.

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::set_imag ( TR  d)
inline

Fills imaginary part.

Sets imaginary part of every element of calling vector to be equal to parameter x and returns reference to the vector changed.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
cvector v(3);
v.set_imag(1.);
std::cout << v;
prints
(0.00e+00,1.00e+00) (0.00e+00,1.00e+00) (0.00e+00,1.00e+00)
Parameters
[in]dtreal number to fill imaginary part with.
Returns
Reference to changed calling vector.

Definition at line 7013 of file cvm.h.

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

Changes size of vector.

Changes size of calling vector to be equal to nNewSize and returns reference to the vector changed. In case of increasing of its size, the vector is filled up with zeroes. Function throws cvmexception in case of negative size passed or memory allocation failure.

Example:
using namespace cvm;
try {
double a[] = {1., 2., 3., 4.};
rvector v (a, 3);
std::cout << v;
v.resize(2);
std::cout << v;
v.resize(4);
std::cout << v;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(1,2) (3,4) (5,6)
(1,2) (3,4)
(1,2) (3,4) (0,0) (0,0)
Parameters
[in]nNewSizeNew size.
Returns
Reference to changed calling vector.

Reimplemented from basic_array< TR, TC >.

Definition at line 7054 of file cvm.h.

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

Vector comparison.

Operator compares calling vector with vector v and returns true if they have same sizes 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., 2., 3., 4.};
cvector v1 ((std::complex<double>*)a, 2);
cvector v2 (2);
v2(1) = std::complex<double>(1.,2.);
v2(2) = std::complex<double>(3.,4.);
std::cout << (v1 == v2) << std::endl;
prints
1
See Also
operator !=()
Parameters
[in]vcvector to compare to.
Returns
bool Result of comparison.

Definition at line 7089 of file cvm.h.

Here is the caller graph for this function:

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

Vector comparison.

Operator compares calling vector with vector v 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., 2., 3., 4.};
cvector v1 ((std::complex<double>*)a, 2);
cvector v2 (2);
std::cout << (v1 != v2) << std::endl;
prints
1
See Also
operator ==()
Parameters
[in]vcvector to compare to.
Returns
bool Result of comparison.

Definition at line 7119 of file cvm.h.

Here is the call graph for this function:

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

Vector replacement.

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

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
try {
cvector v(2);
cvector vc(3);
v(1) = std::complex<double> (1.,2.);
v(2) = std::complex<double> (3.,4.);
std::cout << v << vc << std::endl;
vc << v;
std::cout << vc;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(1.00e+000,2.00e+000) (3.00e+000,4.00e+000)
(0.00e+000,0.00e+000) (0.00e+000,0.00e+000) (0.00e+000,0.00e+000)
(1.00e+000,2.00e+000) (3.00e+000,4.00e+000)
See Also
operator =()
Parameters
[in]vcvector to replace by.
Returns
Reference to changed calling vector.

Definition at line 7160 of file cvm.h.

Here is the call graph for this function:

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

Addition operator.

Creates object of type cvector as sum of calling vector and vector v. Operator throws cvmexception in case of different sizes of the operands or memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
try {
cvector va(3);
cvector vb(3);
va.set(std::complex<double>(1.,1.));
vb.set(std::complex<double>(2.,2.));
std::cout << va + vb;
std::cout << va + va;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(3.00e+000,3.00e+000) (3.00e+000,3.00e+000) (3.00e+000,3.00e+000)
(2.00e+000,2.00e+000) (2.00e+000,2.00e+000) (2.00e+000,2.00e+000)
See Also
sum()
Parameters
[in]vcvector to add to calling one.
Returns
Sum of vectors.

Definition at line 7203 of file cvm.h.

Here is the call graph for this function:

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

Subtraction operator.

Creates object of type cvector as difference of calling vector and vector v. It throws cvmexception in case of different sizes of the operands or memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
try {
cvector va(3);
cvector vb(3);
va.set(std::complex<double> (1.,1.));
vb.set(std::complex<double> (2.,2.));
std::cout << va - vb;
std::cout << va - va;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(-1.00e+000,-1.00e+000) (-1.00e+000,-1.00e+000) (-1.00e+000,-1.00e+000)
(0.00e+000,0.00e+000) (0.00e+000,0.00e+000) (0.00e+000,0.00e+000)
See Also
diff()
Parameters
[in]vcvector to subtract from calling one.
Returns
Difference of vectors.

Definition at line 7246 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::sum ( const basic_cvector< TR, TC > &  v1,
const basic_cvector< TR, TC > &  v2 
) throw (cvmexception)
inline

Sum of vectors.

Assigns sum of vectors v1 and v2 to calling vector and returns reference to the vector changed. It throws cvmexception in case of different sizes of the operands.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
try {
cvector va(3);
cvector vb(3);
cvector v(3);
va.set(std::complex<double> (1.,1.));
vb.set(std::complex<double> (2.,2.));
std::cout << v.sum(va, vb);
std::cout << v.sum(v, va);
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(3.00e+000,3.00e+000) (3.00e+000,3.00e+000) (3.00e+000,3.00e+000)
(4.00e+000,4.00e+000) (4.00e+000,4.00e+000) (4.00e+000,4.00e+000)
See Also
operator +()
Parameters
[in]v1First cvector summand.
[in]v2Second cvector summand.
Returns
Reference to changed calling vector.

Definition at line 7291 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::diff ( const basic_cvector< TR, TC > &  v1,
const basic_cvector< TR, TC > &  v2 
) throw (cvmexception)
inline

Difference of vectors.

Assigns difference of vectors v1 and v2 to calling vector and returns reference to the vector changed. It throws cvmexception in case of different sizes of the operands.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
try {
cvector va(3);
cvector vb(3);
cvector v(3);
va.set(std::complex<double> (1.,1.));
vb.set(std::complex<double> (2.,2.));
std::cout << v.diff(va, vb);
std::cout << v.diff(v, va);
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(-1.00e+000,-1.00e+000) (-1.00e+000,-1.00e+000) (-1.00e+000,-1.00e+000)
(-2.00e+000,-2.00e+000) (-2.00e+000,-2.00e+000) (-2.00e+000,-2.00e+000)
See Also
operator -()
Parameters
[in]v1First cvector subtrahend.
[in]v2Second cvector subtrahend.
Returns
Reference to changed calling vector.

Definition at line 7336 of file cvm.h.

Here is the call graph for this function:

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

Increment operator.

Adds vector v to calling vector and returns reference to the vector changed. It throws cvmexception in case of different sizes of the operands.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
try {
cvector v1(3);
cvector v2(3);
v1.set(std::complex<double> (1.,1.));
v2.set(std::complex<double> (2.,2.));
v1 += v2;
std::cout << v1;
// well, you can do this too, but temporary object would be created
v2 += v2;
std::cout << v2;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(3.00e+000,3.00e+000) (3.00e+000,3.00e+000) (3.00e+000,3.00e+000)
(4.00e+000,4.00e+000) (4.00e+000,4.00e+000) (4.00e+000,4.00e+000)
See Also
operator +()
sum()
Parameters
[in]vcvector to increment by.
Returns
Reference to changed calling vector.

Definition at line 7383 of file cvm.h.

Here is the call graph for this function:

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

Decrement operator.

Subtracts vector v from calling vector and returns reference to the vector changed. It throws cvmexception in case of different sizes of the operands.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
try {
cvector v1(3);
cvector v2(3);
v1.set(std::complex<double> (1.,1.));
v2.set(std::complex<double> (2.,2.));
v1 -= v2;
std::cout << v1;
// well, you can do this too, but temporary object would be created
v2 -= v2;
std::cout << v2;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(-1.00e+000,-1.00e+000) (-1.00e+000,-1.00e+000) (-1.00e+000,-1.00e+000)
(0.00e+000,0.00e+000) (0.00e+000,0.00e+000) (0.00e+000,0.00e+000)
See Also
operator -()
diff()
Parameters
[in]vcvector to decrement by.
Returns
Reference to changed calling vector.

Definition at line 7429 of file cvm.h.

Here is the call graph for this function:

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

Unary minus operator.

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

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
const cvector v ((std::complex<double>*) a, 2);
std::cout << - v;
prints
(-1.00e+000,-2.00e+000) (-3.00e+000,-4.00e+000)
Returns
Result object.

Definition at line 7459 of file cvm.h.

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

Multiply by real number operator.

Creates object of type cvector as product of calling vector and treal number dMult. It throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
const cvector v ((std::complex<double>*) a, 2);
std::cout << v * 5.;
prints
(5.00e+000,1.00e+001) (1.50e+001,2.00e+001)
See Also
operator *=()
Parameters
[in]dMultNumber to multiply by.
Returns
Result object.

Definition at line 7493 of file cvm.h.

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

Divide by real number operator.

Creates object of type cvector as quotient of calling vector and treal 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;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
const cvector v ((std::complex<double>*) a, 2);
std::cout << v / 4.;
prints
(2.50e-001,5.00e-001) (7.50e-001,1.00e+000)
See Also
operator /=()
Parameters
[in]dDivNumber to divide by.
Returns
Result object.

Definition at line 7528 of file cvm.h.

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

Multiply by complex number operator.

Creates object of type cvector as product of calling vector and tcomplex number cMult. It throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
const cvector v ((std::complex<double>*) a, 2);
std::cout << v * std::complex<double>(1.,1.);
prints
(-1.00e+000,3.00e+000) (-1.00e+000,7.00e+000)
See Also
operator *=()
Parameters
[in]cMultNumber to multiply by.
Returns
Result object.

Definition at line 7561 of file cvm.h.

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

Divide by complex number operator.

Creates object of type cvector as quotient of calling vector and tcomplex 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;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
const cvector v ((std::complex<double>*) a, 2);
std::cout << v / std::complex<double>(1.,1.);
prints
(1.50e+000,5.00e-001) (3.50e+000,5.00e-001)
See Also
operator /=()
Parameters
[in]cDivNumber to divide by.
Returns
Result object.

Definition at line 7596 of file cvm.h.

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

Multiply by real number and assign.

Multiplies calling vector by treal number dMult and returns reference to the vector changed.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
cvector v ((std::complex<double>*) a, 2);
std::cout << (v *= 2.);
prints
(2.00e+000,4.00e+000) (6.00e+000,8.00e+000)
See Also
operator *(TR) const
Parameters
[in]dMultNumber to multiply by.
Returns
Reference to changed calling vector.

Definition at line 7627 of file cvm.h.

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

Divide by real number and assign.

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

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
cvector v ((std::complex<double>*) a, 2);
std::cout << (v /= 2.);
prints
(5.00e-001,1.00e+000) (1.50e+000,2.00e+000)
See Also
operator /(TR) const
Parameters
[in]dDivNumber to divide by.
Returns
Reference to changed calling vector.

Definition at line 7659 of file cvm.h.

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::operator*= ( TC  cMult)
inline

Multiply by complex number and assign.

Multiplies calling vector by tcomplex number cMult and returns reference to the vector changed.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
cvector v ((std::complex<double>*) a, 2);
v *= std::complex<double>(1.,1.);
std::cout << v;
prints
(-1.00e+000,3.00e+000) (-1.00e+000,7.00e+000)
See Also
operator *(TC) const
Parameters
[in]cMultNumber to multiply by.
Returns
Reference to changed calling vector.

Definition at line 7690 of file cvm.h.

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::operator/= ( TC  cDiv)
inline

Divide by complex number and assign.

Divides calling vector by tcomplex number cDiv and returns reference to the vector changed. It throws cvmexception if cDiv has absolute value equal or less than cvmMachMin() (the smallest normalized positive number).

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
cvector v ((std::complex<double>*) a, 2);
v /= std::complex<double>(1.,1.);
std::cout << v;
prints
(1.50e+000,5.00e-001) (3.50e+000,5.00e-001)
See Also
operator /(TC) const
Parameters
[in]cDivNumber to divide by.
Returns
Reference to changed calling vector.

Definition at line 7723 of file cvm.h.

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

Vector normalizer.

Normalizes calling vector 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;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
cvector v ((std::complex<double>*) a, 2);
std::cout << v.normalize();
std::cout << v.norm() << std::endl;
prints
(1.83e-01,3.65e-01) (5.48e-01,7.30e-01)
1.00e+00
Returns
Reference to changed calling vector.

Definition at line 7755 of file cvm.h.

template<typename TR, typename TC>
basic_cvector basic_cvector< TR, TC >::operator~ ( ) const throw (cvmexception)
inline

Conjugation operator.

Creates object of type cvector as complex conjugated calling vector It throws cvmexception in case of memory allocation failure.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
const cvector v ((std::complex<double>*) a, 2);
cvector vc(2);
std::cout << ~v;
std::cout << vc.conj(v);
std::cout << vc.conj();
prints
(1.00e+00,-2.00e+00) (3.00e+00,-4.00e+00)
(1.00e+00,-2.00e+00) (3.00e+00,-4.00e+00)
(1.00e+00,2.00e+00) (3.00e+00,4.00e+00)
See Also
conj()
Returns
Result object.

Definition at line 7789 of file cvm.h.

Here is the call graph for this function:

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

Conjugation.

Sets calling vector to be equal to vector v conjugated It throws cvmexception in case of different sizes of the operands)

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
const cvector v ((std::complex<double>*) a, 2);
cvector vc(2);
std::cout << ~v;
std::cout << vc.conj(v);
std::cout << vc.conj();
prints
(1.00e+00,-2.00e+00) (3.00e+00,-4.00e+00)
(1.00e+00,-2.00e+00) (3.00e+00,-4.00e+00)
(1.00e+00,2.00e+00) (3.00e+00,4.00e+00)
See Also
operator ~()
conj()
Returns
Reference to changed calling vector.

Definition at line 7824 of file cvm.h.

Here is the call graph for this function:

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

Conjugation.

Sets calling vector to be equal to conjugated itself.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
const cvector v ((std::complex<double>*) a, 2);
cvector vc(2);
std::cout << ~v;
std::cout << vc.conj(v);
std::cout << vc.conj();
prints
(1.00e+00,-2.00e+00) (3.00e+00,-4.00e+00)
(1.00e+00,-2.00e+00) (3.00e+00,-4.00e+00)
(1.00e+00,2.00e+00) (3.00e+00,4.00e+00)
See Also
operator ~()
conj(const basic_cvector<TR,TC>&)
Returns
Reference to changed calling vector.

Definition at line 7859 of file cvm.h.

Here is the call graph for this function:

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

Scalar product.

Scalar product of calling vector and vector v. It throws cvmexception if the operands have different sizes.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
double b[] = {1., -1., 1., 2.};
const cvector v1((std::complex<double>*) a, 2);
const cvector v2((std::complex<double>*) b, 2);
std::cout << v1 * v2 << std::endl;
prints
(-2.00e+00,1.10e+01)
Parameters
[in]vcvector to compute scalar product with.
Returns
TC Scalar product as tcomplex.

Definition at line 7890 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
TC basic_cvector< TR, TC >::operator% ( const basic_cvector< TR, TC > &  v) const throw (cvmexception)
inline

Scalar product.

Scalar product of complex conjugated calling vector and vector v. It throws cvmexception if the operands have different sizes.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 4.};
double b[] = {1., -1., 1., 2.};
const cvector v1((std::complex<double>*) a, 2);
const cvector v2((std::complex<double>*) b, 2);
std::cout << v1 % v2 << std::endl;
std::cout << ~v1 * v2 << std::endl;
prints
(1.00e+01,-1.00e+00)
(1.00e+01,-1.00e+00)
Parameters
[in]vcvector to compute scalar product with.
Returns
TC Scalar product as tcomplex.

Definition at line 7923 of file cvm.h.

Here is the call graph for this function:

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

Vector-matrix product.

Creates object of type cvector as product of calling vector and matrix m. Use mult() function to avoid new object creation. Operator throws cvmexception if calling vector's size differs from number of rows in matrix m.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., 3., 2., 1.};
double b[] = {1., -1., 1., 2., -2., 1.,
3., -2., 1., 2., -1., 3.};
const cvector v((std::complex<double>*) a, 3);
const cmatrix m((std::complex<double>*) b, 3, 2);
std::cout << v << m << std::endl << v * m;
prints
(1.00e+00,2.00e+00) (3.00e+00,3.00e+00) (2.00e+00,1.00e+00)
(1.00e+00,-1.00e+00) (3.00e+00,-2.00e+00)
(1.00e+00,2.00e+00) (1.00e+00,2.00e+00)
(-2.00e+00,1.00e+00) (-1.00e+00,3.00e+00)
(-5.00e+00,1.00e+01) (-1.00e+00,1.80e+01)
Parameters
[in]mcmatrix to multiply calling vector by.
Returns
Result object.

Definition at line 7963 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::mult ( const basic_cvector< TR, TC > &  v,
const basic_cmatrix< TR, TC > &  m 
) throw (cvmexception)
inline

Vector-matrix product.

Sets calling vector to be equal to product of vector v by matrix m and returns reference to the object changed.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
try {
double a[] = {1., 2., 3., 1., 2., 3.};
double b[] = {1., -1., 1., -1., 1., -1.,
2., -1., 2., -1., 2., -1.};
const cvector v ((std::complex<double>*) a, 3);
const cmatrix m ((std::complex<double>*) b, 3, 2);
const scmatrix sm ((std::complex<double>*) b, 2);
cvector vm (2);
std::cout << vm.mult(v, m) << std::endl;
std::cout << sm << std::endl;
std::cout << vm.mult(vm, sm);
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(1.20e+01,0.00e+00) (1.80e+01,6.00e+00)
(1.00e+00,-1.00e+00) (1.00e+00,-1.00e+00)
(1.00e+00,-1.00e+00) (2.00e+00,-1.00e+00)
(3.60e+01,-2.40e+01) (5.40e+01,-1.80e+01)
See Also
operator *(const basic_cmatrix<TR,TC>&) const
Parameters
[in]vcvector to multiply by matrix m.
[in]mcmatrix to multiply vector v by.
Returns
Reference to changed calling vector.

Definition at line 8012 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::mult ( const basic_cmatrix< TR, TC > &  m,
const basic_cvector< TR, TC > &  v 
) throw (cvmexception)
inline

Matrix-vector product.

Sets calling vector to be equal to product of matrix m by vector v by and returns reference to the object changed.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
try {
double a[] = {1., 2., 3., 1., 2., 3.};
double b[] = {1., -1., 1., -1., 1., -1.,
2., -1., 2., -1., 2., -1.};
const cvector v ((std::complex<double>*) a, 3);
const cmatrix m ((std::complex<double>*) b, 2, 3);
const scmatrix sm ((std::complex<double>*) b, 2);
cvector vm (2);
std::cout << vm.mult(m, v) << std::endl;
std::cout << sm << std::endl;;
std::cout << vm.mult(vm, sm);
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(1.40e+01,3.00e+00) (1.70e+01,4.00e+00)
(1.00e+00,-1.00e+00) (1.00e+00,-1.00e+00)
(1.00e+00,-1.00e+00) (2.00e+00,-1.00e+00)
(3.80e+01,-2.40e+01) (5.50e+01,-2.00e+01)
Parameters
[in]mcmatrix to multiply by vector v.
[in]vcvector to multiply matrix m by.
Returns
Reference to changed calling vector.

Definition at line 8060 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cmatrix<TR,TC> basic_cvector< TR, TC >::rank1update_u ( const basic_cvector< TR, TC > &  v) const
inline

Rank-1 update (unconjugated)

Creates object of type cmatrix as rank-1 update of calling vector and vector v. The rank-1 update of vector-column $x$ of size $m$ and vector-row $y$ of size $n$ is defined as $m\times n$ matrix

\[\begin{pmatrix} x_1 y_1 & x_1 y_2 & \cdots & x_1 y_n \\ x_2 y_1 & x_2 y_2 & \cdots & x_2 y_n \\ \hdotsfor{4} \\ x_m y_1 & x_m y_2 & \cdots & x_m y_n \end{pmatrix} = \begin{pmatrix} x_1 \\ x_2 \\ \vdots \\ x_m \end{pmatrix} \begin{pmatrix} y_1 & y_2 & \cdots & y_n \end{pmatrix}.\]

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., -2., -1., 1.};
double b[] = {4., 5., 3., 2.};
cvector v1((std::complex<double>*) a, 3);
cvector v2((std::complex<double>*) b, 2);
std::cout << v1.rank1update_u (v2);
prints
(-6.00e+00,1.30e+01) (-1.00e+00,8.00e+00)
(2.20e+01,7.00e+00) (1.30e+01,0.00e+00)
(-9.00e+00,-1.00e+00) (-5.00e+00,1.00e+00)
Parameters
[in]vcvector to compute rank-1 update with.
Returns
Result object.

Definition at line 8111 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cmatrix<TR,TC> basic_cvector< TR, TC >::rank1update_c ( const basic_cvector< TR, TC > &  v) const
inline

Rank-1 update (conjugated)

Creates object of type cmatrix as rank-1 update of calling vector and complex conjugated vector v. The rank-1 update of vector-column $x$ of size $m$ and vector-row $y$ of size $n$ is defined as $m\times n$ matrix

\[\begin{pmatrix} x_1 y_1^* & x_1 y_2^* & \cdots & x_1 y_n^* \\ x_2 y_1^* & x_2 y_2^* & \cdots & x_2 y_n^* \\ \hdotsfor{4} \\ x_m y_1^* & x_m y_2^* & \cdots & x_m y_n^* \end{pmatrix} = \begin{pmatrix} x_1 \\ x_2 \\ \vdots \\ x_m \end{pmatrix} \begin{pmatrix} y_1^* & y_2^* & \hdots & y_n^* \end{pmatrix},\]

where $y_i^*$ is $i$-th complex conjugated element of $y$.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
double a[] = {1., 2., 3., -2., -1., 1.};
double b[] = {4., 5., 3., 2.};
cvector v1((std::complex<double>*) a, 3);
cvector v2((std::complex<double>*) b, 2);
std::cout << v1.rank1update_c (v2) << std::endl;
std::cout << v1.rank1update_u (~v2);
prints
(1.40e+01,3.00e+00) (7.00e+00,4.00e+00)
(2.00e+00,-2.30e+01) (5.00e+00,-1.20e+01)
(1.00e+00,9.00e+00) (-1.00e+00,5.00e+00)
(1.40e+01,3.00e+00) (7.00e+00,4.00e+00)
(2.00e+00,-2.30e+01) (5.00e+00,-1.20e+01)
(1.00e+00,9.00e+00) (-1.00e+00,5.00e+00)
Parameters
[in]vcvector to compute rank-1 update with.
Returns
Result object.

Definition at line 8168 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::solve ( const basic_scmatrix< TR, TC > &  mA,
const basic_cvector< TR, TC > &  vB,
TR &  dErr 
) throw (cvmexception)
inline

Linear solver.

Sets calling vector to be equal to solution $x$ of linear equation $A*x=b$ where parameter mA is square matrix $A$ and parameter vB is vector $b$. Function returns reference to the vector changed. It also sets output parameter dErr to be equal to the norm of computation error. Function throws cvmexception in case of inappropriate sizes of the objects or when matrix $A$ is close to singular.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
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);
double dErr = 0.;
std::cout << vx.solve (ma, vb, dErr);
std::cout << dErr << std::endl;
std::cout << ma * vx - vb;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(3.5200000e+00,6.4000000e-01) (2.2400000e+00,-1.3200000e+00)
3.2788531e-15
(-7.7715612e-16,4.4408921e-16) (0.0000000e+00,0.0000000e+00)
See Also
solve(const basic_scmatrix<TR,TC>&,const basic_cvector<TR,TC>&)
Parameters
[in]mAscmatrix $A$.
[in]vBcvector $b$.
[out]dErrNorm of computation error.
Returns
Reference to changed calling vector.

Definition at line 8221 of file cvm.h.

Here is the caller graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::solve_tran ( const basic_scmatrix< TR, TC > &  mA,
const basic_cvector< TR, TC > &  vB,
TR &  dErr 
) throw (cvmexception)
inline

Linear solver (transposed)

Sets calling vector to be equal to solution $x$ of linear equation $A^T*x=b$ (which is equivalent to $x*A=b$) where parameter mA is square matrix $A$ and parameter vB is vector $b$. Function returns reference to the vector changed. It also sets output parameter dErr to be equal to the norm of computation error. Function throws cvmexception in case of inappropriate sizes of the objects 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);
double dErr = 0.;
std::cout << vx.solve_tran (ma, vb, dErr);
std::cout << dErr << std::endl;
std::cout << vx * ma - vb;
std::cout << !ma * vx - vb;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(+1.6000000e-001,-8.8000000e-001) (+1.5600000e+000,-8.0000000e-002)
+3.7480513e-015
(+0.0000000e+000,+0.0000000e+000) (+0.0000000e+000,+4.4408921e-016)
(+0.0000000e+000,+0.0000000e+000) (+0.0000000e+000,+4.4408921e-016)
See Also
solve_tran(const basic_scmatrix<TR,TC>&,const basic_cvector<TR,TC>&)
Parameters
[in]mAscmatrix $A$.
[in]vBcvector $b$.
[out]dErrNorm of computation error.
Returns
Reference to changed calling vector.

Definition at line 8273 of file cvm.h.

Here is the caller graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::solve_conj ( const basic_scmatrix< TR, TC > &  mA,
const basic_cvector< TR, TC > &  vB,
TR &  dErr 
) throw (cvmexception)
inline

Linear solver (conjugated)

Sets calling vector to be equal to solution $x$ of linear equation $A^H*x=b$ (here $A^H$ is conjugated $A$) where parameter mA is square complex matrix $A$ and parameter vB is complex vector $b$. Function returns reference to the vector changed. It also sets output parameter dErr to be equal to the norm of computation error. Functions throws cvmexception in case of inappropriate sizes of the objects or when the 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);
double dErr = 0.;
std::cout << vx.solve_conj (ma, vb, dErr);
std::cout << dErr << std::endl;
std::cout << ~ma * vx - vb;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(+2.3200000e+000,+3.7600000e+000) (+2.1200000e+000,+1.6000000e-001)
+2.1932508e-015
(+0.0000000e+000,-4.4408921e-016) (-8.8817842e-016,-8.8817842e-016)
See Also
solve_conj(const basic_scmatrix<TR,TC>&,const basic_cvector<TR,TC>&)
Parameters
[in]mAscmatrix $A$.
[in]vBcvector $b$.
[out]dErrNorm of computation error.
Returns
Reference to changed calling vector.

Definition at line 8323 of file cvm.h.

Here is the caller graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::solve ( const basic_scmatrix< TR, TC > &  mA,
const basic_cvector< TR, TC > &  vB 
) throw (cvmexception)
inline

Linear solver.

Sets calling vector to be equal to solution $x$ of linear equation $A*x=b$ where parameter mA is square matrix $A$ and parameter vB is vector $b$. Function returns reference to the vector changed. Function throws cvmexception in case of inappropriate sizes of the objects or when matrix $A$ is close to singular.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
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);
std::cout << vx.solve (ma, vb);
std::cout << ma * vx - vb;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(3.5200000e+00,6.4000000e-01) (2.2400000e+00,-1.3200000e+00)
(-7.7715612e-16,4.4408921e-16) (0.0000000e+00,0.0000000e+00)
See Also
solve(const basic_scmatrix<TR,TC>&,const basic_cvector<TR,TC>&,TR&)
Parameters
[in]mAscmatrix $A$.
[in]vBcvector $b$.
Returns
Reference to changed calling vector.

Definition at line 8366 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::solve_tran ( const basic_scmatrix< TR, TC > &  mA,
const basic_cvector< TR, TC > &  vB 
) throw (cvmexception)
inline

Linear solver (transposed)

Sets calling vector to be equal to solution $x$ of linear equation $A^T*x=b$ (which is equivalent to $x*A=b$) where parameter mA is square matrix $A$ and parameter vB is vector $b$. Function returns reference to the vector changed. Function throws cvmexception in case of inappropriate sizes of the objects 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);
std::cout << vx.solve_tran (ma, vb);
std::cout << vx * ma - vb;
std::cout << !ma * vx - vb;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(+1.6000000e-001,-8.8000000e-001) (+1.5600000e+000,-8.0000000e-002)
(+0.0000000e+000,+0.0000000e+000) (+0.0000000e+000,+4.4408921e-016)
(+0.0000000e+000,+0.0000000e+000) (+0.0000000e+000,+4.4408921e-016)
See Also
solve_tran(const basic_scmatrix<TR,TC>&,const basic_cvector<TR,TC>&,TR&)
Parameters
[in]mAscmatrix $A$.
[in]vBcvector $b$.
Returns
Reference to changed calling vector.

Definition at line 8415 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::solve_conj ( const basic_scmatrix< TR, TC > &  mA,
const basic_cvector< TR, TC > &  vB 
) throw (cvmexception)
inline

Linear solver (conjugated)

Sets calling vector to be equal to solution $x$ of linear equation $A^H*x=b$ (here $A^H$ is conjugated $A$) where parameter mA is square complex matrix $A$ and parameter vB is complex vector $b$. Function returns reference to the vector changed. Functions throws cvmexception in case of inappropriate sizes of the objects or when the 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);
std::cout << vx.solve_conj (ma, vb);
std::cout << ~ma * vx - vb;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(+2.3200000e+000,+3.7600000e+000) (+2.1200000e+000,+1.6000000e-001)
(+0.0000000e+000,-4.4408921e-016) (-8.8817842e-016,-8.8817842e-016)
See Also
solve_conj(const basic_scmatrix<TR,TC>&,const basic_cvector<TR,TC>&,TR&)
Parameters
[in]mAscmatrix $A$.
[in]vBcvector $b$.
Returns
Reference to changed calling vector.

Definition at line 8462 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector basic_cvector< TR, TC >::operator/ ( const basic_scmatrix< TR, TC > &  mA) const throw (cvmexception)
inline

Linear solver operator (transposed)

Returns solution $x$ of linear equation $A^T*x=b$ (which is equivalent to $x*A=b$) where parameter mA is square matrix $A$ and calling vector is $b$. This operator throws cvmexception in case of inappropriate sizes of the objects 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 = vb / ma;
std::cout << vx * ma - vb;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(+0.0000000e+000,+0.0000000e+000) (+0.0000000e+000,+4.4408921e-016)
See Also
solve_tran(const basic_scmatrix<TR,TC>&,const basic_cvector<TR,TC>&)
Parameters
[in]mAscmatrix $A$.
Returns
Result object (solution of the equation).

Definition at line 8507 of file cvm.h.

template<typename TR, typename TC>
basic_cvector basic_cvector< TR, TC >::operator% ( const basic_scmatrix< TR, TC > &  mA) const throw (cvmexception)
inline

Linear solver operator.

Returns solution $x$ of linear equation $A*x=b$ where parameter mA is square matrix $A$ and calling vector is $b$. This operator throws cvmexception in case of inappropriate sizes of the objects 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 = vb % ma;
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)
See Also
solve(const basic_scmatrix<TR,TC>&,const basic_cvector<TR,TC>&)
Parameters
[in]mAscmatrix $A$.
Returns
Result object (solution of the equation).

Definition at line 8547 of file cvm.h.

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::solve_lu ( const basic_scmatrix< TR, TC > &  mA,
const basic_scmatrix< TR, TC > &  mLU,
const tint pPivots,
const basic_cvector< TR, TC > &  vB,
TR &  dErr 
) throw (cvmexception)
inline

LU factorization based linear solver.

Sets calling vector to be equal to solution $x$ of linear equation $A*x=b$ where parameter mA is square matrix $A$, parameter mLU is LU factorization (see scmatrix::low_up() ) of matrix $A$, parameter pPivots is an array of pivot numbers created while factorizing matrix $A$ and parameter vB is vector $b$. Function returns reference to the vector changed. It also sets output parameter dErr to be equal to the norm of computation error. This function is useful when you need to solve few linear equations of kind $Ax=b$ with the same matrix $A$ and different vectors $b$. In such case you save on matrix $A$ factorization since it's needed to be performed just one time. Throws cvmexception in case of inappropriate sizes of the objects 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 b1[] = {1., 2., 5., -3.};
double b2[] = {3., -1., 1., 7.};
scmatrix ma((std::complex<double>*) m, 2);
scmatrix mLU(2);
cvector vb1((std::complex<double>*) b1, 2);
cvector vb2((std::complex<double>*) b2, 2);
cvector vx1(2);
cvector vx2(2);
iarray nPivots(2);
double dErr = 0.;
mLU.low_up(ma, nPivots);
std::cout << vx1.solve_lu (ma, mLU, nPivots, vb1, dErr);
std::cout << dErr << std::endl;
std::cout << vx2.solve_lu (ma, mLU, nPivots, vb2, dErr);
std::cout << dErr << std::endl;
std::cout << ma * vx1 - vb1 << ma * vx2 - vb2;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(+3.5200000e+000,+6.4000000e-001) (+2.2400000e+000,-1.3200000e+000)
+3.2191768e-015
(+2.2800000e+000,+1.9600000e+000) (+3.6000000e-001,+5.2000000e-001)
+2.1974410e-015
(-6.6613381e-016,+4.4408921e-016) (+0.0000000e+000,+0.0000000e+000)
(-4.4408921e-016,+0.0000000e+000) (-8.8817842e-016,+8.8817842e-016)
See Also
scmatrix::low_up()
Parameters
[in]mAscmatrix $A$.
[in]mLULU factorization of matrix $A$.
[in]pPivotspivots vector.
[in]vBrvector $b$.
[out]dErrNorm of computation error.
Returns
Reference to changed calling vector.

Definition at line 8614 of file cvm.h.

Here is the call graph for this function:

Here is the caller graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::solve_lu ( const basic_scmatrix< TR, TC > &  mA,
const basic_scmatrix< TR, TC > &  mLU,
const tint pPivots,
const basic_cvector< TR, TC > &  vB 
) throw (cvmexception)
inline

LU factorization based linear solver.

Sets calling vector to be equal to solution $x$ of linear equation $A*x=b$ where parameter mA is square matrix $A$, parameter mLU is LU factorization (see scmatrix::low_up() ) of matrix $A$, parameter pPivots is an array of pivot numbers created while factorizing matrix $A$ and parameter vB is vector $b$. Function returns reference to the vector changed. This function is useful when you need to solve few linear equations of kind $Ax=b$ with the same matrix $A$ and different vectors $b$. In such case you save on matrix $A$ factorization since it's needed to be performed just one time. Throws cvmexception in case of inappropriate sizes of the objects 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 b1[] = {1., 2., 5., -3.};
double b2[] = {3., -1., 1., 7.};
scmatrix ma((std::complex<double>*) m, 2);
scmatrix mLU(2);
cvector vb1((std::complex<double>*) b1, 2);
cvector vb2((std::complex<double>*) b2, 2);
cvector vx1(2);
cvector vx2(2);
iarray nPivots(2);
mLU.low_up(ma, nPivots);
std::cout << vx1.solve_lu (ma, mLU, nPivots, vb1);
std::cout << vx2.solve_lu (ma, mLU, nPivots, vb2);
std::cout << ma * vx1 - vb1 << ma * vx2 - vb2;
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(+3.5200000e+000,+6.4000000e-001) (+2.2400000e+000,-1.3200000e+000)
(+2.2800000e+000,+1.9600000e+000) (+3.6000000e-001,+5.2000000e-001)
(-6.6613381e-016,+4.4408921e-016) (+0.0000000e+000,+0.0000000e+000)
(-4.4408921e-016,+0.0000000e+000) (-8.8817842e-016,+8.8817842e-016)
See Also
scmatrix::low_up()
Parameters
[in]mAscmatrix $A$.
[in]mLULU factorization of matrix $A$.
[in]pPivotspivots vector.
[in]vBrvector $b$.
Returns
Reference to changed calling vector.

Definition at line 8679 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::gels ( bool  conjugate,
const basic_cmatrix< TR, TC > &  mA,
const basic_cvector< TR, TC > &  vB,
TC &  cErr 
) throw (cvmexception)
inline

Overdetermined or underdetermined linear solver.

Solves overdetermined or underdetermined linear system

\[ A*x=b \]

for $m\times n$ matrix $A$ (or conjugateed one) where $b$ is a vector of length $k$ and $k=m$ in non-conjugateed case and $k=n$ otherwise. The algorithm uses QR or LQ factorization of $A$. It is assumed that $A$ has full rank, infinity returned otherwise. Internally function uses ZGELS LAPACK routines. If $m>n$ and conjugate=false or $m<n$ and conjugate=true, then the system is overdetermined, thus the algorithm tries to find the least squares solution $x$ of the problem

\[ \|A*x-b\|_2\to\min\quad\text{or}\quad\|A^H*x-b\|_2\to\min\ \]

respectively. Complex number cErr returns residuals. The system is underdetermined otherwise, and the algorithm finds its minimum norm solution. In this case cErr is set to zero. In both cases the solution computed satisfies $\newcommand{\pinv}{\operatorname{pinv}} x=\pinv(A)*b$, but this algorithm is faster than pseudo inversion. Function sets calling object to be the solution and returns reference to it. It 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 (7);
cvm::cmatrix a(7, 5);
cvm::cvector b(7), bt(5);
tcomplex cErr, cErrc;
a.randomize_real(-10., 10.);
a.randomize_imag(-10., 10.);
b.randomize_real(-10., 10.);
b.randomize_imag(-10., 10.);
bt.randomize_real(-10., 10.);
bt.randomize_imag(-10., 10.);
x.gels(false, a, b, cErr);
xt.gels(true, a, bt, cErrc);
std::cout << (a.pinv()*b - x).norm() << std::endl;
std::cout << cErr << std::endl;
std::cout << (~a.pinv()*bt - xt).norm() << std::endl;
std::cout << cErrc << std::endl;
prints
+1.0444134e-015
(-3.2469414e+001,+5.8375852e+000)
+1.9786436e-015
(+0.0000000e+000,+0.0000000e+000)
See Also
cmatrix::pinv()
Parameters
[in]conjugateTrue to compute for conjugated matrix $A$.
[in]mAcmatrix $A$.
[in]vBcvector $b$.
[out]cErrNorm of computation error.
Returns
Reference to changed calling vector.

Definition at line 8750 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::gelsy ( const basic_cmatrix< TR, TC > &  mA,
const basic_cvector< TR, TC > &  vB,
tint rank,
TR  tol = basic_cvmMachSp<TR>() 
) throw (cvmexception)
inline

Linear least squares problem.

Computes minimum-norm solution to linear least squares problem

\[ \|A*x-b\|_2\to\min \]

using complete orthogonal factorization of $m\times n$ matrix $A$. Here $b$ is a vector of length $m$. Matrix $A$ may be rank-deficient, function returns its effective rank in rank output parameter using tol tolerance. Internally function uses ZGELSY LAPACK routines, see more details about the algorithm in those routines' documentation. Matrix $A$ is passed as argument mA. Function sets calling object to be the solution and returns reference to it. It throws cvmexception in case of inappropriate sizes of the operands.

Example:
using namespace cvm;
cmatrix a(4, 5);
cvector b(4);
a.randomize_real(-10., 10.);
a.randomize_imag(-10., 10.);
b.randomize_real(-10., 10.);
b.randomize_imag(-10., 10.);
cvector x(5);
tint rank;
x.gelsy(a, b, rank);
std::cout << (a*x - b).norm() << std::endl;
std::cout << rank << " " << a.rank() << std::endl;
prints
+6.0116296e-015
+4 +4
See Also
gelss()
gelsd()
cmatrix::gelsy()
cmatrix::rank()
Parameters
[in]mAcmatrix $A$.
[in]vBcvector $b$.
[out]rankRank of matrix $A$.
[in]tolRank computation tolerance.
Returns
Reference to changed calling vector.

Definition at line 8813 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::gelss ( const basic_cmatrix< TR, TC > &  mA,
const basic_cvector< TR, TC > &  vB,
basic_rvector< TR > &  sv,
tint rank,
TR  tol = basic_cvmMachSp<TR>() 
) throw (cvmexception)
inline

Linear least squares problem.

Computes minimum-norm solution to linear least squares problem

\[ \|A*x-b\|_2\to\min \]

using singular value decomposition of $m\times n$ matrix $A$. Here $b$ is a vector of length $m$. Matrix $A$ may be rank-deficient, function returns its effective rank in rank output parameter using tol tolerance. This function also computes singular values of $A$ in decreasing order and returns them in sv output parameter having $\min(m,n)$ size. Internally function uses ZGELSS LAPACK routines, see more details about the algorithm in that routine's documentation. Matrix $A$ is passed as argument mA. Function sets calling object to be the solution and returns reference to it. It throws cvmexception in case of inappropriate sizes of the operands.

Example:
using namespace cvm;
cvm::cmatrix a(4, 5);
a.randomize_real(-10., 10.);
a.randomize_imag(-10., 10.);
b.randomize_real(-10., 10.);
b.randomize_imag(-10., 10.);
tint rank;
x.gelss(a, b, sv, rank);
std::cout << (a*x - b).norm() << std::endl;
std::cout << (sv - a.svd()).norm() << std::endl;
std::cout << rank << " " << a.rank() << std::endl;
prints
+9.4127010e-015
+4.4408921e-015
+4 +4
See Also
gelsy()
gelsd()
cmatrix::gelsy()
cmatrix::rank()
Parameters
[in]mAcmatrix $A$.
[in]vBcvector $b$.
[out]svSingular values of matrix $A$.
[out]rankRank of matrix $A$.
[in]tolRank computation tolerance.
Returns
Reference to changed calling vector.

Definition at line 8880 of file cvm.h.

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::gelsd ( const basic_cmatrix< TR, TC > &  mA,
const basic_cvector< TR, TC > &  vB,
basic_rvector< TR > &  sv,
tint rank,
TR  tol = basic_cvmMachSp<TR>() 
) throw (cvmexception)
inline

Linear least squares problem.

Computes minimum-norm solution to linear least squares problem

\[ \|A*x-b\|_2\to\min \]

using singular value decomposition of $m\times n$ matrix $A$ and divide and conquer method. Here $b$ is a vector of length $m$. Matrix $A$ may be rank-deficient, function returns its effective rank in rank output parameter using tol tolerance. This function also computes singular values of $A$ in decreasing order and returns them in sv output parameter having $\min(m,n)$ size. Internally function uses ZGELSD LAPACK routines, see more details about the algorithm in that routine's documentation. Matrix $A$ is passed as argument mA. Function sets calling object to be the solution and returns reference to it. It throws cvmexception in case of inappropriate sizes of the operands.

Example:
using namespace cvm;
cvm::cmatrix a(4, 5);
a.randomize_real(-10., 10.);
a.randomize_imag(-10., 10.);
b.randomize_real(-10., 10.);
b.randomize_imag(-10., 10.);
tint rank;
x.gelsd(a, b, sv, rank);
std::cout << (a*x - b).norm() << std::endl;
std::cout << (sv - a.svd()).norm() << std::endl;
std::cout << rank << " " << a.rank() << std::endl;
prints
+1.3491837e-014
+5.3290705e-015
+4 +4
See Also
gelss()
gelsy()
cmatrix::gelsy()
cmatrix::rank()
Parameters
[in]mAcmatrix $A$.
[in]vBcvector $b$.
[out]svSingular values of matrix $A$.
[out]rankRank of matrix $A$.
[in]tolRank computation tolerance.
Returns
Reference to changed calling vector.

Definition at line 8942 of file cvm.h.

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::eig ( const basic_srmatrix< TR > &  mA) throw (cvmexception)
inline

Eigenvalues.

Solves eigenvalue problem and sets calling vector to be equal to eigenvalues of real square matrix mA. The eigenvalue problem is defined as follows: given square matrix $A$, find the eigenvalues $\lambda$ and the corresponding eigenvectors $z$ that satisfy the equation

\[ Az = \lambda z. \]

Some eigenvalues may be complex even for real matrix $A$. Moreover, if real nonsymmetric matrix has a complex eigenvalue $a+bi$ corresponding to an eigenvector $z$, then $a-bi$ is also an eigenvalue. The eigenvalue $a-bi$ corresponds to the eigenvector whose elements are complex conjugate to the elements of $z$. Function returns reference to the vector changed and 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::showpos);
std::cout.precision (2);
try {
srmatrix m(3);
cvector vl(3);
m(1,1) = 0.1; m(1,2) = 0.2; m(1,3) = 0.1;
m(2,1) = 0.11; m(2,2) = 2.9; m(2,3) = -8.4;
m(3,1) = 0.; m(3,2) = 2.91; m(3,3) = 8.2;
std::cout << vl.eig (m);
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(+9.69e-002,+0.00e+000) (+5.55e+000,+4.17e+000) (+5.55e+000,-4.17e+000)
See Also
cvector::eig()
Parameters
[in]mAsrmatrix $A$.
Returns
Reference to changed calling vector.

Definition at line 8994 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::eig ( const basic_srmatrix< TR > &  mA,
basic_scmatrix< TR, TC > &  mEigVect,
bool  bRightVect = true 
) throw (cvmexception)
inline

Eigenvalues anf eigenvectors.

Solves eigenvalue problem and sets calling vector to be equal to eigenvalues of real square matrix mA. The eigenvalue problem is defined as follows: given square matrix $A$, find the eigenvalues $\lambda$ and the corresponding eigenvectors $z$ that satisfy the equation

\[ Az = \lambda z. \]

Some eigenvalues may be complex even for real matrix $A$. Moreover, if real nonsymmetric matrix has a complex eigenvalue $a+bi$ corresponding to an eigenvector $z$, then $a-bi$ is also an eigenvalue. The eigenvalue $a-bi$ corresponds to the eigenvector whose elements are complex conjugate to the elements of $z$. Function sets output parameter mEigVect to be equal to square matrix containing eigenvectors as columns. Function also computes "left" eigenvectors if parameter bRightVect is set to false. Left eigencectors satisfy

\[ z^HA = \lambda z^H. \]

Function returns reference to the vector changed and 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::showpos);
std::cout.precision (2);
try {
srmatrix m(3);
scmatrix me(3);
cvector vl(3);
m(1,1) = 0.1; m(1,2) = 0.2; m(1,3) = 0.1;
m(2,1) = 0.11; m(2,2) = 2.9; m(2,3) = -8.4;
m(3,1) = 0.; m(3,2) = 2.91; m(3,3) = 8.2;
scmatrix mc(m);
std::cout << vl.eig (m, me) << std::endl;
std::cout << mc * me(1) - me(1) * vl(1);
std::cout << mc * me(2) - me(2) * vl(2);
std::cout << mc * me(3) - me(3) * vl(3) << std::endl;
std::cout << vl.eig (m, me, false) << std::endl;
std::cout << ~(me(1)) * mc - ~(me(1)) * vl(1);
std::cout << ~(me(2)) * mc - ~(me(2)) * vl(2);
std::cout << ~(me(3)) * mc - ~(me(3)) * vl(3);
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(+9.69e-002,+0.00e+000) (+5.55e+000,+4.17e+000) (+5.55e+000,-4.17e+000)
(+2.78e-017,+0.00e+000) (+7.24e-017,+0.00e+000) (+4.23e-017,+0.00e+000)
(-2.22e-016,-1.94e-016) (-7.11e-015,-4.88e-015) (+0.00e+000,+1.78e-015)
(-2.22e-016,+1.94e-016) (-7.11e-015,+4.88e-015) (+0.00e+000,-1.78e-015)
(+9.69e-002,+0.00e+000) (+5.55e+000,+4.17e+000) (+5.55e+000,-4.17e+000)
(+0.00e+000,+0.00e+000) (+4.19e-017,+0.00e+000) (+1.13e-017,+0.00e+000)
(+0.00e+000,-4.16e-017) (+2.22e-016,-1.78e-015) (-8.88e-016,+2.44e-015)
(+0.00e+000,+4.16e-017) (+2.22e-016,+1.78e-015) (-8.88e-016,-2.44e-015)
See Also
cvector::eig()
Parameters
[in]mAsrmatrix $A$.
[out]mEigVectEigenvectors of mA.
[in]bRightVecttrue (default) to compute right eigenvectors.
Returns
Reference to changed calling vector.

Definition at line 9074 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::eig ( const basic_scmatrix< TR, TC > &  mA) throw (cvmexception)
inline

Eigenvalues.

Solves eigenvalue problem and sets calling vector to be equal to eigenvalues of complex square matrix mA. The eigenvalue problem is defined as follows: given square matrix $A$, find the eigenvalues $\lambda$ and the corresponding eigenvectors $z$ that satisfy the equation

\[ Az = \lambda z. \]

Function returns reference to the vector changed and 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::showpos);
std::cout.precision (2);
try {
double re[] = {1., 2., 3., 4., 5., 6., 7., 8., 9.};
double im[] = {9., 8., 7., 6., 5., 4., 3., 2., 1.};
scmatrix m(re, im, 3);
cvector vl(3);
std::cout << vl.eig (m);
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(+1.37e+001,+1.37e+001) (+1.32e+000,+1.32e+000) (+1.10e-015,+1.92e-015)
See Also
cvector::eig()
Parameters
[in]mAscmatrix $A$.
Returns
Reference to changed calling vector.

Definition at line 9120 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::eig ( const basic_scmatrix< TR, TC > &  mA,
basic_scmatrix< TR, TC > &  mEigVect,
bool  bRightVect = true 
) throw (cvmexception)
inline

Eigenvalues.

Solves eigenvalue problem and sets calling vector to be equal to eigenvalues of complex square matrix mA. The eigenvalue problem is defined as follows: given square matrix $A$, find the eigenvalues $\lambda$ and the corresponding eigenvectors $z$ that satisfy the equation

\[ Az = \lambda z. \]

Function sets output parameter mEigVect to be equal to square matrix containing eigenvectors as columns. Function also computes "left" eigenvectors if parameter bRightVect is set to false. Left eigencectors satisfy

\[ z^HA = \lambda z^H. \]

Function returns reference to the vector changed and 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::showpos);
std::cout.precision (2);
try {
double re[] = {1., 2., 3., 4., 5., 6., 7., 8., 9.};
double im[] = {9., 8., 7., 6., 5., 4., 3., 2., 1.};
scmatrix m(re, im, 3);
scmatrix me(3);
cvector vl(3);
std::cout << vl.eig (m, me) << std::endl;
std::cout << m * me(1) - me(1) * vl(1);
std::cout << m * me(2) - me(2) * vl(2);
std::cout << m * me(3) - me(3) * vl(3) << std::endl;
std::cout << vl.eig (m, me, false) << std::endl;
std::cout << ~(me(1)) * m - ~(me(1)) * vl(1);
std::cout << ~(me(2)) * m - ~(me(2)) * vl(2);
std::cout << ~(me(3)) * m - ~(me(3)) * vl(3);
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(+1.37e+001,+1.37e+001) (+1.32e+000,+1.32e+000) (+1.10e-015,+1.92e-015)
(+1.55e-015,+3.55e-015) (-8.44e-015,+1.07e-014) (-7.99e-015,+7.11e-015)
(-2.89e-015,+5.11e-015) (-1.39e-015,+5.00e-016) (+3.33e-016,+2.00e-015)
(-2.91e-016,+6.23e-016) (+3.25e-015,+1.86e-015) (+1.82e-015,-4.87e-016)
(+1.37e+001,+1.37e+001) (+1.32e+000,+1.32e+000) (+1.10e-015,+1.92e-015)
(+3.55e-015,+2.31e-014) (+1.24e-014,+3.55e-015) (+5.33e-015,-7.11e-015)
(+1.33e-015,+0.00e+000) (-2.33e-015,-4.33e-015) (-1.55e-015,-6.44e-015)
(-9.67e-016,-1.00e-019) (-6.44e-017,+4.33e-015) (+4.43e-018,+3.66e-015)
See Also
cvector::eig()
Parameters
[in]mAscmatrix $A$.
[out]mEigVectEigenvectors of mA.
[in]bRightVecttrue (default) to compute right eigenvectors.
Returns
Reference to changed calling vector.

Definition at line 9192 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::geneig ( const basic_srmatrix< TR > &  mA,
const basic_srmatrix< TR > &  mB,
basic_rvector< TR > &  vBeta 
) throw (cvmexception)
inline

Generalized eigenvalues of real matrices.

Solves generalized eigenvalue problem and sets calling complex vector to be equal to generalized eigenvalues of real square matrices mA and mB. A generalized eigenvalue for a pair of matrices $(A,B)$ is a scalar $\lambda$ or a ratio $\alpha / \beta = \lambda$, such that $A - \lambda B$ is singular. It is usually represented as the pair $(\alpha, \beta)$, as there is a reasonable interpretation for $\beta =0$ and even for both being zero. Function returns reference to the vector changed and throws cvmexception in case of inappropriate calling object sizes or in case of convergence error.

Example:
using namespace cvm;
try {
srmatrix a(5), b(5);
a.randomize(-10., 10.);
b.randomize(-10., 10.);
cvector alpha(5);
rvector beta(5);
alpha.geneig(a, b, beta);
for (int i = 1; i <= 5; ++i) {
std::cout << (scmatrix(a) - alpha[i] / beta[i] * scmatrix(b)).svd();
}
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
99.7613 88.2526 59.4809 39.999 2.58055e-015
99.7613 88.2526 59.4809 39.999 2.58055e-015
26.8437 17.4907 9.30819 6.88691 1.44137e-015
25.8221 19.8697 11.2967 4.08824 3.35041e-015
25.8221 19.8697 11.2967 4.08824 1.43823e-015
See Also
cvector::eig()
Parameters
[in]mAsrmatrix $A$.
[in]mBsrmatrix $B$.
[out]vBetacomputed $\beta$ values.
Returns
Reference to changed calling vector with computed $\alpha$ values

Definition at line 9244 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::geneig ( const basic_srmatrix< TR > &  mA,
const basic_srmatrix< TR > &  mB,
basic_rvector< TR > &  vBeta,
basic_scmatrix< TR, TC > &  mEigVectLeft,
basic_scmatrix< TR, TC > &  mEigVectRight 
) throw (cvmexception)
inline

Generalized eigenvalues of real matrices.

Solves generalized eigenvalue problem and sets calling complex vector to be equal to generalized eigenvalues of real square matrices mA and mB. A generalized eigenvalue for a pair of matrices $(A,B)$ is a scalar $\lambda$ or a ratio $\alpha / \beta = \lambda$, such that $A - \lambda B$ is singular. It is usually represented as the pair $(\alpha, \beta)$, as there is a reasonable interpretation for $\beta =0$ and even for both being zero. The right generalized eigenvector $v_j$ corresponding to the generalized eigenvalue $\lambda_j$ of $(A,B)$ satisfies

\[ A v_j = \lambda_j B v_j. \]

The left generalized eigenvector $u_j$ corresponding to the generalized eigenvalue $\lambda_j$ of $(A,B)$ satisfies

\[ u_j^H A = \lambda_j u_j^H B. \]

Function returns reference to the vector changed and throws cvmexception in case of inappropriate calling object sizes or in case of convergence error.

Example:
using namespace cvm;
try {
srmatrix a(5), b(5);
a.randomize(-10., 10.);
b.randomize(-10., 10.);
cvector alpha(5);
rvector beta(5);
scmatrix eigVectLeft(5), eigVectRight(5);
int i;
alpha.geneig(a, b, beta, eigVectLeft, eigVectRight);
for (i = 1; i <= 5; ++i) {
std::cout << (scmatrix(a) - alpha[i] / beta[i] * scmatrix(b)).svd();
}
for (i = 1; i <= 5; ++i) {
std::cout << (~(eigVectLeft(i)) * scmatrix(a) - (alpha[i] / beta[i]) * ~(eigVectLeft(i)) * scmatrix(b)).norm() << std::endl;
}
for (i = 1; i <= 5; ++i) {
std::cout << (scmatrix(a) * eigVectRight(i) - (alpha[i] / beta[i]) * scmatrix(b) * eigVectRight(i)).norm() << std::endl;
}
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
99.7613 88.2526 59.4809 39.999 2.58055e-015
99.7613 88.2526 59.4809 39.999 2.58055e-015
26.8437 17.4907 9.30819 6.88691 1.44137e-015
25.8221 19.8697 11.2967 4.08824 3.35041e-015
25.8221 19.8697 11.2967 4.08824 1.43823e-015
5.91177e-014
5.91177e-014
9.42906e-015
5.56442e-015
5.65233e-015
4.10838e-014
4.10838e-014
5.82212e-015
9.46038e-015
9.15244e-015
See Also
cvector::eig()
Parameters
[in]mAsrmatrix $A$.
[in]mBsrmatrix $B$.
[out]vBetacomputed $\beta$ values.
[out]mEigVectLeftcomputed left eigenvectors $v_j$ stored by columns.
[out]mEigVectRightcomputed right eigenvectors $u_j$ stored by columns.
Returns
Reference to changed calling vector with computed $\alpha$ values

Definition at line 9328 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::geneig ( const basic_srmatrix< TR > &  mA,
const basic_srmatrix< TR > &  mB,
basic_rvector< TR > &  vBeta,
basic_scmatrix< TR, TC > &  mEigVect,
bool  bRightVect = true 
) throw (cvmexception)
inline

Generalized eigenvalues of real matrices.

Solves generalized eigenvalue problem and sets calling complex vector to be equal to generalized eigenvalues of real square matrices mA and mB. A generalized eigenvalue for a pair of matrices $(A,B)$ is a scalar $\lambda$ or a ratio $\alpha / \beta = \lambda$, such that $A - \lambda B$ is singular. It is usually represented as the pair $(\alpha, \beta)$, as there is a reasonable interpretation for $\beta =0$ and even for both being zero. The right generalized eigenvector $v_j$ corresponding to the generalized eigenvalue $\lambda_j$ of $(A,B)$ satisfies

\[ A v_j = \lambda_j B v_j. \]

The left generalized eigenvector $u_j$ corresponding to the generalized eigenvalue $\lambda_j$ of $(A,B)$ satisfies

\[ u_j^H A = \lambda_j u_j^H B. \]

Function returns reference to the vector changed and throws cvmexception in case of inappropriate calling object sizes or in case of convergence error.

Example:
using namespace cvm;
try {
srmatrix a(5), b(5);
a.randomize(-10., 10.);
b.randomize(-10., 10.);
cvector alpha(5);
rvector beta(5);
scmatrix eigVectLeft(5), eigVectRight(5);
int i;
alpha.geneig(a, b, beta, eigVectLeft, false);
for (i = 1; i <= 5; ++i) {
std::cout << (scmatrix(a) - alpha[i] / beta[i] * scmatrix(b)).svd();
}
for (i = 1; i <= 5; ++i) {
std::cout << (~(eigVectLeft(i)) * scmatrix(a) - (alpha[i] / beta[i]) * ~(eigVectLeft(i)) * scmatrix(b)).norm() << std::endl;
}
alpha.geneig(a, b, beta, eigVectRight);
for (i = 1; i <= 5; ++i) {
std::cout << (scmatrix(a) * eigVectRight(i) - (alpha[i] / beta[i]) * scmatrix(b) * eigVectRight(i)).norm() << std::endl;
}
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
240.642 220.334 163.141 92.5529 3.32163e-014
25.6052 17.872 12.2668 8.4116 1.95196e-015
24.2267 22.5764 12.2631 8.3846 7.13458e-016
24.2267 22.5764 12.2631 8.3846 8.2214e-016
20.9593 15.6649 14.9024 9.8735 2.25102e-016
1.3069e-013
1.39588e-014
1.07709e-014
1.10778e-014
1.07411e-014
8.41382e-014
1.46503e-014
1.78632e-014
1.80007e-014
4.56813e-015
See Also
cvector::eig()
Parameters
[in]mAsrmatrix $A$.
[in]mBsrmatrix $B$.
[out]vBetacomputed $\beta$ values.
[out]mEigVectcomputed left or right eigenvectors stored by columns.
[in]bRightVecttrue to compute right eigenvectors, false to compute left ones.
Returns
Reference to changed calling vector with computed $\alpha$ values

Definition at line 9418 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::geneig ( const basic_scmatrix< TR, TC > &  mA,
const basic_scmatrix< TR, TC > &  mB,
basic_cvector< TR, TC > &  vBeta 
) throw (cvmexception)
inline

Generalized eigenvalues of complex matrices.

Solves generalized eigenvalue problem and sets calling complex vector to be equal to generalized eigenvalues of complex square matrices mA and mB. A generalized eigenvalue for a pair of matrices $(A,B)$ is a scalar $\lambda$ or a ratio $\alpha / \beta = \lambda$, such that $A - \lambda B$ is singular. It is usually represented as the pair $(\alpha, \beta)$, as there is a reasonable interpretation for $\beta =0$ and even for both being zero. Function returns reference to the vector changed and throws cvmexception in case of inappropriate calling object sizes or in case of convergence error.

Example:
using namespace cvm;
try {
scmatrix a(5);
scbmatrix b(5, 2, 1);
a.randomize_real(-10., 10.);
a.randomize_imag(-10., 10.);
b.randomize_real(-10., 10.);
b.randomize_imag(-10., 10.);
cvector alpha(5);
cvector beta(5);
alpha.geneig(a, b, beta);
for (int i = 1; i <= 5; ++i) {
std::cout << (a - alpha[i] / beta[i] * b).svd();
}
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
54.1185 36.7609 17.7991 12.6586 2.06711e-015
37.0357 23.4901 14.3073 7.11309 6.83749e-015
43.3951 31.7299 26.8143 12.1402 9.2184e-015
30.6039 21.3084 14.9713 6.77913 1.25657e-015
29.5941 22.7112 18.7702 4.08204 4.24703e-015
See Also
cvector::eig()
Parameters
[in]mAscmatrix $A$.
[in]mBscmatrix $B$.
[out]vBetacomputed $\beta$ values.
Returns
Reference to changed calling vector with computed $\alpha$ values

Definition at line 9479 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::geneig ( const basic_scmatrix< TR, TC > &  mA,
const basic_scmatrix< TR, TC > &  mB,
basic_cvector< TR, TC > &  vBeta,
basic_scmatrix< TR, TC > &  mEigVectLeft,
basic_scmatrix< TR, TC > &  mEigVectRight 
) throw (cvmexception)
inline

Generalized eigenvalues of complex matrices.

Solves generalized eigenvalue problem and sets calling complex vector to be equal to generalized eigenvalues of complex square matrices mA and mB. A generalized eigenvalue for a pair of matrices $(A,B)$ is a scalar $\lambda$ or a ratio $\alpha / \beta = \lambda$, such that $A - \lambda B$ is singular. It is usually represented as the pair $(\alpha, \beta)$, as there is a reasonable interpretation for $\beta =0$ and even for both being zero. The right generalized eigenvector $v_j$ corresponding to the generalized eigenvalue $\lambda_j$ of $(A,B)$ satisfies

\[ A v_j = \lambda_j B v_j. \]

The left generalized eigenvector $u_j$ corresponding to the generalized eigenvalue $\lambda_j$ of $(A,B)$ satisfies

\[ u_j^H A = \lambda_j u_j^H B. \]

Function returns reference to the vector changed and throws cvmexception in case of inappropriate calling object sizes or in case of convergence error.

Example:
using namespace cvm;
try {
scmatrix a(5);
scbmatrix b(5, 2, 1);
a.randomize_real(-10., 10.);
a.randomize_imag(-10., 10.);
b.randomize_real(-10., 10.);
b.randomize_imag(-10., 10.);
cvector alpha(5);
cvector beta(5);
scmatrix eigVectLeft(5), eigVectRight(5);
int i;
alpha.geneig(a, b, beta, eigVectLeft, eigVectRight);
for (i = 1; i <= 5; ++i) {
std::cout << (a - alpha[i] / beta[i] * b).svd();
}
for (i = 1; i <= 5; ++i) {
std::cout << (~(eigVectLeft(i)) * a - (alpha[i] / beta[i]) * ~(eigVectLeft(i)) * b).norm() << std::endl;
}
for (i = 1; i <= 5; ++i) {
std::cout << (a * eigVectRight(i) - (alpha[i] / beta[i]) * b * eigVectRight(i)).norm() << std::endl;
}
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
143.311 87.1305 55.3041 33.3983 2.72025e-015
62.2127 52.6279 28.8434 16.4205 3.40966e-015
63.6394 45.3556 35.0669 11.7007 6.85443e-015
32.1333 23.8343 16.2436 9.18922 6.82323e-016
38.7917 31.4473 18.0545 10.806 7.54392e-015
4.65447e-014
3.99705e-014
2.64046e-014
1.54875e-014
1.3413e-014
3.93314e-014
1.06772e-014
2.37629e-014
1.14562e-014
1.99e-014
See Also
cvector::eig()
Parameters
[in]mAscmatrix $A$.
[in]mBscmatrix $B$.
[out]vBetacomputed $\beta$ values.
[out]mEigVectLeftcomputed left eigenvectors $v_j$ stored by columns.
[out]mEigVectRightcomputed right eigenvectors $u_j$ stored by columns.
Returns
Reference to changed calling vector with computed $\alpha$ values

Definition at line 9566 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::geneig ( const basic_scmatrix< TR, TC > &  mA,
const basic_scmatrix< TR, TC > &  mB,
basic_cvector< TR, TC > &  vBeta,
basic_scmatrix< TR, TC > &  mEigVect,
bool  bRightVect = true 
) throw (cvmexception)
inline

Generalized eigenvalues of complex matrices.

Solves generalized eigenvalue problem and sets calling complex vector to be equal to generalized eigenvalues of complex square matrices mA and mB. A generalized eigenvalue for a pair of matrices $(A,B)$ is a scalar $\lambda$ or a ratio $\alpha / \beta = \lambda$, such that $A - \lambda B$ is singular. It is usually represented as the pair $(\alpha, \beta)$, as there is a reasonable interpretation for $\beta =0$ and even for both being zero. The right generalized eigenvector $v_j$ corresponding to the generalized eigenvalue $\lambda_j$ of $(A,B)$ satisfies

\[ A v_j = \lambda_j B v_j. \]

The left generalized eigenvector $u_j$ corresponding to the generalized eigenvalue $\lambda_j$ of $(A,B)$ satisfies

\[ u_j^H A = \lambda_j u_j^H B. \]

Function returns reference to the vector changed and throws cvmexception in case of inappropriate calling object sizes or in case of convergence error.

Example:
using namespace cvm;
try {
scmatrix a(5);
scbmatrix b(5, 2, 1);
a.randomize_real(-10., 10.);
a.randomize_imag(-10., 10.);
b.randomize_real(-10., 10.);
b.randomize_imag(-10., 10.);
cvector alpha(5);
cvector beta(5);
scmatrix eigVectLeft(5), eigVectRight(5);
int i;
alpha.geneig(a, b, beta, eigVectLeft, false);
for (i = 1; i <= 5; ++i) {
std::cout << (a - alpha[i] / beta[i] * b).svd();
}
for (i = 1; i <= 5; ++i) {
std::cout << (~(eigVectLeft(i)) * a - (alpha[i] / beta[i]) * ~(eigVectLeft(i)) * b).norm() << std::endl;
}
alpha.geneig(a, b, beta, eigVectRight);
for (i = 1; i <= 5; ++i) {
std::cout << (a * eigVectRight(i) - (alpha[i] / beta[i]) * b * eigVectRight(i)).norm() << std::endl;
}
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
58.952 50.1947 28.6464 11.6468 2.4024e-015
51.8029 34.5341 26.584 13.8821 8.32486e-015
41.2836 33.7583 20.6533 11.6938 6.15973e-015
32.6075 20.7858 16.4413 10.8306 3.05618e-015
31.678 23.5424 14.7232 9.93361 2.41315e-015
2.13093e-014
3.15991e-014
1.97457e-014
7.83081e-015
8.88248e-015
3.70688e-014
1.89195e-014
2.79176e-014
5.49027e-015
9.92253e-015
See Also
cvector::eig()
Parameters
[in]mAscmatrix $A$.
[in]mBscmatrix $B$.
[out]vBetacomputed $\beta$ values.
[out]mEigVectcomputed left or right eigenvectors stored by columns.
[in]bRightVecttrue to compute right eigenvectors, false to compute left ones.
Returns
Reference to changed calling vector with computed $\alpha$ values

Definition at line 9659 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::gemv ( bool  bLeft,
const basic_cmatrix< TR, TC > &  m,
TC  dAlpha,
const basic_cvector< TR, TC > &  v,
TC  dBeta 
) throw (cvmexception)
inline

Generic matrix-vector operation.

Calls one of ZGEMV routines of the BLAS Library performing matrix-vector operation defined as

\[ c=\alpha\,M\cdot v + \beta\,c\quad\text{or}\quad c=\alpha\,v\cdot M + \beta\, c, \]

where $\alpha$ and $\beta$ are complex numbers (parameters dAlpha and dBeta), $M$ is matrix (parameter m) and $v$ and $c$ are vectors (parameter v and calling vector respectively). First operation is performed if bLeft passed is false and second one otherwise. Function returns reference to the vector changed and throws cvmexception in case of inappropriate sizes of the operands.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
try {
std::complex<double> alpha = std::complex<double>(1.3,-0.7);
std::complex<double> beta = std::complex<double>(0.15,-1.09);
cmatrix m(3,2);
cvector c(3);
cvector v(2);
m.randomize_real(-1., 2.); m.randomize_imag(0., 1.);
v.randomize_real(-1., 3.); v.randomize_imag(2., 4.);
c.randomize_real(0., 2.); c.randomize_imag(3., 7.);
std::cout << m * v * alpha + c * beta;
std::cout << c.gemv(false, m, alpha, v, beta);
std::cout << c * m * alpha + v * beta;
std::cout << v.gemv(true, m, alpha, c, beta);
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(2.71e-01,2.44e+00) (2.20e+01,7.16e+00) (-7.89e-01,2.45e+00)
(2.71e-01,2.44e+00) (2.20e+01,7.16e+00) (-7.89e-01,2.45e+00)
(5.92e+01,-1.47e+01) (3.54e+01,-3.14e+00)
(5.92e+01,-1.47e+01) (3.54e+01,-3.14e+00)
See Also
http://www.netlib.org/blas
Parameters
[in]bLeftIs left-sided multiplication.
[in]mcmatrix $M$.
[in]dAlphaNumber $\alpha$.
[in]vcvector $v$.
[in]dBetaNumber $\beta$.
Returns
Reference to changed calling vector.

Definition at line 9725 of file cvm.h.

Here is the call graph for this function:

template<typename TR, typename TC>
basic_cvector& basic_cvector< TR, TC >::gbmv ( bool  bLeft,
const basic_scbmatrix< TR, TC > &  m,
TC  dAlpha,
const basic_cvector< TR, TC > &  v,
TC  dBeta 
) throw (cvmexception)
inline

Generic band matrix-vector operation.

Calls one of ZGBMV routines of the BLAS Library performing matrix-vector operation defined as

\[ c=\alpha\,M\cdot v + \beta\,c\quad\text{or}\quad c=\alpha\,v\cdot M + \beta\, c, \]

where $\alpha$ and $\beta$ are complex numbers (parameters dAlpha and dBeta), $M$ is band matrix (parameter m) and $v$ and $c$ are vectors (parameter v and calling vector respectively). First operation is performed if bLeft passed is false and second one otherwise. Function returns reference to the vector changed and throws cvmexception in case of inappropriate sizes of the operands.

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
try {
std::complex<double> alpha = std::complex<double>(1.3,-0.7);
std::complex<double> beta = std::complex<double>(0.15,-1.09);
scbmatrix m(3,1,0);
cvector c(3);
cvector v(3);
m.randomize_real(-1., 2.); m.randomize_imag(0., 1.);
v.randomize_real(-1., 3.); v.randomize_imag(2., 4.);
c.randomize_real(0., 2.); c.randomize_imag(3., 7.);
std::cout << m * v * alpha + c * beta;
std::cout << c.gbmv(false, m, alpha, v, beta);
std::cout << c * m * alpha + v * beta;
std::cout << v.gbmv(true, m, alpha, c, beta);
}
catch (std::exception& e) {
std::cout << "Exception " << e.what () << std::endl;
}
prints
(3.73e+00,7.96e-01) (6.89e+00,1.07e+01) (2.16e+00,3.28e+00)
(3.73e+00,7.96e-01) (6.89e+00,1.07e+01) (2.16e+00,3.28e+00)
(3.11e+01,2.51e+01) (-4.93e+00,1.34e+01) (1.70e+00,3.93e+00)
(3.11e+01,2.51e+01) (-4.93e+00,1.34e+01) (1.70e+00,3.93e+00)
See Also
http://www.netlib.org/blas
Parameters
[in]bLeftIs left-sided multiplication.
[in]mscbmatrix $M$.
[in]dAlphaNumber $\alpha$.
[in]vcvector $v$.
[in]dBetaNumber $\beta$.
Returns
Reference to changed calling vector.

Definition at line 9786 of file cvm.h.

Here is the call graph for this function:

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

Randomizer (real part)

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

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
cvector v(3);
v.randomize_real(-2.,3.);
std::cout << v;
prints
(-4.93e-01,0.00e+00) (1.37e+00,0.00e+00) (-1.49e-01,0.00e+00)
Parameters
[in]dFromFirst limit.
[in]dToSecond limit.
Returns
Reference to changed calling vector.

Definition at line 9817 of file cvm.h.

Here is the call graph for this function:

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

Randomizer (imaginary part)

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

Example:
using namespace cvm;
std::cout.setf (std::ios::scientific | std::ios::left);
std::cout.precision (2);
cvector v(3);
v.randomize_imag(-2.,3.);
std::cout << v;
prints
(0.00e+00,-4.37e-01) (0.00e+00,-1.59e+00) (0.00e+00,2.42e+00)
Parameters
[in]dFromFirst limit.
[in]dToSecond limit.
Returns
Reference to changed calling vector.

Definition at line 9846 of file cvm.h.

Here is the call graph for this function:


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