Documentation for module almo_scf_methods

Subroutines for ALMO SCF

source: almo_scf_methods.F
Loading...

public Subroutines/Functions:

computes ALMOs by diagonalizing the projected blocked KS matrix uses the diagonalization code for blocks computes both the occupied and virtual orbitals
computes the projected KS from the total KS matrix also computes the DIIS error vector as a by-product
builds projected KS matrices for the overlapping domains also computes the DIIS error vector as a by-product
ALMOs by diagonalizing the KS domain submatrices computes both the occupied and virtual orbitals
computes occupied ALMOs from the superimposed atomic density blocks
computes the idempotent density matrix from ALMOs
orthogonalize ALMOs within a domain (obsolete, use orthogonalize_mos)
computes the idempotent density matrix from MOs MOs can be either orthogonal or non-orthogonal
Parallel code for domain specific operations (my_action) 0. out = op1 * in 1. out = in - op2 * op1 * in
applies projector to the orbitals |psi_out> = P |psi_in> OR |psi_out> = (1-P) |psi_in>, where P = |psi_proj> ()^{-1}
Constructs preconditioners for each domain -1. projected preconditioner 0. simple preconditioner
Constructs subblocks of the covariant-covariant DM
Constructs S_inv block for each domain
Constructs S^(+1/2) and S^(-1/2) submatrices for each domain
Tests construction and release of domain submatrices
Load balancing of the submatrix computations
computes a unitary matrix from an arbitrary "generator" matrix U = ( 1 - X + tr(X) ) ( 1 + X - tr(X) )^(-1)
Computes the overlap matrix of MO orbitals
computes the step matrix from the gradient and Hessian using the Newton-Raphson method
orthogonalize MOs
inverts block-diagonal blocks of a cp_dbcsr_matrix
create the initial guess for XALMOs

SUBROUTINEalmo_scf_ks_blk_to_tv_blk(almo_scf_env)

computes ALMOs by diagonalizing the projected blocked KS matrix uses the diagonalization code for blocks computes both the occupied and virtual orbitals

Arguments:
INTENT(inout)
:: almo_scf_env ...

SUBROUTINEalmo_scf_ks_to_ks_blk(almo_scf_env)

computes the projected KS from the total KS matrix also computes the DIIS error vector as a by-product

Arguments:
INTENT(inout)
:: almo_scf_env ...

SUBROUTINEalmo_scf_ks_to_ks_xx(almo_scf_env)

builds projected KS matrices for the overlapping domains also computes the DIIS error vector as a by-product

Arguments:
INTENT(inout)
:: almo_scf_env ...

SUBROUTINEalmo_scf_ks_xx_to_tv_xx(almo_scf_env)

ALMOs by diagonalizing the KS domain submatrices computes both the occupied and virtual orbitals

Arguments:
INTENT(inout)
:: almo_scf_env ...

SUBROUTINEalmo_scf_p_blk_to_t_blk(almo_scf_env, ionic)

computes occupied ALMOs from the superimposed atomic density blocks

Arguments:
INTENT(inout)
:: almo_scf_env ...
LOGICAL,
INTENT(in)
:: ionic ...

SUBROUTINEalmo_scf_t_blk_to_p(almo_scf_env, use_sigma_inv_guess)

computes the idempotent density matrix from ALMOs

Arguments:
INTENT(inout)
:: almo_scf_env ...
LOGICAL,
INTENT(in),
OPTIONAL
:: use_sigma_inv_guess ...

SUBROUTINEalmo_scf_t_blk_to_t_blk_orthonormal(almo_scf_env)

orthogonalize ALMOs within a domain (obsolete, use orthogonalize_mos)

Arguments:
INTENT(inout)
:: almo_scf_env ...

SUBROUTINEalmo_scf_t_to_p(t, p, eps_filter, orthog_orbs, s, sigma, sigma_inv, use_guess)

computes the idempotent density matrix from MOs MOs can be either orthogonal or non-orthogonal

Arguments:
TYPE(dbcsr_type),
INTENT(in)
:: t ...
TYPE(dbcsr_type),
INTENT(inout)
:: p ...
REAL(dp),
INTENT(in)
:: eps_filter ...
LOGICAL,
INTENT(in)
:: orthog_orbs ...
TYPE(dbcsr_type),
INTENT(in),
OPTIONAL
:: s ...
TYPE(dbcsr_type),
INTENT(inout),
OPTIONAL
:: sigma ...
TYPE(dbcsr_type),
INTENT(inout),
OPTIONAL
:: sigma_inv ...
LOGICAL,
INTENT(in),
OPTIONAL
:: use_guess ...

SUBROUTINEapply_domain_operators(matrix_in, matrix_out, operator1, operator2, dpattern, map, node_of_domain, my_action, filter_eps, matrix_trimmer, use_trimmer)

Parallel code for domain specific operations (my_action) 0. out = op1 * in 1. out = in - op2 * op1 * in

