Those functions are available both for dense and sparse vectors. In the case of dense vectors, Blas subroutines are called if SELDON_WITH_BLAS
(or, for backward compatibility, SELDON_WITH_CBLAS
) is defined.
Mlt | multiplies the elements of the vector/matrix by a scalar, or performs a matrix-vector product |
MltAdd | performs a matrix-vector product or a matrix-matrix product |
Add | adds two vectors or two matrices |
Copy | copies a vector into another one, or a matrix into another matrix |
Swap | exchanges two vectors |
SwapPointer | exchanges two vectors (only adresses are exchanged) |
GenRot | computes givens rotation |
GenModifRot | computes givens roation |
ApplyRot | applies givens rotation |
ApplyModifRot | applies givens rotation |
DotProd | scalar product between two vectors |
DotProdConj | scalar product between two vectors, first vector being conjugated |
Conjugate | conjugates all elements of a vector |
GetMaxAbsIndex | returns index where highest absolute value is reached |
Norm1 | returns 1-norm of a vector |
Norm2 | returns 2-norm of a vector |
Rank1Update | Adds a contribution X.Y' to a matrix |
Rank2Update | Adds a contribution X.Y' + Y.X' to a symmetric matrix |
Solve | solves a triangular system |
void Mlt(T0, Vector&); void Mlt(T0, Matrix&); void Mlt(T0, Array3D&);
void Mlt(const Matrix&, const Vector&, Vector&); void Mlt(SeldonTrans, const Matrix&, const Matrix&, Matrix&); void Mlt(const Matrix&, const Matrix&, Matrix&);
This function multiplies all elements of vector/matrix/array by a constant. It can also perform a matrix-vector product. This function is implemented both for dense and sparse vectors/matrices.
Vector<double> X(3); X.Fill(); // we multiply all elements by 2 Mlt(2.0, X); // similar method for matrix Matrix<double> A(3, 3); A.Fill(); // we multiply all elements by -3.4 Mlt(-3.4, A); // similar stuff for a 3-D array Array3D<double> C(3, 4, 6); C.Fill(); Mlt(-2.5, C); // and sparse matrix as well Matrix<double, General, ArrayRowSparse> Asp(3,3), Bsp(3, 3), Csp(3, 3); Asp(0, 1) = 1.4; Asp(1, 1) = -0.3; Asp(2, 0) = 2.1; Mlt(2.1, Asp); // with Mlt, you can compute Y = Asp*X or Y = A*X Vector<double> Y(3); Mlt(A, X, Y); Mlt(Asp, X, Y); // matrix-matrix product // available for dense and sparse matrices Mlt(Asp, Bsp, Csp); // you can multiply with the transpose matrix // Y = A^T x Mlt(SeldonTrans, A, X, Y); // for triangular matrices, the result is overwritten in X // X = T X Matrix<double, General, RowUpTriangPacked> T(3, 3); Mlt(T, X); // multiplication with transpose Mlt(SeldonTrans, SeldonNonUnit, T, X); // multiplication with T, assuming that diagonal of T is equal to 1 X = transpose(T) X Mlt(SeldonNoTrans, SeldonUnit, T, X); // multiplication triangular/ dense matrices // M = T M Matrix<double> M(3, 5); Mlt(T, M);
Vector operators
Matrix operators
MltAdd
Functions_Vector.cxx
Functions_MatVect.cxx
Functions_Matrix.cxx
Functions_MatrixArray.cxx
Blas_1.cxx
Blas_2.cxx
Blas_3.cxx
void MltAdd(T0, const Matrix&, const Vector&, T0, Vector&); void MltAdd(T0, SeldonTrans, const Matrix&, const Vector&, T0, Vector&); void MltAdd(T0, const Matrix&, const Matrix&, T0, Matrix&); void MltAdd(T0, SeldonTrans, const Matrix&, SeldonTrans, const Matrix&, T0, Matrix&);
This function performs a matrix-vector product or a matrix-matrix product. You can specify a product with the transpose of the matrix, or the conjugate of the transpose. The size of matrices need to be compatible in order to complete the matrix-matrix product.
// matrix-vector product Matrix<double> A(3,3); Vector<double> X(3), B(3); A.Fill(); X.Fill(); B.Fill(1.5); double beta = 2, alpha = 3; // computation of B = beta*B + alpha*A*X MltAdd(alpha, A, X, beta, B); // same function for sparse matrices Matrix<double, General, RowSparse> Asp(3, 3, nnz, values, ptr, ind); // computation of B = beta*B + alpha*Asp*X MltAdd(alpha, Asp, X, beta, B); // you can compute B = beta*B + alpha*transpose(Asp)*X MltAdd(alpha, SeldonTrans, Asp, X, beta, B); // similar method for matrix-matrix product // available for sparse and dense matrices Matrix<double> M(3,4), N(4,3); // computation of A = beta*A + alpha*M*N MltAdd(alpha, M, N, beta, A); // and you can use SeldonTrans, SeldonNoTrans, SeldonConjTrans // computation of A = beta*A + alpha*transpose(M)*N MltAdd(alpha, SeldonTrans, M, SeldonNoTrans, N, beta, A); // SeldonConjTrans is meaningfull for complex matrices Matrix<complex<double> > Ac(3, 3), Bc(3, 3), Cc(3, 3); // computation of Ac = beta * Ac // + alpha * conjugate(transpose(Bc)) * transpose(Cc) MltAdd(alpha, SeldonConjTrans, Bc, SeldonTrans, Cc, beta, Ac);
Vector operators
Matrix operators
Mlt
Functions_MatVect.cxx
Functions_Matrix.cxx
Functions_MatrixArray.cxx
Blas_2.cxx
Blas_3.cxx
void Add(T0, const Vector&, Vector&); void Add(T0, const Matrix&, Matrix&);
This function adds two vectors or two matrices. This function is available both for sparse or dense vectors/matrices. The size of the vectors or the matrices to add should be the same.
Vector<double> X(3), Y(3); X.Fill(); Y.Fill(1); double alpha = 2; // computation of Y = Y + alpha*X Add(alpha, X, Y); // similar method for matrix Matrix<double> A(3,3), B(3,3); A.Fill(); B.SetIdentity(); Add(alpha, A, B); // and sparse matrix as well Matrix<double, General, ArrayRowSparse> Asp(3,3), Bsp(3,3); Asp(0, 1) = 1.4; Asp(1, 1) = -0.3; Asp(2, 0) = 2.1; Bsp.SetIdentity(); Add(alpha, Asp, Bsp);
Functions_Vector.cxx
Functions_Matrix.cxx
Functions_MatrixArray.cxx
Blas_1.cxx
void Copy(const Vector&, Vector&); void Copy(const Matrix& Matrix&);
A vector is copied into another one. If BLAS is enabled, the two dense vectors have to be of the same size. The operator = is more flexible. For sparse matrices, the function Copy
is overloaded to enable the conversion between the different formats.
// vectors need to have the same size Vector<float> V(3), W(3); V.Fill(); Copy(V, W); // it is equivalent to write W = V; W = V; // for sparse matrices, you can use Copy as a convert tool Matrix<double, General, ArrayRowSparse> A(3, 3); A.AddInteraction(0, 2, 2.5); A.AddInteraction(1, 1, -1.0); Matrix<double, General, RowSparse> B; Copy(A, B);
Functions_Vector.cxx
Blas_1.cxx
Matrix_Conversions.cxx
void Swap(Vector& Vector&);
The two vectors are exchanged. If BLAS is enabled, the two dense vectors have to be of the same size.
// vectors need to have the same size Vector<float> V(3), W(3); V.Fill(); W.Fill(1.0); Swap(V, W);
void SwapPointer(Vector& Vector&);
The two vectors are exchanged (only pointers are exchanged).
// vectors do not need to have the same size Vector<float> V(3), W(5); V.Fill(); W.Fill(1.0); Swap(V, W);
Functions_Vector.cxx
Blas_1.cxx
void GenRot(T&, T&, T&, T&); void ApplyRot(Vector&, Vector&, T, T); void GenModifRot(T&, T&, T&, T&, T*); void ApplyModifRot(Vector&, Vector&, T*);
This function constructs the Givens rotation G = [cos_, sin_; -sin_, cos_] that zeros the second entry of the two vector (a,b).
double a, b, cos_, sin_; a = 2.0; b = 1.0; // we compute cos_ and sin_ so that // G = [cos_, sin_; -sin_, cos_] zeros second entry of vector [a;b] // i.e. G*[a;b] = [r;0] GenRot(a, b, cos_, sin_); // then we can apply rotation G to vectors X, Y // G*[x_i, y_i] for all i Vector<double> X(10), Y(10); X.FillRand(); Y.FillRand(); ApplyRot(X, Y, cos_, sin_);
Functions_Vector.cxx
Blas_1.cxx
T DotProd(const Vector&, const Vector&); T DotProdConj(const Vector&, const Vector&);
This function returns the scalar product between two vectors. For DotProdConj
, the first vector is conjugated. In the case of dense vectors, they have to be of the same size.
// we construct two vectors of same size Vector<double> X(10), Y(10); X.FillRand(); Y.FillRand(); // computation of X*Y double scal = DotProd(X, Y); // for complex numbers Vector<complex<double> > Xc(10), Yc(10); Xc.Fill(); Yc.Fill(); // computation of conj(X)*Y complex<double> scalc = DotProdConj(X, Y);
Functions_Vector.cxx
Blas_1.cxx
void Conjugate(Vector&); void Conjugate(Matrix&);
Each component of the vector (or matrix) is conjugated. In the case of real numbers, the vector (or matrix) is not modified.
// complex vector Vector<complex<double> > X(10); X.FillRand(); // we take the conjugate Conjugate(X); // complex matrix Matrix<complex<double> > A(10, 8); A.FillRand(); // we take the conjugate Conjugate(A);
Functions_Vector.cxx
Blas_1.cxx
int GetMaxAbsIndex(const Vector&);
This function returns the index for which the vector reached its highest absolute value.
// complex vector Vector<complex<double> > X(10); X.FillRand(); int imax = GetMaxAbsIndex(X); // maximum ? double maximum = abs(X(imax));
Functions_Vector.cxx
Blas_1.cxx
T Norm1(const Vector&);
This function returns the 1-norm of the vector.
// complex vector Vector<complex<double> > X(10); X.Fill(); // 1-norm double norme = Norm1(X);
Functions_Vector.cxx
Blas_1.cxx
T Norm2(const Vector&);
This function returns the 2-norm of the vector.
// complex vector Vector<complex<double> > X(10); X.Fill(); // 2-norm double norme = Norm2(X);
Functions_Vector.cxx
Blas_1.cxx
void Rank1Update(T, const Vector&, const Vector&, Matrix&); void Rank1Update(T, const Vector&, SeldonConj, const Vector&, Matrix&); void Rank1Update(T, const Vector&, Matrix&);
This function adds a contribution of rank 1 to matrix.
// two complex vectors Vector<complex<double> > X(10), Y(10); X.FillRand(); Y.FillRand(); // complex matrix Matrix<complex<double> > A(10, 10); A.SetIdentity(); complex<double> alpha(1, 1); // we compute A = A + alpha*X*transpose(Y) Rank1Update(alpha, X, Y, A); // we can also compute A = A + alpha*X*conjugate(transpose(Y)) Rank1Update(alpha, X, SeldonConj, Y, A); // for hermitian matrices, use of X only Matrix<complex<double>, General, RowHermPacked> B(10, 10); // we compute B = B + alpha*X*conjugate(transpose(X)) Rank1Update(alpha, X, B); // for real numbers, use of symmetric matrices Vector<double> Xr(10); Xr.FillRand(); Matrix<double, Symmetric, RowSymPacked> Br(10, 10); Br.Zero(); double beta = 2; // computation of Br = Br + beta*Xr*tranpose(Xr) Rank1Update(beta, Xr, Br);
void Rank2Update(T, const Vector&, const Vector&, Matrix&);
This function adds a contribution of rank 2 to a symmetric or Hermitian matrix.
// two complex vectors Vector<complex<double> > X(10), Y(10); X.FillRand(); Y.FillRand(); // hermitian matrix Matrix<complex<double>, General, RowHermPacked> A(10, 10); A.SetIdentity(); complex<double> alpha(1, 1); // we compute A = A + alpha * (X * conjugate(transpose(Y)) // + Y * conjugate(transpose(X))) Rank2Update(alpha, X, Y, A); // for real numbers, use of symmetric matrices Vector<double> Xr(10), Yr(10); Xr.FillRand(); Yr.FillRand(); Matrix<double, Symmetric, RowSymPacked> Br(10, 10); Br.Zero(); double beta = 2; // computation of Br = Br + beta*(Xr*transpose(Yr) + Yr*transpose(Xr)) Rank2Update(beta, Xr, Yr, Br);
void Solve(const Matrix&, Vector&); void Solve(SeldonTrans, SeldonNonUnit, const Matrix&, Vector&);
This function solves a triangular linear system.
// right hand side Vector<double> B(3); B.Fill(); // triangular matrix Matrix<double, General, RowUpTriangPacked> T(3, 3); T.Fill(); // now you can solve T*X = B, the result overwrites B Solve(T, B); // if you want to solve transpose(T)*X = B Solve(SeldonTrans, SeldonNonUnit, T, B); // SeldonUnit can be used if you are assuming that T has // a unitary diagonal (with 1) // we force unitary diagonal for (int i = 0; i < T.GetM(); i++) T(i, i) = 1.0; // then we can call Solve with SeldonUnit to solve T*X = B // B can be also a dense matrix Solve(SeldonNoTrans, SeldonUnit, T, B);