Documentation for module cp_cfm_basic_linalg

Basic linear algebra operations for complex full matrices.

source: cp_cfm_basic_linalg.F
Loading...

Generic procedures:

cp_cfm_scale

public Subroutines/Functions:

Used to replace a symmetric positive definite matrix M with its Cholesky decomposition U: M = U^T * U, with U upper triangular.
Scales columns of the full matrix by corresponding factors.
Performs one of the matrix-matrix operations: matrix_c = alpha * op1( matrix_a ) * op2( matrix_b ) + beta*matrix_c.
Computes LU decomposition of a given matrix.
Inverts a matrix using LU decomposition. The input matrix will be overwritten.
FUNCTION
REAL(dp)
cp_cfm_norm (matrix, mode)
Norm of matrix using (p)zlange.
Scale and add two BLACS matrices (a = alpha*a + beta*b).
Scale and add two BLACS matrices (a = alpha*a + beta*b). where b is a real matrix (adapted from cp_cfm_scale_and_add).
Computes the element-wise (Schur) product of two matrices: C = A \circ B .
Solve the system of linear equations A*b=A_general using LU decomposition. Pay attention that both matrices are overwritten on exit and that the result is stored into the matrix 'general_a'.
Returns the trace of matrix_a^T matrix_b, i.e sum_{i,j}(matrix_a(i,j)*matrix_b(i,j)) .
Transposes a BLACS distributed complex matrix.
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 transa="N" and invert_tr=.false.) triangular_matrix^T (if transa="T" and invert_tr=.false.) triangular_matrix^H (if transa="C" and invert_tr=.false.) triangular_matrix^(-1) (if transa="N" and invert_tr=.true.) triangular_matrix^(-T) (if transa="T" and invert_tr=.true.) triangular_matrix^(-H) (if transa="C" and invert_tr=.true.)

Generic procedure cp_cfm_scale

...

REAL(dp),
INTENT(in)
:: alpha ×
COMPLEX(dp),
INTENT(in)
:: alpha ×
POINTER
:: matrix_a × ×

SUBROUTINEcp_cfm_cholesky_decompose(matrix, n, info_out)

Used to replace a symmetric positive definite matrix M with its Cholesky decomposition U: M = U^T * U, with U upper triangular.

Arguments:
POINTER
:: matrix the matrix to replace with its Cholesky decomposition
INTEGER,
INTENT(in),
OPTIONAL
:: n the number of row (and columns) of the matrix & (defaults to the min(size(matrix)))
INTEGER,
INTENT(out),
OPTIONAL
:: info_out if present, outputs info from (p)zpotrf

SUBROUTINEcp_cfm_column_scale(matrix_a, scaling)

Scales columns of the full matrix by corresponding factors.

Arguments:
POINTER
:: matrix_a matrix to scale
COMPLEX(dp),
INTENT(in)
:: scaling(:) scale factors for every column. The actual number of scaled columns is limited by the number of scale factors given or by the actual number of columns whichever is smaller.

SUBROUTINEcp_cfm_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)

Performs one of the matrix-matrix operations: matrix_c = alpha * op1( matrix_a ) * op2( matrix_b ) + beta*matrix_c.

Arguments:
CHARACTER(1),
INTENT(in)
:: transa form of op1( matrix_a ): op1( matrix_a ) = matrix_a, when transa == 'N' , op1( matrix_a ) = matrix_a^T, when transa == 'T' , op1( matrix_a ) = matrix_a^H, when transa == 'C' ,
CHARACTER(1),
INTENT(in)
:: transb form of op2( matrix_b )
INTEGER,
INTENT(in)
:: m number of rows of the matrix op1( matrix_a )
INTEGER,
INTENT(in)
:: n number of columns of the matrix op2( matrix_b )
INTEGER,
INTENT(in)
:: k number of columns of the matrix op1( matrix_a ) as well as number of rows of the matrix op2( matrix_b )
COMPLEX(dp),
INTENT(in)
:: alpha scale factor
POINTER
:: matrix_a matrix A
POINTER
:: matrix_b matrix B
COMPLEX(dp),
INTENT(in)
:: beta scale factor
POINTER
:: matrix_c matrix C
INTEGER,
INTENT(in),
OPTIONAL
:: a_first_col the first column of the matrix_a to multiply
INTEGER,
INTENT(in),
OPTIONAL
:: a_first_row the first row of the matrix_a to multiply
INTEGER,
INTENT(in),
OPTIONAL
:: b_first_col the first column of the matrix_b to multiply
INTEGER,
INTENT(in),
OPTIONAL
:: b_first_row the first row of the matrix_b to multiply
INTEGER,
INTENT(in),
OPTIONAL
:: c_first_col the first column of the matrix_c
INTEGER,
INTENT(in),
OPTIONAL
:: c_first_row the first row of the matrix_c

SUBROUTINEcp_cfm_lu_decompose(matrix_a, determinant)

Computes LU decomposition of a given matrix.

Arguments:
POINTER
:: matrix_a full matrix
COMPLEX(dp),
INTENT(out)
:: determinant determinant

SUBROUTINEcp_cfm_lu_invert(matrix, info_out)

Inverts a matrix using LU decomposition. The input matrix will be overwritten.

Arguments:
POINTER
:: matrix input a general square non-singular matrix, outputs its inverse
INTEGER,
INTENT(out),
OPTIONAL
:: info_out optional, if present outputs the info from (p)zgetri

FUNCTIONcp_cfm_norm(matrix, mode)

Norm of matrix using (p)zlange.

Return Value ::
REAL(dp)
Arguments:
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_cfm_scale_and_add(alpha, matrix_a, beta, matrix_b)

Scale and add two BLACS matrices (a = alpha*a + beta*b).

Arguments:
COMPLEX(dp),
INTENT(in)
:: alpha ...
POINTER
:: matrix_a ...
COMPLEX(dp),
INTENT(in),
OPTIONAL
:: beta ...
OPTIONAL, POINTER
:: matrix_b ...

SUBROUTINEcp_cfm_scale_and_add_fm(alpha, matrix_a, beta, matrix_b)

Scale and add two BLACS matrices (a = alpha*a + beta*b). where b is a real matrix (adapted from cp_cfm_scale_and_add).

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

SUBROUTINEcp_cfm_schur_product(matrix_a, matrix_b, matrix_c)

Computes the element-wise (Schur) product of two matrices: C = A \circ B .

Arguments:
POINTER
:: matrix_a the first input matrix
POINTER
:: matrix_b the second input matrix
POINTER
:: matrix_c matrix to store the result

SUBROUTINEcp_cfm_solve(matrix_a, general_a, determinant)

Solve the system of linear equations A*b=A_general using LU decomposition. Pay attention that both matrices are overwritten on exit and that the result is stored into the matrix 'general_a'.

Arguments:
POINTER
:: matrix_a matrix A (overwritten on exit)
POINTER
:: general_a (input) matrix A_general, (output) matrix B
COMPLEX(dp),
OPTIONAL
:: determinant determinant

SUBROUTINEcp_cfm_trace(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:
POINTER
:: matrix_a a complex matrix
POINTER
:: matrix_b another complex matrix
COMPLEX(dp),
INTENT(out)
:: trace value of the trace operator

SUBROUTINEcp_cfm_transpose(matrix, trans, matrixt)

Transposes a BLACS distributed complex matrix.

Arguments:
POINTER
:: matrix input matrix
CHARACTER,
INTENT(in)
:: trans 'T' for transpose, 'C' for Hermitian conjugate
POINTER
:: matrixt output matrix

SUBROUTINEcp_cfm_triangular_invert(matrix_a, uplo, info_out)

Inverts a triangular matrix.

Arguments:
POINTER
:: matrix_a ...
CHARACTER,
INTENT(in),
OPTIONAL
:: uplo ...
INTEGER,
INTENT(out),
OPTIONAL
:: info_out ...

SUBROUTINEcp_cfm_triangular_multiply(triangular_matrix, matrix_b, side, transa_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 transa="N" and invert_tr=.false.) triangular_matrix^T (if transa="T" and invert_tr=.false.) triangular_matrix^H (if transa="C" and invert_tr=.false.) triangular_matrix^(-1) (if transa="N" and invert_tr=.true.) triangular_matrix^(-T) (if transa="T" and invert_tr=.true.) triangular_matrix^(-H) (if transa="C" and invert_tr=.true.)

Arguments:
POINTER
:: triangular_matrix the triangular matrix that multiplies the other
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')
CHARACTER,
INTENT(in),
OPTIONAL
:: transa_tr ...
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))
COMPLEX(dp),
INTENT(in),
OPTIONAL
:: alpha ...

SUBROUTINEcp_cfm_dscale(alpha, matrix_a)

Scales a complex matrix by a real number. matrix_a = alpha * matrix_b

Arguments:
REAL(dp),
INTENT(in)
:: alpha scale factor
POINTER
:: matrix_a complex matrix to scale

SUBROUTINEcp_cfm_zscale(alpha, matrix_a)

Scales a complex matrix by a complex number. matrix_a = alpha * matrix_b

Arguments:
COMPLEX(dp),
INTENT(in)
:: alpha scale factor
POINTER
:: matrix_a complex matrix to scale