Arguments:
TYPE(dbcsr_type),
INTENT(in)
:: matrix_in ...
TYPE(dbcsr_type),
INTENT(inout)
:: matrix_out ...
INTENT(in)
:: operator1(:) ...
INTENT(in),
OPTIONAL
:: operator2(:) ...
TYPE(dbcsr_type),
INTENT(in)
:: dpattern ...
INTENT(in)
:: map ...
INTEGER,
INTENT(in)
:: node_of_domain(:) ...
INTEGER,
INTENT(in)
:: my_action ...
REAL(dp)
:: filter_eps ...
TYPE(dbcsr_type),
INTENT(in),
OPTIONAL
:: matrix_trimmer ...
LOGICAL,
INTENT(in),
OPTIONAL
:: use_trimmer ...

SUBROUTINEapply_projector(psi_in, psi_out, psi_projector, metric, project_out, psi_projector_orthogonal, proj_in_template, eps_filter, sig_inv_projector, sig_inv_template)

applies projector to the orbitals |psi_out> = P |psi_in> OR |psi_out> = (1-P) |psi_in>, where P = |psi_proj> ()^{-1}

Arguments:
TYPE(dbcsr_type),
INTENT(in)
:: psi_in ...
TYPE(dbcsr_type),
INTENT(inout)
:: psi_out ...
TYPE(dbcsr_type),
INTENT(in)
:: psi_projector ...
TYPE(dbcsr_type),
INTENT(in)
:: metric ...
LOGICAL,
INTENT(in)
:: project_out ...
LOGICAL,
INTENT(in)
:: psi_projector_orthogonal ...
TYPE(dbcsr_type),
INTENT(in)
:: proj_in_template ...
REAL(dp),
INTENT(in)
:: eps_filter ...
TYPE(dbcsr_type),
INTENT(in),
OPTIONAL
:: sig_inv_projector ...
TYPE(dbcsr_type),
INTENT(in),
OPTIONAL
:: sig_inv_template ...

SUBROUTINEconstruct_domain_preconditioner(matrix_main, subm_s_inv, subm_r_down, matrix_trimmer, dpattern, map, node_of_domain, preconditioner, use_trimmer, my_action)

Constructs preconditioners for each domain -1. projected preconditioner 0. simple preconditioner

Arguments:
TYPE(dbcsr_type),
INTENT(inout)
:: matrix_main ...
INTENT(in),
OPTIONAL
:: subm_s_inv(:) ...
INTENT(in),
OPTIONAL
:: subm_r_down(:) ...
TYPE(dbcsr_type),
INTENT(inout),
OPTIONAL
:: matrix_trimmer ...
TYPE(dbcsr_type),
INTENT(in)
:: dpattern ...
INTENT(in)
:: map ...
INTEGER,
INTENT(in)
:: node_of_domain(:) ...
INTENT(inout)
:: preconditioner(:) ...
LOGICAL,
INTENT(in),
OPTIONAL
:: use_trimmer ...
INTEGER,
INTENT(in)
:: my_action ...

SUBROUTINEconstruct_domain_r_down(matrix_t, matrix_sigma_inv, matrix_s, subm_r_down, dpattern, map, node_of_domain, filter_eps)

Constructs subblocks of the covariant-covariant DM

Arguments:
TYPE(dbcsr_type),
INTENT(in)
:: matrix_t ...
TYPE(dbcsr_type),
INTENT(in)
:: matrix_sigma_inv ...
TYPE(dbcsr_type),
INTENT(in)
:: matrix_s ...
INTENT(inout)
:: subm_r_down(:) ...
TYPE(dbcsr_type),
INTENT(in)
:: dpattern ...
INTENT(in)
:: map ...
INTEGER,
INTENT(in)
:: node_of_domain(:) ...
REAL(dp)
:: filter_eps ...

SUBROUTINEconstruct_domain_s_inv(matrix_s, subm_s_inv, dpattern, map, node_of_domain)

Constructs S_inv block for each domain

Arguments:
TYPE(dbcsr_type),
INTENT(in)
:: matrix_s ...
INTENT(inout)
:: subm_s_inv(:) ...
TYPE(dbcsr_type),
INTENT(in)
:: dpattern ...
INTENT(in)
:: map ...
INTEGER,
INTENT(in)
:: node_of_domain(:) ...

SUBROUTINEconstruct_domain_s_sqrt(matrix_s, subm_s_sqrt, subm_s_sqrt_inv, dpattern, map, node_of_domain)

Constructs S^(+1/2) and S^(-1/2) submatrices for each domain

Arguments:
TYPE(dbcsr_type),
INTENT(in)
:: matrix_s ...
INTENT(inout)
:: subm_s_sqrt(:) ...
INTENT(inout)
:: subm_s_sqrt_inv(:) ...
TYPE(dbcsr_type),
INTENT(in)
:: dpattern ...
INTENT(in)
:: map ...
INTEGER,
INTENT(in)
:: node_of_domain(:) ...

SUBROUTINEconstruct_test(matrix_no, dpattern, map, node_of_domain)

Tests construction and release of domain submatrices

Arguments:
TYPE(dbcsr_type),
INTENT(in)
:: matrix_no ...
TYPE(dbcsr_type),
INTENT(in)
:: dpattern ...
INTENT(in)
:: map ...
INTEGER,
INTENT(in)
:: node_of_domain(:) ...

