Class for SU(N) matrix. More...
#include <sun_matrix.h>
Public Member Functions | |
| const SU & | make_unitary () |
| Makes the matrix unitary by orthogonalizing the rows. | |
| const SU & | fix_det () |
| Fix determinant. | |
| const SU & | reunitarize () |
| Reunitarize SU(N) matrix. | |
| const SU & | random (int nhits=16) |
| Generate random SU(N) matrix. | |
| Algebra< SU< N, T > > | project_to_algebra () const |
| Complex< T > | e (const int i, const int j) const |
| Standard array indexing operation for matrices and vectors. | |
| Complex< T > | operator[] (const int i) const |
| Indexing operation [] defined only for vectors. | |
| RowVector< m, Complex< T > > | row (int r) const |
| Return row of a matrix as a RowVector. | |
| void | set_row (int r, const RowVector< m, S > &v) |
| Set row of Matrix with RowVector if types are assignable. | |
| Vector< n, Complex< T > > | column (int c) const |
| Returns column vector as value at index c. | |
| void | set_column (int c, const Vector< n, S > &v) |
| get column of a matrix | |
| DiagonalMatrix< n, Complex< T > > | diagonal () |
| Return diagonal of square matrix. | |
| void | set_diagonal (const Vector< n, S > &v) |
| Set diagonal of square matrix to Vector which is passed to the method. | |
| const Array< n, m, Complex< T > > & | asArray () const |
| Cast Matrix to Array. | |
| const DiagonalMatrix< n, Complex< T > > & | asDiagonalMatrix () const |
| Cast Vector to DiagonalMatrix. | |
| SU< N, T > | operator- () const |
| Unary - operator. | |
| const auto & | operator+ () const |
| Unary + operator. | |
| bool | operator== (const Matrix< n2, m2, S > &rhs) const |
| Boolean operator == to determine if two matrices are exactly the same. | |
| bool | operator!= (const Matrix< n, m, S > &rhs) const |
| Boolean operator != to check if matrices are exactly different. | |
| auto & | operator+= (const Matrix_t< n, m, S, MT > &rhs) & |
| Add assign operator with matrix or scalar. | |
| auto & | operator+= (const DiagonalMatrix< n, S > &rhs) & |
| add and sub assign a DiagonalMatrix | |
| auto & | operator-= (const Matrix_t< n, m, S, MT > &rhs) & |
| Subtract assign operator with matrix or scalar. | |
| auto & | operator*= (const Matrix_t< m, p, S, MT > &rhs) & |
| Multiply assign scalar or matrix. | |
| auto & | operator*= (const DiagonalMatrix< m, S > &rhs) & |
| mult and divide assign a diagonal - cols must match diagonal matrix rows | |
| auto & | operator/= (const S rhs) & |
| Divide assign scalar. | |
| const auto & | fill (const S rhs) |
| Matrix fill. | |
| Rtype | transpose () const |
| Transpose of matrix. | |
| const RowVector< n, Complex< T > > & | transpose () const |
| Transpose of vector. | |
| SU< N, T > | conj () const |
| Returns complex conjugate of Matrix. | |
| Rtype | dagger () const |
| Hermitian conjugate of matrix. | |
| Rtype | adjoint () const |
| Adjoint of matrix. | |
| auto | abs () const |
| Returns absolute value of Matrix. | |
| Matrix< n, m, hila::arithmetic_type< Complex< T > > > | real () const |
| Returns real part of Matrix or Vector. | |
| Matrix< n, m, hila::arithmetic_type< Complex< T > > > | imag () const |
| Returns imaginary part of Matrix or Vector. | |
| Complex< T > | trace () const |
| Computes Trace for Matrix. | |
| hila::type_mul< Complex< T >, S > | mul_trace (const Matrix_t< p, q, S, MT > &rm) const |
| Multiply with given matrix and compute trace of result. | |
| hila::arithmetic_type< Complex< T > > | squarenorm () const |
| Calculate square norm - sum of squared elements. | |
| hila::arithmetic_type< Complex< T > > | norm () const |
| Calculate vector norm - sqrt of squarenorm. | |
| Complex< T > | max () const |
| Find max or min value - only for arithmetic types. | |
| R | dot (const Matrix< p, q, S > &rhs) const |
| Dot product. | |
| Matrix< n, p, R > | outer_product (const Matrix< p, q, S > &rhs) const |
| Outer product. | |
| SU< N, T > & | random () |
| dot with matrix - matrix | |
| SU< N, T > & | gaussian_random (base_type width=1.0) |
| Fills Matrix with gaussian random elements. | |
| SU< N, T > | permute_columns (const Vector< m, int > &permutation) const |
| permute columns of Matrix | |
| SU< N, T > | permute_rows (const Vector< n, int > &permutation) const |
| permute rows of Matrix | |
| SU< N, T > | permute (const Vector< N, int > &permutation) const |
| Permute elements of vector. | |
| SU< N, T > | sort (Vector< N, int > &permutation, hila::sort order=hila::sort::ascending) const |
| Sort method for Vector. | |
| void | mult_by_2x2_left (int p, int q, const Matrix_t< 2, 2, R, Mt > &M) |
| Multiply \( n \times m \)-matrix from the left by \( n \times m \) matrix defined by \( 2 \times 2 \) sub matrix. | |
| void | mult_by_2x2_right (int p, int q, const Matrix_t< 2, 2, R, Mt > &M) |
| Multiply \( n \times m \)-matrix from the right by \( n \times m \) matrix defined by \( 2 \times 2 \) sub matrix. | |
| Complex< T > | det_lu () const |
| following calculate the determinant of a square matrix det() is the generic interface, using laplace for small matrices and LU for large | |
| Complex< T > | det_laplace () const |
| Determinant of matrix, using Laplace method. | |
| Complex< T > | det () const |
| determinant function - if matrix size is < 5, use Laplace, otherwise LU | |
| int | eigen_hermitean (DiagonalMatrix< n, Et > &eigenvalues, Matrix_t< n, n, Mt, MT > &eigenvectors, enum hila::sort sorted=hila::sort::unsorted) const |
| Calculate eigenvalues and -vectors of an hermitean (or symmetric) matrix. | |
| eigen_result< SU< N, T > > | eigen_hermitean (enum hila::sort sorted=hila::sort::unsorted) const |
| eigenvalues and -vectors of hermitean/symmetric matrix, alternative interface | |
| int | svd (Matrix_t< n, n, Mt, MT > &_U, DiagonalMatrix< n, Et > &_D, Matrix_t< n, n, Mt, MT > &_V, enum hila::sort sorted=hila::sort::unsorted) const |
| Singular value decomposition: divide matrix A = U S V*, where U,V unitary and S diagonal matrix of real singular values. Unpivoted Jacobi rotations. | |
| svd_result< SU< N, T > > | svd (enum hila::sort sorted=hila::sort::unsorted) const |
| svd and svd_pivot - alternative interface | |
| int | svd_pivot (Matrix_t< n, n, Mt, MT > &_U, DiagonalMatrix< n, Et > &_D, Matrix_t< n, n, Mt, MT > &_V, enum hila::sort sorted=hila::sort::unsorted) const |
| Singular value decomposition: divide matrix A = U S V*, where U,V unitary and S diagonal matrix of real singular values. Fully pivoted Jacobi rotations. | |
| std::string | str (int prec=8, char separator=' ') const |
Static Public Member Functions | |
| static constexpr bool | is_vector () |
| Returns true if Matrix is a vector. | |
| static constexpr bool | is_square () |
| Returns true if matrix is a square matrix. | |
| static constexpr int | rows () |
| Define constant methods rows(), columns() - may be useful in template code. | |
| static constexpr int | columns () |
| Returns column length. | |
| static constexpr int | size () |
| Returns size of Vector or square Matrix. | |
Public Attributes | |
| Complex< T > | c [n *m] |
| The data as a one dimensional array. | |
Detailed Description
class SU< N, T >
Class for SU(N) matrix.
Class for Special unitary group SU(N).
SU class is a special case inherited from Matrix_t class. SU specific or overloaded methods are:
- Template Parameters
-
N Dimensionality of SU(N) matrix T Arithmetic type of Complex<T> number which SU(N) matrix elements consist of
Definition at line 110 of file sun_matrix.h.
Member Function Documentation
◆ abs()
◆ adjoint()
◆ asArray()
◆ asDiagonalMatrix()
|
inlineinherited |
◆ column()
◆ columns()
◆ conj()
◆ dagger()
◆ det_laplace()
Determinant of matrix, using Laplace method.
Defined only for square matrices
For perfomance overloads exist for \( 1 \times 1 \), \( 2 \times 2 \) and \( 3 \times 3 \) matrices.
- Parameters
-
mat matrix to compute determinant for
- Returns
- T result determinant
Definition at line 1530 of file matrix_linalg.h.
◆ det_lu()
following calculate the determinant of a square matrix det() is the generic interface, using laplace for small matrices and LU for large
Matrix determinant with LU decomposition.
Algorithm: numerical Recipes, 2nd ed. p. 47 ff Works for Real and Complex matrices Defined only for square matrices
- Returns
- Complex<radix> Determinant
Definition at line 1528 of file matrix_linalg.h.
◆ diagonal()
|
inlineinherited |
◆ dot()
|
inlineinherited |
Dot product.
Only works between two Vector objects
- Template Parameters
-
p Row length for rhs q Column length for rhs S Type for rhs R Gives resulting type of lhs and rhs multiplication
- Parameters
-
rhs Vector to compute dot product with
- Returns
- R Value of dot product
◆ e()
|
inlineinherited |
Standard array indexing operation for matrices and vectors.
Accessing singular elements is insufficient, but matrix elements are often quite small.
Exammple for matrix:
Example for vector:
- Parameters
-
i row index j column index
- Returns
- T matrix element type
◆ eigen_hermitean() [1/2]
|
inherited |
Calculate eigenvalues and -vectors of an hermitean (or symmetric) matrix.
Calculate eigenvalues and vectors of an hermitean or real symmetric matrix.
Returns the number of Jacobi iterations, or -1 if did not converge. Arguments will contain eigenvalues and eigenvectors in columns of the "eigenvectors" matrix. Computation is done in double precision despite the input matrix types.
- Parameters
-
eigenvaluevec Vector of computed eigenvalue eigenvectors Eigenvectors as columns of \( n \times n \) Matrix
Algorithm uses fully pivoted Jacobi rotations.
Two interfaces:
H.eigen_hermitean(E, U, [optional: sort]); E: output is DiagnoalMatrix containing real eigenvalues U: nxn unitary matrix, columns are normalized eigenvectors
auto res = H.eigen_hermitean([optional: sort]); This saves the trouble of defining E and U (see below)
Thus, H = U E U^* and H U = U E
Both interfaces allow optional sorting according to eigenvalues: hila::sort::unsorted [default] hila::sort::ascending / descending
Example:
- Note
- The matrix has to be hermitean/symmetric
- Parameters
-
E diagonal matrix of real eigenvalues U Unitary nxn matrix of eigenvectors sorted sorting of eigenvalues (default:unsorted)
- Returns
- int number of jacobi rotations
Definition at line 1546 of file matrix_linalg.h.
◆ eigen_hermitean() [2/2]
|
inherited |
eigenvalues and -vectors of hermitean/symmetric matrix, alternative interface
result = eigen_hermitean(hila::sort [optional]) returns struct eigen_result<Mtype>, with fields eigenvalues: DiagonalMatrix of eigenvalues eigenvectors: nxn unitary matrix of eigenvectors
Example:
This interface saves the trouble of defining the eigenvalue and -vector variables.
Definition at line 1550 of file matrix_linalg.h.
◆ fill()
Matrix fill.
Fills the matrix with element if it is assignable to matrix type T
Works as follows:
- Template Parameters
-
S Element type to be assigned
- Parameters
-
rhs Element to fill matrix with
- Returns
- const Mtype&
◆ fix_det()
Fix determinant.
Set the determinant of the SU(N) matrix to 1
- Returns
- const SU&
Definition at line 174 of file sun_matrix.h.
◆ gaussian_random()
◆ imag()
◆ is_square()
◆ is_vector()
◆ make_unitary()
Makes the matrix unitary by orthogonalizing the rows.
This is not the most optimal method, but it is simple:
Let SU<N> H be a special unitary matrix and define the indicies \( i \in
\{0,...,N-1\} \). The method is as follows:
- Normalize
H.row(i) - Make rows
H.row(i+1)...H.row(n-1)orthogonal with respect to rowH.row(i)
- Returns
- const SU&
Definition at line 137 of file sun_matrix.h.
◆ mul_trace()
|
inlineinherited |
Multiply with given matrix and compute trace of result.
Slightly cheaper operation than
- Template Parameters
-
p q S MT
- Parameters
-
rm
- Returns
- hila::type_mul<T, S>
◆ mult_by_2x2_left()
|
inlineinherited |
Multiply \( n \times m \)-matrix from the left by \( n \times m \) matrix defined by \( 2 \times 2 \) sub matrix.
The multiplication is defined as follows, let \(M\) as the \( 2 \times 2 \) input matrix and \(B\) be (this) matrix, being the matrix stored in the object this method is called for. Let \(A = I\) be a \( n \times m \) unit matrix. We then set the values of A to be:
\begin{align} A_{p,p} = M_{0,0}, \hspace{5px} A_{p,q} = M_{0,1}, \hspace{5px} A_{q,p} = M_{1,0}, \hspace{5px} A_{q,q} = M_{1,1}. \end{align}
Then the resulting matrix will be:
\begin{align} B = A \cdot B \end{align}
- Template Parameters
-
R Element type of M Mt Matrix type of M
- Parameters
-
p First row and column q Second row and column M \( 2 \times 2\) Matrix to multiply with
◆ mult_by_2x2_right()
|
inlineinherited |
Multiply \( n \times m \)-matrix from the right by \( n \times m \) matrix defined by \( 2 \times 2 \) sub matrix.
See Matrix::mult_by_2x2_left, only difference being that the multiplication is from the right.
- Template Parameters
-
R Element type of M Mt Matrix type of M
- Parameters
-
p First row and column q Second row and column M \( 2 \times 2\) Matrix to multiply with
◆ norm()
◆ operator!=()
|
inlineinherited |
Boolean operator != to check if matrices are exactly different.
if matrices are exactly the same then this will return false
- Template Parameters
-
S Type for MAtrix which is being compared to
- Parameters
-
rhs right hand side Matrix which we are comparing
- Returns
- true
- false
◆ operator*=()
|
inlineinherited |
Multiply assign scalar or matrix.
Multiplication works as defined for matrix multiplication and scalar matrix multiplication.
Matrix multiply assign only defined for square matrices, since the matrix dimensions would change otherwise.
Multiply assign operator can be used in the following ways
Multiply assign matrix:
Multiply assign scalar:
- Parameters
-
rhs Matrix to multiply with
- Returns
- template <int p, typename S, typename MT, std::enable_if_t<hila::is_assignable<T &, hila::type_mul<T, S>>::value, int> = 0>&
◆ operator+()
◆ operator+=() [1/2]
|
inlineinherited |
add and sub assign a DiagonalMatrix
This is possible only for square matrices
◆ operator+=() [2/2]
|
inlineinherited |
Add assign operator with matrix or scalar.
Add assign operator can be used in the following ways
Add assign matrix:
Add assign scalar:
Adds scalar \( a \) to square matrix as \( M + a\cdot\mathbb{1} \)
- Template Parameters
-
S Element type of rhs MT Matrix type of rhs
- Parameters
-
rhs Matrix to add
- Returns
- Mtype&
◆ operator-()
Unary - operator.
casting from one Matrix (number) type to another: do not do this automatically. but require an explicit cast operator. This makes it easier to write code. or should it be automatic? keep/remove explicit? TODO: CHECK AVX CONVERSIONS
Returns matrix with the signs of all the elements in the Matrix flipped.
- Returns
- Mtype
◆ operator-=()
|
inlineinherited |
Subtract assign operator with matrix or scalar.
Subtract assign operator can be used in the following ways
Subtract assign matrix:
Subtract assign scalar:
Subtract scalar \( a \) to square matrix as \( M - a\cdot\mathbb{1} \)
- Parameters
-
rhs Matrix to subtract with
- Returns
- template <typename S, typename MT, std::enable_if_t<hila::is_assignable<T &, S>::value, int> = 0>&
◆ operator/=()
Divide assign scalar.
Divide works as defined for scalar matrix division.
Division assign operator can be used in the following ways
Divide assign scalar:
- Parameters
-
rhs Matrix to divide with
- Returns
- template <int p, typename S, typename MT, std::enable_if_t<hila::is_assignable<T &, hila::type_mul<T, S>>::value, int> = 0>&
◆ operator==()
|
inlineinherited |
◆ operator[]()
|
inlineinherited |
Indexing operation [] defined only for vectors.
Example:
- Template Parameters
-
q row size n p column size m
- Parameters
-
i row or vector index depending on which is being indexed
- Returns
- T
◆ outer_product()
|
inlineinherited |
Outer product.
Only works between two Vector objects
- Template Parameters
-
p Row length for rhs q Column length for rhs S Element type for rhs R Type between lhs and rhs multiplication
- Parameters
-
rhs Vector to compute outer product with
- Returns
- Matrix<n, p, R>
◆ permute()
◆ permute_columns()
◆ permute_rows()
◆ project_to_algebra()
Project matrix to antihermitean and traceless algebra of the group. suN generators, normalized as Tr(\lambda_n \lambda_m) = -1/2 \delta_nm , or equivalently: Tr(\lambda^{\dagger}_n \lambda_m) = 1/2 \delta_nm
off-diagonal are just that: \lambda^od_nm = i/2 for elements nm and mn \lambda^od_nm = 1/2 for nm; -1/2 for mn
diagonals: su(N) has N-1 diag generators parametrize these recursively: \lambda_1 = diag(i,-i,0,0,..)/sqrt(1)/2 \lambda_2 = diag(i,i,-2i,0,..)/sqrt(3)/2 \lambda_3 = diag(i,i,i,-3i,0,..)/sqrt(6)/2 \lambda_k = diag(i,.. i,-ki,0,..)/sqrt(k(k+1)/2)/2 .. \lambda_N-1 = diag(i,.. i,-(N-1)i)/sqrt(N(N-1)/2)/2
Define \lambda's so that diagonals come first
Dividing U = U_ah + U_h, U_ah = 1/2 (U - U^+) = a_k \lambda_k + tr.im I/N => Tr(\lambda_k U_ah) = -1/2 a_k = 1/2 (\lambda_k)_lm (u_ml - u_lm^*) => a_k = - (\lambda_k)_lm (u_ml - u_lm^*)
Thus, for diags, a_k = (u_00 + ..u_(k-1)(k-1) - k*u_kk).im 2/(sqrt(2k(k+1)))
and off-diags: symm: a_i = -i (u_kj - u_kj^* + u_jk - u_jk^*)/2 = -i i(u_kj.im + u_jk.im) = (u_kj.im + u_jk.im) antisymm: a_i = -i (i u_kj - i u_jk^* - i u_jk + i u_kj^*)/2 = (u_kj.re - u_jk.re)
Definition at line 281 of file sun_matrix.h.
◆ random() [1/2]
◆ random() [2/2]
Generate random SU(N) matrix.
If N=2 random SU(N) matrix is generated by using Pauli matrix representation.
If N > 2 then first we generate a random SU(2) matrix, and multiply the Identity matrix I(N) by this matrix using SU::mult_by_2x2_left. After this we SU::reunitarize the matrix due to numerical errors in the multiplication.
- Parameters
-
nhits Number of times I(N) matrix is multiplied to generate random SU(N) matrix, only relevant for N > 2s
- Returns
- const SU&
Definition at line 211 of file sun_matrix.h.
◆ real()
◆ reunitarize()
Reunitarize SU(N) matrix.
Steps to reunitarize are:
- Make SU(N) matrix unitary
- Fix determinant
- Returns
- const SU&
Definition at line 193 of file sun_matrix.h.
◆ row()
◆ rows()
◆ set_column()
◆ set_diagonal()
|
inlineinherited |
Set diagonal of square matrix to Vector which is passed to the method.
If called for non square matrix the program will throw an error.
Example:
- Template Parameters
-
S type vector to assign values to
- Parameters
-
v Vector to assign to diagonal
◆ set_row()
◆ size()
◆ sort()
|
inlineinherited |
Sort method for Vector.
Two interfaces: first returns permutation vector, which can be used to permute other vectors/matrices second does only sort
Direct sort:
Permutation vector:
◆ squarenorm()
◆ str()
◆ svd() [1/2]
|
inherited |
svd and svd_pivot - alternative interface
svd(hila::sort [optional]) returns struct svd_result<Mtype>, with fields U: nxn unitary matrix singularvalues: DiagonalMatrix of singular values (real, > 0) V: nxn unitary matrix
auto res = M.svd(); now res.U : nxn unitary res.singularvalues : DiagonalMatrix of singular values res.V : nxn unitary
result satifies M = res.U res.singularvalues res.V.dagger()
Definition at line 1571 of file matrix_linalg.h.
◆ svd() [2/2]
|
inherited |
Singular value decomposition: divide matrix A = U S V*, where U,V unitary and S diagonal matrix of real singular values. Unpivoted Jacobi rotations.
Unpivoted rotation is generally faster than pivoted (esp. on gpu), but a bit less accurate
Use: M.svd(U, S, V, [sorted]);
Result satisfies M = U S V.dagger()
The algorithm works by one-sided Jacobi rotations.
- U = V = 1
- Loop
- for i=0..(n-2); j=(i+1)..(n-1)
- calculate B_ii, B_ij, B_jj, (B_ij = (A^* A)_ij)
- If B_ij > 0
- diagonalize J^* [ B_ii B_ij ] J = diag(D_ii, D_jj) [ B_ji B_jj ]
- A <- A J , Columns i and j change
- V <- V J
- endfor
- if nothing changed, break Loop
- for i=0..(n-2); j=(i+1)..(n-1)
- Endloop Now A = U S, A^* A = S* S = S^2
- S_ii = A_ii / |A_i| (norm of column i)
- U = A / S
Definition at line 1567 of file matrix_linalg.h.
◆ svd_pivot()
|
inherited |
Singular value decomposition: divide matrix A = U S V*, where U,V unitary and S diagonal matrix of real singular values. Fully pivoted Jacobi rotations.
Use: M.svd_pivot(U, S, V, [sorted]);
The algorithm works by one-sided Jacobi rotations.
- U = V = 1
- Compute Gram matrix B = A^* A. B is Hermitean, with B_ii >= 0.
- Loop:
- Find largest B_ij, j > i - full pivoting if |B_ij|^2 <= tol * B_ii B_jj we're done
- find J to diagonalize J^* [ B_ii B_ij ] J = diag(D_ii, D_jj) [ B_ji B_jj ]
- A <- A J , Columns i and j change
- V <- V J
- Recompute B_ik and B_jk, k /= i,j (B_ii = D_ii, B_jj = D_jj, B_ij = 0) (n^2)
- Endloop Now A = U S, A^* A = S* S = S^2
- S_ii = A_ii / |A_i| (norm of column i)
- U = A / S
- Parameters
-
U
Definition at line 1576 of file matrix_linalg.h.
◆ trace()
◆ transpose() [1/2]
◆ transpose() [2/2]
The documentation for this class was generated from the following file:
