Documentation for module cp_cfm_basic_linalg

basic linear algebra operations for complex full matrixes

source: cp_cfm_basic_linalg.F
Loading...

public Subroutines/Functions:

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_add)
used to replace a symmetric positive def. matrix M with its cholesky decomposition U: M = U^T * U, with U upper triangular
scales column i of matrix a with scaling(i)
BLACS interface to the BLAS routine zgemm.
Computes the LU decomposition of a given matrix the actual purpose right now is to compute the determinant of a given matrix which is most efficiently done this way, but, indeed, destroys the matrix SERIOUS WARNING (KNOWN BUG) : the sign of the determinant depends on ipivot one should be able to find out if ipivot is an even or an odd permutation...
scales a matrix matrix_a = alpha * matrix_b
...
computs the the solution to A*b=A_general using lu decomposition pay attention, both matrices are overwritten, a_general contais the result
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.)

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

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

SUBROUTINEcp_cfm_cholesky_decompose(matrix, n, info_out)

used to replace a symmetric positive def. 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(matrixa, scaling)

scales column i of matrix a with scaling(i)

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

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)

BLACS interface to the BLAS routine zgemm.

Arguments:
CHARACTER(1),
INTENT(in)
:: transa ...
CHARACTER(1),
INTENT(in)
:: transb ...
INTEGER,
INTENT(in)
:: m ...
INTEGER,
INTENT(in)
:: n ...
INTEGER,
INTENT(in)
:: k ...
COMPLEX(dp),
INTENT(in)
:: alpha ...
POINTER
:: matrix_a ...
POINTER
:: matrix_b ...
COMPLEX(dp),
INTENT(in)
:: beta ...
POINTER
:: matrix_c ...
INTEGER,
INTENT(in),
OPTIONAL
:: a_first_col ...
INTEGER,
INTENT(in),
OPTIONAL
:: a_first_row ...
INTEGER,
INTENT(in),
OPTIONAL
:: b_first_col ...
INTEGER,
INTENT(in),
OPTIONAL
:: b_first_row ...
INTEGER,
INTENT(in),
OPTIONAL
:: c_first_col ...
INTEGER,
INTENT(in),
OPTIONAL
:: c_first_row ...

SUBROUTINEcp_cfm_lu_decompose(matrix_a, almost_determinant)

Computes the LU decomposition of a given matrix the actual purpose right now is to compute the determinant of a given matrix which is most efficiently done this way, but, indeed, destroys the matrix SERIOUS WARNING (KNOWN BUG) : the sign of the determinant depends on ipivot one should be able to find out if ipivot is an even or an odd permutation...

Arguments:
POINTER
:: matrix_a ...
COMPLEX(dp),
INTENT(out)
:: almost_determinant ...

SUBROUTINEcp_cfm_scale(alpha, matrix_a)

scales a matrix matrix_a = alpha * matrix_b

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

SUBROUTINEcp_cfm_schur_product(matrix_a, matrix_b, matrix_c)

...

Arguments:
POINTER
:: matrix_a ...
POINTER
:: matrix_b ...
POINTER
:: matrix_c ...

SUBROUTINEcp_cfm_solve(matrix_a, general_a, determinant)

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

Arguments:
POINTER
:: matrix_a ...
POINTER
:: general_a ...
COMPLEX(dp),
OPTIONAL
:: determinant ...

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 ...