Documentation for module cp_fm_basic_linalg

basic linear algebra operations for full matrices

source: cp_fm_basic_linalg.F
Loading...

Generic procedures:

cp_fm_contracted_trace
cp_fm_trace

public Subroutines/Functions:

Convenience function. Computes the matrix multiplications needed for the multiplication of complex matrices. C = beta * C + alpha * ( A ** transa ) * ( B ** transb )
scales column i of matrix a with scaling(i)
computes the Frobenius norm of matrix_a
computes matrix_c = beta * matrix_c + alpha * ( matrix_a ** transa ) * ( matrix_b ** transb )
Inverts a cp_fm_type matrix, optionally returning the determinant of the input matrix
FUNCTION
REAL(dp)
cp_fm_norm (matrix, mode)
norm of matrix using (p)dlange
perfoms a QR factorization of the input rectangular matrix A or of a submatrix of A the computed upper triangular matrix R is in output in the submatrix sub(A) of size NxN M and M give the dimention of the submatrix that has to be factorized (MxN) with M>N
SUBROUTINE
cp_fm_scale (alpha, matrix_a)
scales a matrix matrix_a = alpha * matrix_b
calc A <- alpha*A + beta*B optimized for alpha == 1.0 (just add beta*B) and beta == 0.0 (just scale A)
computes the schur product of two matrices c_ij = a_ij * b_ij
computs the the solution to A*b=A_general using lu decomposition pay attention, both matrices are overwritten, a_general contais the result
computes matrix_c = beta * matrix_c + alpha * matrix_a * matrix_b computes matrix_c = beta * matrix_c + alpha * matrix_b * matrix_a where matrix_a is symmetric
performs a rank-k update of a symmetric matrix_c matrix_c = beta * matrix_c + alpha * matrix_a * transpose ( matrix_a )
transposes a matrix matrixt = matrix ^ T
inverts a triangular matrix
multiplies in place by a triangular matrix: matrix_b = alpha op(triangular_matrix) matrix_b or (if side='R') matrix_b = alpha matrix_b op(triangular_matrix) op(triangular_matrix) is: triangular_matrix (if transpose_tr=.false. and invert_tr=.false.) triangular_matrix^T (if transpose_tr=.true. and invert_tr=.false.) triangular_matrix^(-1) (if transpose_tr=.false. and invert_tr=.true.) triangular_matrix^(-T) (if transpose_tr=.true. and invert_tr=.true.)
given an upper triangular matrix computes the corresponding full matrix

Generic procedure cp_fm_contracted_trace

Compute trace(i,j) = \sum_k Tr (matrix_a(k,i)^T matrix_b(k,j)).

INTENT(in)
:: matrix_a(:,:) ×
INTENT(in)
:: matrix_b(:,:) ×
REAL(dp),
INTENT(out)
:: trace(:,:) ×

Generic procedure cp_fm_trace

...

TYPE(cp_fm_type),
POINTER
:: matrix_a ×
INTENT(in)
:: matrix_a(:) ×
TYPE(cp_fm_type),
POINTER
:: matrix_b ×
INTENT(in)
:: matrix_b(:) ×
REAL(dp),
INTENT(out)
:: trace ×
REAL(dp),
INTENT(out)
:: trace(:) ×

SUBROUTINEcp_complex_fm_gemm(transa, transb, m, n, k, alpha, a_re, a_im, b_re, b_im, beta, c_re, c_im, a_first_col, a_first_row, b_first_col, b_first_row, c_first_col, c_first_row)

Convenience function. Computes the matrix multiplications needed for the multiplication of complex matrices. C = beta * C + alpha * ( A ** transa ) * ( B ** transb )