SUBROUTINEdistribute_domains(almo_scf_env)

Load balancing of the submatrix computations

Arguments:
INTENT(inout)
:: almo_scf_env ...

SUBROUTINEgenerator_to_unitary(x, u, eps_filter)

computes a unitary matrix from an arbitrary "generator" matrix U = ( 1 - X + tr(X) ) ( 1 + X - tr(X) )^(-1)

Arguments:
TYPE(dbcsr_type),
INTENT(in)
:: x ...
TYPE(dbcsr_type),
INTENT(inout)
:: u ...
REAL(dp),
INTENT(in)
:: eps_filter ...

SUBROUTINEget_overlap(bra, ket, overlap, metric, retain_overlap_sparsity, eps_filter)

Computes the overlap matrix of MO orbitals

Arguments:
TYPE(dbcsr_type),
INTENT(in)
:: bra ...
TYPE(dbcsr_type),
INTENT(in)
:: ket ...
TYPE(dbcsr_type),
INTENT(inout)
:: overlap ...
TYPE(dbcsr_type),
INTENT(in)
:: metric ...
LOGICAL,
INTENT(in),
OPTIONAL
:: retain_overlap_sparsity ...
REAL(dp)
:: eps_filter ...

SUBROUTINEnewton_grad_to_step(matrix_grad, matrix_step, matrix_s, matrix_ks, matrix_t, matrix_sigma_inv, quench_t, spin_factor, eps_filter)

computes the step matrix from the gradient and Hessian using the Newton-Raphson method

Arguments:
TYPE(dbcsr_type),
INTENT(inout)
:: matrix_grad ...
TYPE(dbcsr_type),
INTENT(inout)
:: matrix_step ...
TYPE(dbcsr_type),
INTENT(inout)
:: matrix_s ...
TYPE(dbcsr_type),
INTENT(in)
:: matrix_ks ...
TYPE(dbcsr_type),
INTENT(in)
:: matrix_t ...
TYPE(dbcsr_type),
INTENT(inout)
:: matrix_sigma_inv ...
TYPE(dbcsr_type),
INTENT(inout)
:: quench_t ...
REAL(dp),
INTENT(in)
:: spin_factor ...
REAL(dp),
INTENT(in)
:: eps_filter ...

SUBROUTINEorthogonalize_mos(ket, overlap, metric, retain_locality, only_normalize, eps_filter, order_lanczos, eps_lanczos, max_iter_lanczos)

orthogonalize MOs

Arguments:
TYPE(dbcsr_type),
INTENT(inout)
:: ket ...
TYPE(dbcsr_type),
INTENT(inout)
:: overlap ...
TYPE(dbcsr_type),
INTENT(in)
:: metric ...
LOGICAL,
INTENT(in)
:: retain_locality ...
LOGICAL,
INTENT(in)
:: only_normalize ...
REAL(dp)
:: eps_filter ...
INTEGER,
INTENT(in)
:: order_lanczos ...
REAL(dp),
INTENT(in)
:: eps_lanczos ...
INTEGER,
INTENT(in)
:: max_iter_lanczos ...

SUBROUTINEpseudo_invert_diagonal_blk(matrix_in, matrix_out, nocc)

inverts block-diagonal blocks of a cp_dbcsr_matrix

Arguments:
TYPE(dbcsr_type),
INTENT(in)
:: matrix_in ...
TYPE(dbcsr_type),
INTENT(inout)
:: matrix_out ...
INTEGER
:: nocc(:) ...

SUBROUTINExalmo_initial_guess(m_guess, m_t_in, m_t0, m_quench_t, m_overlap, m_sigma_tmpl, nspins, xalmo_history, assume_t0_q0x, optimize_theta, envelope_amplitude, eps_filter, order_lanczos, eps_lanczos, max_iter_lanczos)

create the initial guess for XALMOs

Arguments:
TYPE(dbcsr_type),
INTENT(inout)
:: m_guess(:) ...
TYPE(dbcsr_type),
INTENT(in)
:: m_t_in(:) ...
TYPE(dbcsr_type),
INTENT(in)
:: m_t0(:) ...
TYPE(dbcsr_type),
INTENT(in)
:: m_quench_t(:) ...
TYPE(dbcsr_type),
INTENT(in)
:: m_overlap ...
TYPE(dbcsr_type),
INTENT(in)
:: m_sigma_tmpl(:) ...
INTEGER,
INTENT(in)
:: nspins ...
INTENT(in)
:: xalmo_history ...
LOGICAL,
INTENT(in)
:: assume_t0_q0x ...
LOGICAL,
INTENT(in)
:: optimize_theta ...
REAL(dp),
INTENT(in)
:: envelope_amplitude ...
REAL(dp),
INTENT(in)
:: eps_filter ...
INTEGER,
INTENT(in)
:: order_lanczos ...
REAL(dp),
INTENT(in)
:: eps_lanczos ...
INTEGER,
INTENT(in)
:: max_iter_lanczos ...