Arguments:
CHARACTER(1),
INTENT(in)
:: transa 'N' -> normal 'T' -> transpose alpha,beta :: can be 0.0_dp and 1.0_dp
CHARACTER(1),
INTENT(in)
:: transb ...
INTEGER,
INTENT(in)
:: m ...
INTEGER,
INTENT(in)
:: n ...
INTEGER,
INTENT(in)
:: k ...
REAL(dp),
INTENT(in)
:: alpha ...
TYPE(cp_fm_type),
POINTER
:: a_re m x k matrix ( ! for transa = 'N'), real part
TYPE(cp_fm_type),
POINTER
:: a_im m x k matrix ( ! for transa = 'N'), imaginary part
TYPE(cp_fm_type),
POINTER
:: b_re k x n matrix ( ! for transa = 'N'), real part
TYPE(cp_fm_type),
POINTER
:: b_im k x n matrix ( ! for transa = 'N'), imaginary part
REAL(dp),
INTENT(in)
:: beta ...
TYPE(cp_fm_type),
POINTER
:: c_re m x n matrix, real part
TYPE(cp_fm_type),
POINTER
:: c_im m x n matrix, imaginary part
INTEGER,
INTENT(in),
OPTIONAL
:: a_first_col ...
INTEGER,
INTENT(in),
OPTIONAL
:: a_first_row ...
INTEGER,
INTENT(in),
OPTIONAL
:: b_first_col the k x n matrix starts at col b_first_col of matrix_b (avoid usage)
INTEGER,
INTENT(in),
OPTIONAL
:: b_first_row ...
INTEGER,
INTENT(in),
OPTIONAL
:: c_first_col ...
INTEGER,
INTENT(in),
OPTIONAL
:: c_first_row ...

SUBROUTINEcp_fm_column_scale(matrixa, scaling)

scales column i of matrix a with scaling(i)

Arguments:
TYPE(cp_fm_type),
POINTER
:: matrixa ...
REAL(dp),
INTENT(in)
:: scaling(:) an array used for scaling the columns, SIZE(scaling) determines the number of columns to be scaled

SUBROUTINEcp_fm_frobenius_norm(matrix_a, norm)

computes the Frobenius norm of matrix_a

Arguments:
TYPE(cp_fm_type),
POINTER
:: matrix_a m x n matrix
REAL(dp),
INTENT(inout)
:: norm ...

SUBROUTINEcp_fm_gemm(transa, transb, m, n, k, alpha, matrix_a, matrix_b, beta, matrix_c, a_first_col, a_first_row, b_first_col, b_first_row, c_first_col, c_first_row)

computes matrix_c = beta * matrix_c + alpha * ( matrix_a ** transa ) * ( matrix_b ** transb )

Arguments:
CHARACTER(1),
INTENT(in)
:: transa 'N' -> normal 'T' -> transpose alpha,beta :: can be 0.0_dp and 1.0_dp
CHARACTER(1),
INTENT(in)
:: transb ...
INTEGER,
INTENT(in)
:: m ...
INTEGER,
INTENT(in)
:: n ...
INTEGER,
INTENT(in)
:: k ...
REAL(dp),
INTENT(in)
:: alpha ...
TYPE(cp_fm_type),
POINTER
:: matrix_a m x k matrix ( ! for transa = 'N')
TYPE(cp_fm_type),
POINTER
:: matrix_b k x n matrix ( ! for transb = 'N')
REAL(dp),
INTENT(in)
:: beta ...
TYPE(cp_fm_type),
POINTER
:: matrix_c m x n matrix
INTEGER,
INTENT(in),
OPTIONAL
:: a_first_col ...
INTEGER,
INTENT(in),
OPTIONAL
:: a_first_row ...
INTEGER,
INTENT(in),
OPTIONAL
:: b_first_col the k x n matrix starts at col b_first_col of matrix_b (avoid usage)
INTEGER,
INTENT(in),
OPTIONAL
:: b_first_row ...
INTEGER,
INTENT(in),
OPTIONAL
:: c_first_col ...
INTEGER,
INTENT(in),
OPTIONAL
:: c_first_row ...

SUBROUTINEcp_fm_invert(matrix_a, matrix_inverse, det_a, eps_svd)

Inverts a cp_fm_type matrix, optionally returning the determinant of the input matrix

Arguments:
TYPE(cp_fm_type),
POINTER
:: matrix_a the matrix to invert
TYPE(cp_fm_type),
POINTER
:: matrix_inverse the inverse of matrix_a
REAL(dp),
INTENT(out),
OPTIONAL
:: det_a the determinant of matrix_a
REAL(dp),
INTENT(in),
OPTIONAL
:: eps_svd optional parameter to active SVD based inversion, singular values below eps_svd are screened

FUNCTIONcp_fm_norm(matrix, mode)

norm of matrix using (p)dlange

Return Value ::
REAL(dp)
the norm according to mode
Arguments:
TYPE(cp_fm_type),
POINTER
:: matrix input a general matrix
CHARACTER,
INTENT(in)
:: mode 'M' max abs element value, '1' or 'O' one norm, i.e. maximum column sum 'I' infinity norm, i.e. maximum row sum 'F' or 'E' Frobenius norm, i.e. sqrt of sum of all squares of elements

SUBROUTINEcp_fm_qr_factorization(matrix_a, matrix_r, nrow_fact, ncol_fact, first_row, first_col)

perfoms a QR factorization of the input rectangular matrix A or of a submatrix of A the computed upper triangular matrix R is in output in the submatrix sub(A) of size NxN M and M give the dimention of the submatrix that has to be factorized (MxN) with M>N

Arguments:
TYPE(cp_fm_type),
POINTER
:: matrix_a ...
TYPE(cp_fm_type),
POINTER
:: matrix_r ...
INTEGER,
INTENT(in),
OPTIONAL
:: nrow_fact ...
INTEGER,
INTENT(in),
OPTIONAL
:: ncol_fact ...
INTEGER,
INTENT(in),
OPTIONAL
:: first_row ...
INTEGER,
INTENT(in),
OPTIONAL
:: first_col ...

SUBROUTINEcp_fm_scale(alpha, matrix_a)

scales a matrix matrix_a = alpha * matrix_b

Arguments:
REAL(dp),
INTENT(in)
:: alpha ...
TYPE(cp_fm_type),
POINTER
:: matrix_a ...

SUBROUTINEcp_fm_scale_and_add(alpha, matrix_a, beta, matrix_b)

calc A <- alpha*A + beta*B optimized for alpha == 1.0 (just add beta*B) and beta == 0.0 (just scale A)

Arguments:
REAL(dp),
INTENT(in)
:: alpha ...
TYPE(cp_fm_type),
POINTER
:: matrix_a ...
REAL(dp),
INTENT(in),
OPTIONAL
:: beta ...
TYPE(cp_fm_type),
OPTIONAL, POINTER
:: matrix_b ...

SUBROUTINEcp_fm_schur_product(matrix_a, matrix_b, matrix_c)

computes the schur product of two matrices c_ij = a_ij * b_ij

Arguments:
TYPE(cp_fm_type),
POINTER
:: matrix_a ...
TYPE(cp_fm_type),
POINTER
:: matrix_b ...
TYPE(cp_fm_type),
POINTER
:: matrix_c ...

SUBROUTINEcp_fm_solve(matrix_a, general_a)

computs the the solution to A*b=A_general using lu decomposition pay attention, both matrices are overwritten, a_general contais the result

Arguments:
TYPE(cp_fm_type),
POINTER
:: matrix_a ...
TYPE(cp_fm_type),
POINTER
:: general_a ...

SUBROUTINEcp_fm_symm(side, uplo, m, n, alpha, matrix_a, matrix_b, beta, matrix_c)

computes matrix_c = beta * matrix_c + alpha * matrix_a * matrix_b computes matrix_c = beta * matrix_c + alpha * matrix_b * matrix_a where matrix_a is symmetric

Arguments:
CHARACTER(1),
INTENT(in)
:: side 'L' -> matrix_a is on the left 'R' -> matrix_a is on the right alpha,beta :: can be 0.0_dp and 1.0_dp
CHARACTER(1),
INTENT(in)
:: uplo ...
INTEGER,
INTENT(in)
:: m ...
INTEGER,
INTENT(in)
:: n ...
REAL(dp),
INTENT(in)
:: alpha ...
TYPE(cp_fm_type),
POINTER
:: matrix_a m x m matrix
TYPE(cp_fm_type),
POINTER
:: matrix_b m x n matrix
REAL(dp),
INTENT(in)
:: beta ...
TYPE(cp_fm_type),
POINTER
:: matrix_c m x n matrix

SUBROUTINEcp_fm_syrk(uplo, trans, k, alpha, matrix_a, ia, ja, beta, matrix_c)

performs a rank-k update of a symmetric matrix_c matrix_c = beta * matrix_c + alpha * matrix_a * transpose ( matrix_a )

Arguments:
CHARACTER(1),
INTENT(in)
:: uplo 'U' ('L')
CHARACTER(1),
INTENT(in)
:: trans 'N' ('T')
INTEGER,
INTENT(in)
:: k number of cols to use in matrix_a ia,ja :: 1,1 (could be used for selecting subblock of a)
REAL(dp),
INTENT(in)
:: alpha ...
TYPE(cp_fm_type),
POINTER
:: matrix_a ...
INTEGER,
INTENT(in)
:: ia ...
INTEGER,
INTENT(in)
:: ja ...
REAL(dp),
INTENT(in)
:: beta ...
TYPE(cp_fm_type),
POINTER
:: matrix_c ...

SUBROUTINEcp_fm_transpose(matrix, matrixt)

transposes a matrix matrixt = matrix ^ T

Arguments:
TYPE(cp_fm_type),
POINTER
:: matrix ...
TYPE(cp_fm_type),
POINTER
:: matrixt ...

SUBROUTINEcp_fm_triangular_invert(matrix_a, uplo_tr)

inverts a triangular matrix

Arguments:
TYPE(cp_fm_type),
POINTER
:: matrix_a ...
CHARACTER,
INTENT(in),
OPTIONAL
:: uplo_tr ...

SUBROUTINEcp_fm_triangular_multiply(triangular_matrix, matrix_b, side, transpose_tr, invert_tr, uplo_tr, unit_diag_tr, n_rows, n_cols, alpha)

multiplies in place by a triangular matrix: matrix_b = alpha op(triangular_matrix) matrix_b or (if side='R') matrix_b = alpha matrix_b op(triangular_matrix) op(triangular_matrix) is: triangular_matrix (if transpose_tr=.false. and invert_tr=.false.) triangular_matrix^T (if transpose_tr=.true. and invert_tr=.false.) triangular_matrix^(-1) (if transpose_tr=.false. and invert_tr=.true.) triangular_matrix^(-T) (if transpose_tr=.true. and invert_tr=.true.)

Arguments:
TYPE(cp_fm_type),
POINTER
:: triangular_matrix the triangular matrix that multiplies the other
TYPE(cp_fm_type),
POINTER
:: matrix_b the matrix that gets multiplied and stores the result
CHARACTER,
INTENT(in),
OPTIONAL
:: side on which side of matrix_b stays op(triangular_matrix) (defaults to 'L')
LOGICAL,
INTENT(in),
OPTIONAL
:: transpose_tr if the triangular matrix should be transposed (defaults to false)
LOGICAL,
INTENT(in),
OPTIONAL
:: invert_tr if the triangular matrix should be inverted (defaults to false)
CHARACTER,
INTENT(in),
OPTIONAL
:: uplo_tr if triangular_matrix is stored in the upper ('U') or lower ('L') triangle (defaults to 'U')
LOGICAL,
INTENT(in),
OPTIONAL
:: unit_diag_tr if the diagonal elements of triangular_matrix should be assumed to be 1 (defaults to false)
INTEGER,
INTENT(in),
OPTIONAL
:: n_rows the number of rows of the result (defaults to size(matrix_b,1))
INTEGER,
INTENT(in),
OPTIONAL
:: n_cols the number of columns of the result (defaults to size(matrix_b,2))
REAL(dp),
INTENT(in),
OPTIONAL
:: alpha ...

SUBROUTINEcp_fm_upper_to_full(matrix, work)

given an upper triangular matrix computes the corresponding full matrix

Arguments:
TYPE(cp_fm_type),
POINTER
:: matrix the upper triangular matrix as input, the full matrix as output
TYPE(cp_fm_type),
POINTER
:: work a matrix of the same size as matrix

SUBROUTINEcp_fm_contracted_trace_a2b2t2(matrix_a, matrix_b, trace)

Compute trace(i,j) = \sum_k Tr (matrix_a(k,i)^T matrix_b(k,j)).

Arguments:
INTENT(in)
:: matrix_a(:,:) list of A matrices
INTENT(in)
:: matrix_b(:,:) list of B matrices
REAL(dp),
INTENT(out)
:: trace(:,:) computed traces

SUBROUTINEcp_fm_trace_a0b0t0(matrix_a, matrix_b, trace)

returns the trace of matrix_a^T matrix_b, i.e sum_{i,j}(matrix_a(i,j)*matrix_b(i,j))

Arguments:
TYPE(cp_fm_type),
POINTER
:: matrix_a a matrix
TYPE(cp_fm_type),
POINTER
:: matrix_b another matrix
REAL(dp),
INTENT(out)
:: trace ...

SUBROUTINEcp_fm_trace_a1b1t1(matrix_a, matrix_b, trace)

Compute trace(k) = Tr (matrix_a(k)^T matrix_b(k)) for each pair of matrices A_k and B_k.

Arguments:
INTENT(in)
:: matrix_a(:) list of A matrices
INTENT(in)
:: matrix_b(:) list of B matrices
REAL(dp),
INTENT(out)
:: trace(:) computed traces