determinants

Contains everything for the computation of the Hamiltonian matrix elements in the basis of orthogonal Slater determinants built on a restricted spin-orbitals basis.

The main providers for this module are:

  • determinants n_states: number of states to be computed

  • psi_det: list of determinants in the wave function used in many routines/providers of the Quantum Package.

  • psi_coef: list of coefficients, for all determinants n_states states, and all determinants.

The main routines for this module are:

  • i_H_j(): computes the Hamiltonian matrix element between two arbitrary Slater determinants.

  • i_H_j_s2(): computes the Hamiltonian and (\(\widehat{S^2}\)) matrix element between two arbitrary Slater determinants.

  • i_H_j_verbose(): returns the decomposition in terms of one- and two-body components of the Hamiltonian matrix elements between two arbitrary Slater determinants. Also return the fermionic phase factor.

  • i_H_psi(): computes the Hamiltonian matrix element between an arbitrary Slater determinant and a wave function composed of a sum of arbitrary Slater determinants.

For an example of how to use these routines and providers, take a look at example.irp.f.

EZFIO parameters

n_det_max

Maximum number of determinants in the wave function

Default: 1000000

n_det_print_wf

Maximum number of determinants to be printed with the program print_wf

Default: 10000

n_states

Number of states to consider

Default: 1

read_wf

If true, read the wave function from the EZFIO file

Default: False

pruning

If p>0., remove p*Ndet determinants at every iteration

Default: 0.

s2_eig

Force the wave function to be an eigenfunction of \(\widehat{S^2}\)

Default: True

weight_one_e_dm

Weight used in the calculation of the one-electron density matrix. 0: 1./(c_0^2), 1: 1/N_states, 2: input state-average weight, 3: 1/(Norm_L3(Psi))

Default: 2

weight_selection

Weight used in the selection. 0: input state-average weight, 1: 1./(c_0^2), 2: rPT2 matching, 3: variance matching, 4: variance and rPT2 matching, 5: variance minimization and matching, 6: CI coefficients

Default: 2

threshold_generators

Thresholds on generators (fraction of the square of the norm)

Default: 0.999

n_int

Number of integers required to represent bitstrings (set in module bitmask module)

bit_kind

(set in module bitmask module)

mo_label

Label of the MOs on which the determinants are expressed

n_det

Number of determinants in the current wave function

n_det_qp_edit

Number of determinants to print in qp_edit

psi_coef

Coefficients of the wave function

psi_det

Determinants of the variational space

psi_coef_qp_edit

Coefficients of the wave function

psi_det_qp_edit

Determinants of the variational space

expected_s2

Expected value of \(\widehat{S^2}\)

target_energy

Energy that should be obtained when truncating the wave function (optional)

Default: 0.

state_average_weight

Weight of the states in state-average calculations.

selection_factor

f such that the number of determinants to add is f * N_det * sqrt(N_states)

Default: 1.

thresh_sym

Thresholds to check if a determinant is connected with HF

Default: 1.e-15

pseudo_sym

If true, discard any Slater determinants with an interaction smaller than thresh_sym with HF.

Default: False

Providers

abs_psi_coef_max

File : determinants/determinants.irp.f

double precision, allocatable   :: psi_coef_max (N_states)
double precision, allocatable   :: psi_coef_min (N_states)
double precision, allocatable   :: abs_psi_coef_max     (N_states)
double precision, allocatable   :: abs_psi_coef_min     (N_states)

Max and min values of the coefficients

Needs:

  • mpi_master

  • n_states

  • psi_coef

abs_psi_coef_min

File : determinants/determinants.irp.f

double precision, allocatable   :: psi_coef_max (N_states)
double precision, allocatable   :: psi_coef_min (N_states)
double precision, allocatable   :: abs_psi_coef_max     (N_states)
double precision, allocatable   :: abs_psi_coef_min     (N_states)

Max and min values of the coefficients

Needs:

  • mpi_master

  • n_states

  • psi_coef

barycentric_electronic_energy

File : determinants/energy.irp.f

double precision, allocatable   :: barycentric_electronic_energy        (N_states)

\(E_n = \sum_i {c_i^{(n)}}^2 H_{ii}\)

Needs:

  • diagonal_h_matrix_on_psi_det

  • n_det

  • n_states

  • psi_coef

c0_weight

File : determinants/density_matrix.irp.f

double precision, allocatable   :: c0_weight    (N_states)

Weight of the states in the selection : \(\frac{1}{c_0^2}\) .

Needs:

  • n_states

  • psi_coef

Needed by:

  • state_average_weight

det_alpha_norm

File : determinants/spindeterminants.irp.f

double precision, allocatable   :: det_alpha_norm       (N_det_alpha_unique)
double precision, allocatable   :: det_beta_norm        (N_det_beta_unique)

Norm of the \(\alpha\) and \(\beta\) spin determinants in the wave function:

\(||D_\alpha||_i = \sum_j C_{ij}^2\)

Needs:

  • n_det

  • n_states

  • psi_bilinear_matrix_values

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • state_average_weight

det_beta_norm

File : determinants/spindeterminants.irp.f

double precision, allocatable   :: det_alpha_norm       (N_det_alpha_unique)
double precision, allocatable   :: det_beta_norm        (N_det_beta_unique)

Norm of the \(\alpha\) and \(\beta\) spin determinants in the wave function:

\(||D_\alpha||_i = \sum_j C_{ij}^2\)

Needs:

  • n_det

  • n_states

  • psi_bilinear_matrix_values

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • state_average_weight

det_to_occ_pattern

File : determinants/occ_pattern.irp.f

integer, allocatable    :: det_to_occ_pattern   (N_det)

Returns the index of the occupation pattern for each determinant

Needs:

  • elec_alpha_num

  • n_det

  • n_int

  • psi_det

  • psi_occ_pattern

Needed by:

  • pruned

  • psi_occ_pattern_hii

  • weight_occ_pattern

  • weight_occ_pattern_average

diagonal_h_matrix_on_psi_det

File : determinants/energy.irp.f

double precision, allocatable   :: diagonal_h_matrix_on_psi_det (N_det)

Diagonal of the Hamiltonian ordered as psi_det

Needs:

  • elec_alpha_num

  • elec_beta_num

  • elec_num

  • n_det

  • n_int

  • psi_det

  • ref_bitmask

  • ref_bitmask_energy

Needed by:

  • barycentric_electronic_energy

double_exc_bitmask

File : determinants/determinants_bitmasks.irp.f

integer(bit_kind), allocatable  :: double_exc_bitmask   (N_int,4,N_double_exc_bitmasks)

double_exc_bitmask(:,1,i) is the bitmask for holes of excitation 1

double_exc_bitmask(:,2,i) is the bitmask for particles of excitation 1

double_exc_bitmask(:,3,i) is the bitmask for holes of excitation 2

double_exc_bitmask(:,4,i) is the bitmask for particles of excitation 2

for a given couple of hole/particle excitations i.

Needs:

  • hf_bitmask

  • n_double_exc_bitmasks

  • n_int

expected_s2

File : determinants/s2.irp.f

double precision        :: expected_s2

Expected value of \(\widehat{S^2}\) : S*(S+1)

Needs:

  • elec_alpha_num

  • elec_beta_num

Needed by:

  • ci_electronic_energy

fock_operator_closed_shell_ref_bitmask

File : determinants/single_excitations.irp.f

double precision, allocatable   :: fock_operator_closed_shell_ref_bitmask       (mo_num,mo_num)

Needs:

  • full_ijkl_bitmask

  • mo_integrals_map

  • mo_num

  • mo_one_e_integrals

  • mo_two_e_integrals_in_map

  • n_int

  • ref_closed_shell_bitmask

fock_wee_closed_shell

File : determinants/single_excitation_two_e.irp.f

double precision, allocatable   :: fock_wee_closed_shell        (mo_num,mo_num)

Needs:

  • full_ijkl_bitmask

  • mo_integrals_map

  • mo_num

  • mo_two_e_integrals_in_map

  • n_int

  • ref_closed_shell_bitmask

h_apply_buffer_allocated

File : determinants/h_apply.irp.f

logical :: h_apply_buffer_allocated
integer(omp_lock_kind), allocatable     :: h_apply_buffer_lock  (64,0:nproc-1)

Buffer of determinants/coefficients/perturbative energy for H_apply. Uninitialized. Filled by H_apply subroutines.

Needs:

  • n_det

  • n_int

  • n_states

  • nproc

h_apply_buffer_lock

File : determinants/h_apply.irp.f

logical :: h_apply_buffer_allocated
integer(omp_lock_kind), allocatable     :: h_apply_buffer_lock  (64,0:nproc-1)

Buffer of determinants/coefficients/perturbative energy for H_apply. Uninitialized. Filled by H_apply subroutines.

Needs:

  • n_det

  • n_int

  • n_states

  • nproc

h_matrix_all_dets

File : determinants/utils.irp.f

double precision, allocatable   :: h_matrix_all_dets    (N_det,N_det)

\(\hat H\) matrix on the basis of the Slater determinants defined by psi_det

Needs:

  • big_array_coulomb_integrals

  • big_array_coulomb_integrals

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_map

  • mo_two_e_integrals_in_map

  • n_det

  • n_int

  • psi_det

Needed by:

  • ci_electronic_energy

  • psi_energy

h_matrix_cas

File : determinants/psi_cas.irp.f

double precision, allocatable   :: h_matrix_cas (N_det_cas,N_det_cas)

Needs:

  • big_array_coulomb_integrals

  • big_array_coulomb_integrals

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_map

  • mo_two_e_integrals_in_map

  • n_int

  • psi_cas

Needed by:

  • psi_cas_energy

  • psi_coef_cas_diagonalized

idx_cas

File : determinants/psi_cas.irp.f

integer(bit_kind), allocatable  :: psi_cas      (N_int,2,psi_det_size)
double precision, allocatable   :: psi_cas_coef (psi_det_size,n_states)
integer, allocatable    :: idx_cas      (psi_det_size)
integer :: n_det_cas

CAS wave function, defined from the application of the CAS bitmask on the determinants. idx_cas gives the indice of the CAS determinant in psi_det.

Needs:

  • act_bitmask

  • hf_bitmask

  • mpi_master

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • h_matrix_cas

  • psi_cas_energy

  • psi_cas_sorted_bit

  • psi_coef_cas_diagonalized

  • psi_non_cas

  • psi_non_cas_sorted_bit

idx_non_cas

File : determinants/psi_cas.irp.f

integer(bit_kind), allocatable  :: psi_non_cas  (N_int,2,psi_det_size)
double precision, allocatable   :: psi_non_cas_coef     (psi_det_size,n_states)
integer, allocatable    :: idx_non_cas  (psi_det_size)
integer :: n_det_non_cas

Set of determinants which are not part of the CAS, defined from the application of the CAS bitmask on the determinants. idx_non_cas gives the indice of the determinant in psi_det.

Needs:

  • n_det

  • n_int

  • n_states

  • psi_cas

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • psi_non_cas_sorted_bit

max_degree_exc

File : determinants/determinants.irp.f

integer :: max_degree_exc

Maximum degree of excitation in the wave function with respect to the Hartree-Fock determinant.

Needs:

  • hf_bitmask

  • n_det

  • n_int

  • psi_det

n_det

File : determinants/determinants.irp.f

integer :: n_det

Number of determinants in the wave function

Needs:

  • ezfio_filename

  • mo_label

  • mpi_master

  • nproc

  • read_wf

Needed by:

  • act_2_rdm_aa_mo

  • act_2_rdm_ab_mo

  • act_2_rdm_bb_mo

  • act_2_rdm_spin_trace_mo

  • barycentric_electronic_energy

  • ci_electronic_energy

  • ci_energy

  • det_alpha_norm

  • det_to_occ_pattern

  • diag_algorithm

  • diagonal_h_matrix_on_psi_det

  • dressed_column_idx

  • dressing_column_h

  • h_apply_buffer_allocated

  • h_matrix_all_dets

  • max_degree_exc

  • n_det_qp_edit

  • one_e_dm_mo_alpha

  • pruned

  • psi_average_norm_contrib

  • psi_bilinear_matrix

  • psi_bilinear_matrix_columns_loc

  • psi_bilinear_matrix_order_reverse

  • psi_bilinear_matrix_order_transp_reverse

  • psi_bilinear_matrix_transp_rows_loc

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_values

  • psi_cas

  • psi_coef

  • psi_det

  • psi_det_alpha

  • psi_det_alpha_unique

  • psi_det_beta

  • psi_det_beta_unique

  • psi_det_hii

  • psi_det_sorted

  • psi_det_sorted_bit

  • psi_energy

  • psi_energy_two_e

  • psi_non_cas

  • psi_occ_pattern

  • psi_occ_pattern_hii

  • s2_matrix_all_dets

  • s2_values

  • state_av_act_2_rdm_aa_mo

  • state_av_act_2_rdm_ab_mo

  • state_av_act_2_rdm_bb_mo

  • state_av_act_2_rdm_spin_trace_mo

  • weight_occ_pattern

  • weight_occ_pattern_average

n_det_alpha_unique

File : determinants/spindeterminants.irp.f_template_144

integer(bit_kind), allocatable  :: psi_det_alpha_unique (N_int,psi_det_size)
integer :: n_det_alpha_unique

Unique \(\alpha\) determinants

Needs:

  • mpi_master

  • n_det

  • n_int

  • psi_det_alpha

  • psi_det_size

Needed by:

  • det_alpha_norm

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix

  • psi_bilinear_matrix_transp_rows_loc

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_values

  • singles_alpha_csc

  • singles_alpha_csc_idx

n_det_beta_unique

File : determinants/spindeterminants.irp.f_template_144

integer(bit_kind), allocatable  :: psi_det_beta_unique  (N_int,psi_det_size)
integer :: n_det_beta_unique

Unique \(\beta\) determinants

Needs:

  • mpi_master

  • n_det

  • n_int

  • psi_det_beta

  • psi_det_size

Needed by:

  • det_alpha_norm

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix

  • psi_bilinear_matrix_columns_loc

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_values

  • singles_beta_csc

  • singles_beta_csc_idx

n_det_cas

File : determinants/psi_cas.irp.f

integer(bit_kind), allocatable  :: psi_cas      (N_int,2,psi_det_size)
double precision, allocatable   :: psi_cas_coef (psi_det_size,n_states)
integer, allocatable    :: idx_cas      (psi_det_size)
integer :: n_det_cas

CAS wave function, defined from the application of the CAS bitmask on the determinants. idx_cas gives the indice of the CAS determinant in psi_det.

Needs:

  • act_bitmask

  • hf_bitmask

  • mpi_master

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • h_matrix_cas

  • psi_cas_energy

  • psi_cas_sorted_bit

  • psi_coef_cas_diagonalized

  • psi_non_cas

  • psi_non_cas_sorted_bit

n_det_non_cas

File : determinants/psi_cas.irp.f

integer(bit_kind), allocatable  :: psi_non_cas  (N_int,2,psi_det_size)
double precision, allocatable   :: psi_non_cas_coef     (psi_det_size,n_states)
integer, allocatable    :: idx_non_cas  (psi_det_size)
integer :: n_det_non_cas

Set of determinants which are not part of the CAS, defined from the application of the CAS bitmask on the determinants. idx_non_cas gives the indice of the determinant in psi_det.

Needs:

  • n_det

  • n_int

  • n_states

  • psi_cas

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • psi_non_cas_sorted_bit

n_det_qp_edit

File : determinants/determinants.irp.f

integer :: n_det_qp_edit

Number of determinants to print in qp_edit

Needs:

  • n_det

n_double_exc_bitmasks

File : determinants/determinants_bitmasks.irp.f

integer :: n_double_exc_bitmasks

Number of double excitation bitmasks

Needed by:

  • double_exc_bitmask

n_occ_pattern

File : determinants/occ_pattern.irp.f

integer(bit_kind), allocatable  :: psi_occ_pattern      (N_int,2,psi_det_size)
integer :: n_occ_pattern

Array of the occ_patterns present in the wave function.

psi_occ_pattern(:,1,j) = j-th occ_pattern of the wave function : represents all the single occupations

psi_occ_pattern(:,2,j) = j-th occ_pattern of the wave function : represents all the double occupations

The occ patterns are sorted by occ_pattern_search_key()

Needs:

  • elec_alpha_num

  • n_det

  • n_int

  • psi_det

  • psi_det_size

Needed by:

  • det_to_occ_pattern

  • pruned

  • psi_occ_pattern_hii

  • psi_occ_pattern_sorted

  • weight_occ_pattern

  • weight_occ_pattern_average

n_single_exc_bitmasks

File : determinants/determinants_bitmasks.irp.f

integer :: n_single_exc_bitmasks

Number of single excitation bitmasks

Needed by:

  • single_exc_bitmask

one_e_dm_ao_alpha

File : determinants/density_matrix.irp.f

double precision, allocatable   :: one_e_dm_ao_alpha    (ao_num,ao_num)
double precision, allocatable   :: one_e_dm_ao_beta     (ao_num,ao_num)

One body density matrix on the AO basis : \(\rho_{AO}(\alpha), \rho_{AO}(\beta)\) .

Needs:

  • ao_num

  • mo_coef

  • mo_num

  • one_e_dm_mo_alpha_average

one_e_dm_ao_beta

File : determinants/density_matrix.irp.f

double precision, allocatable   :: one_e_dm_ao_alpha    (ao_num,ao_num)
double precision, allocatable   :: one_e_dm_ao_beta     (ao_num,ao_num)

One body density matrix on the AO basis : \(\rho_{AO}(\alpha), \rho_{AO}(\beta)\) .

Needs:

  • ao_num

  • mo_coef

  • mo_num

  • one_e_dm_mo_alpha_average

one_e_dm_dagger_mo_spin_index

File : determinants/density_matrix.irp.f

double precision, allocatable   :: one_e_dm_dagger_mo_spin_index        (mo_num,mo_num,N_states,2)

Needs:

  • mo_num

  • n_states

  • one_e_dm_mo_alpha

one_e_dm_mo

File : determinants/density_matrix.irp.f

double precision, allocatable   :: one_e_dm_mo  (mo_num,mo_num)

One-body density matrix

Needs:

  • mo_num

  • one_e_dm_mo_alpha_average

one_e_dm_mo_alpha

File : determinants/density_matrix.irp.f

double precision, allocatable   :: one_e_dm_mo_alpha    (mo_num,mo_num,N_states)
double precision, allocatable   :: one_e_dm_mo_beta     (mo_num,mo_num,N_states)

\(\alpha\) and \(\beta\) one-body density matrix for each state

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_num

  • n_det

  • n_int

  • n_states

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_values

  • psi_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

Needed by:

  • full_occ_2_rdm_aa_mo

  • full_occ_2_rdm_ab_mo

  • full_occ_2_rdm_bb_mo

  • full_occ_2_rdm_spin_trace_mo

  • one_e_dm_dagger_mo_spin_index

  • one_e_dm_mo_alpha_average

  • one_e_dm_mo_alpha_for_dft

  • one_e_dm_mo_beta_for_dft

  • one_e_dm_mo_diff

  • one_e_dm_mo_spin_index

  • psi_energy_h_core

one_e_dm_mo_alpha_average

File : determinants/density_matrix.irp.f

double precision, allocatable   :: one_e_dm_mo_alpha_average    (mo_num,mo_num)
double precision, allocatable   :: one_e_dm_mo_beta_average     (mo_num,mo_num)

\(\alpha\) and \(\beta\) one-body density matrix for each state

Needs:

  • mo_num

  • n_states

  • one_e_dm_mo_alpha

  • state_average_weight

Needed by:

  • one_e_dm_ao_alpha

  • one_e_dm_mo

  • one_e_dm_mo_alpha_for_dft

  • one_e_dm_mo_beta_for_dft

  • one_e_spin_density_mo

  • state_av_full_occ_2_rdm_aa_mo

  • state_av_full_occ_2_rdm_ab_mo

  • state_av_full_occ_2_rdm_bb_mo

  • state_av_full_occ_2_rdm_spin_trace_mo

one_e_dm_mo_beta

File : determinants/density_matrix.irp.f

double precision, allocatable   :: one_e_dm_mo_alpha    (mo_num,mo_num,N_states)
double precision, allocatable   :: one_e_dm_mo_beta     (mo_num,mo_num,N_states)

\(\alpha\) and \(\beta\) one-body density matrix for each state

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_num

  • n_det

  • n_int

  • n_states

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_values

  • psi_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

Needed by:

  • full_occ_2_rdm_aa_mo

  • full_occ_2_rdm_ab_mo

  • full_occ_2_rdm_bb_mo

  • full_occ_2_rdm_spin_trace_mo

  • one_e_dm_dagger_mo_spin_index

  • one_e_dm_mo_alpha_average

  • one_e_dm_mo_alpha_for_dft

  • one_e_dm_mo_beta_for_dft

  • one_e_dm_mo_diff

  • one_e_dm_mo_spin_index

  • psi_energy_h_core

one_e_dm_mo_beta_average

File : determinants/density_matrix.irp.f

double precision, allocatable   :: one_e_dm_mo_alpha_average    (mo_num,mo_num)
double precision, allocatable   :: one_e_dm_mo_beta_average     (mo_num,mo_num)

\(\alpha\) and \(\beta\) one-body density matrix for each state

Needs:

  • mo_num

  • n_states

  • one_e_dm_mo_alpha

  • state_average_weight

Needed by:

  • one_e_dm_ao_alpha

  • one_e_dm_mo

  • one_e_dm_mo_alpha_for_dft

  • one_e_dm_mo_beta_for_dft

  • one_e_spin_density_mo

  • state_av_full_occ_2_rdm_aa_mo

  • state_av_full_occ_2_rdm_ab_mo

  • state_av_full_occ_2_rdm_bb_mo

  • state_av_full_occ_2_rdm_spin_trace_mo

one_e_dm_mo_diff

File : determinants/density_matrix.irp.f

double precision, allocatable   :: one_e_dm_mo_diff     (mo_num,mo_num,2:N_states)

Difference of the one-body density matrix with respect to the ground state

Needs:

  • mo_num

  • n_states

  • one_e_dm_mo_alpha

one_e_dm_mo_spin_index

File : determinants/density_matrix.irp.f

double precision, allocatable   :: one_e_dm_mo_spin_index       (mo_num,mo_num,N_states,2)

Needs:

  • mo_num

  • n_states

  • one_e_dm_mo_alpha

one_e_spin_density_ao

File : determinants/density_matrix.irp.f

double precision, allocatable   :: one_e_spin_density_ao        (ao_num,ao_num)

One body spin density matrix on the AO basis : \(\rho_{AO}(\alpha) - \rho_{AO}(\beta)\)

Needs:

  • ao_num

  • mo_coef

  • mo_num

  • one_e_spin_density_mo

one_e_spin_density_mo

File : determinants/density_matrix.irp.f

double precision, allocatable   :: one_e_spin_density_mo        (mo_num,mo_num)

\(\rho(\alpha) - \rho(\beta)\)

Needs:

  • mo_num

  • one_e_dm_mo_alpha_average

Needed by:

  • one_e_spin_density_ao

pruned

File : determinants/prune_wf.irp.f

logical, allocatable    :: pruned       (N_det)

True if determinant is removed by pruning

Needs:

  • det_to_occ_pattern

  • n_det

  • pruning

  • psi_average_norm_contrib

  • psi_det_sorted

  • psi_occ_pattern

  • psi_occ_pattern_sorted

  • s2_eig

psi_average_norm_contrib

File : determinants/determinants.irp.f

double precision, allocatable   :: psi_average_norm_contrib     (psi_det_size)

Contribution of determinants to the state-averaged density.

Needs:

  • n_det

  • n_states

  • psi_coef

  • psi_det_size

  • state_average_weight

Needed by:

  • pruned

  • psi_det_sorted

psi_average_norm_contrib_sorted

File : determinants/determinants.irp.f

integer(bit_kind), allocatable  :: psi_det_sorted       (N_int,2,psi_det_size)
double precision, allocatable   :: psi_coef_sorted      (psi_det_size,N_states)
double precision, allocatable   :: psi_average_norm_contrib_sorted      (psi_det_size)
integer, allocatable    :: psi_det_sorted_order (psi_det_size)

Wave function sorted by determinants contribution to the norm (state-averaged)

psi_det_sorted_order(i) -> k : index in psi_det

Needs:

  • n_det

  • n_int

  • n_states

  • psi_average_norm_contrib

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • pruned

psi_bilinear_matrix

File : determinants/spindeterminants.irp.f

double precision, allocatable   :: psi_bilinear_matrix  (N_det_alpha_unique,N_det_beta_unique,N_states)

Coefficient matrix if the wave function is expressed in a bilinear form :

\(D_\alpha^\dagger.C.D_\beta\)

Needs:

  • n_det

  • n_states

  • psi_bilinear_matrix_values

  • psi_det_alpha_unique

  • psi_det_beta_unique

psi_bilinear_matrix_columns

File : determinants/spindeterminants.irp.f

double precision, allocatable   :: psi_bilinear_matrix_values   (N_det,N_states)
integer, allocatable    :: psi_bilinear_matrix_rows     (N_det)
integer, allocatable    :: psi_bilinear_matrix_columns  (N_det)
integer, allocatable    :: psi_bilinear_matrix_order    (N_det)
Sparse coefficient matrix if the wave function is expressed in a bilinear form :

\(D_\alpha^\dagger.C.D_\beta\)

Rows are \(\alpha\) determinants and columns are \(\beta\) .

Order refers to psi_det

Needs:

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det_sorted_bit

  • psi_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

Needed by:

  • act_2_rdm_aa_mo

  • act_2_rdm_ab_mo

  • act_2_rdm_bb_mo

  • act_2_rdm_spin_trace_mo

  • det_alpha_norm

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix

  • psi_bilinear_matrix_columns_loc

  • psi_bilinear_matrix_order_reverse

  • psi_bilinear_matrix_transp_values

  • state_av_act_2_rdm_aa_mo

  • state_av_act_2_rdm_ab_mo

  • state_av_act_2_rdm_bb_mo

  • state_av_act_2_rdm_spin_trace_mo

psi_bilinear_matrix_columns_loc

File : determinants/spindeterminants.irp.f

integer, allocatable    :: psi_bilinear_matrix_columns_loc      (N_det_beta_unique+1)

Sparse coefficient matrix if the wave function is expressed in a bilinear form :

\(D_\alpha^\dagger.C.D_\beta\)

Rows are \(\alpha\) determinants and columns are \(\beta\) .

Order refers to psi_det

Needs:

  • n_det

  • psi_bilinear_matrix_values

  • psi_det_beta_unique

psi_bilinear_matrix_order

File : determinants/spindeterminants.irp.f

double precision, allocatable   :: psi_bilinear_matrix_values   (N_det,N_states)
integer, allocatable    :: psi_bilinear_matrix_rows     (N_det)
integer, allocatable    :: psi_bilinear_matrix_columns  (N_det)
integer, allocatable    :: psi_bilinear_matrix_order    (N_det)
Sparse coefficient matrix if the wave function is expressed in a bilinear form :

\(D_\alpha^\dagger.C.D_\beta\)

Rows are \(\alpha\) determinants and columns are \(\beta\) .

Order refers to psi_det

Needs:

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det_sorted_bit

  • psi_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

Needed by:

  • act_2_rdm_aa_mo

  • act_2_rdm_ab_mo

  • act_2_rdm_bb_mo

  • act_2_rdm_spin_trace_mo

  • det_alpha_norm

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix

  • psi_bilinear_matrix_columns_loc

  • psi_bilinear_matrix_order_reverse

  • psi_bilinear_matrix_transp_values

  • state_av_act_2_rdm_aa_mo

  • state_av_act_2_rdm_ab_mo

  • state_av_act_2_rdm_bb_mo

  • state_av_act_2_rdm_spin_trace_mo

psi_bilinear_matrix_order_reverse

File : determinants/spindeterminants.irp.f

integer, allocatable    :: psi_bilinear_matrix_order_reverse    (N_det)

Order which allows to go from psi_bilinear_matrix to psi_det

Needs:

  • n_det

  • psi_bilinear_matrix_values

Needed by:

  • act_2_rdm_aa_mo

  • act_2_rdm_ab_mo

  • act_2_rdm_bb_mo

  • act_2_rdm_spin_trace_mo

  • state_av_act_2_rdm_aa_mo

  • state_av_act_2_rdm_ab_mo

  • state_av_act_2_rdm_bb_mo

  • state_av_act_2_rdm_spin_trace_mo

psi_bilinear_matrix_order_transp_reverse

File : determinants/spindeterminants.irp.f

integer, allocatable    :: psi_bilinear_matrix_order_transp_reverse     (N_det)

Order which allows to go from psi_bilinear_matrix_order_transp to psi_bilinear_matrix

Needs:

  • n_det

  • psi_bilinear_matrix_transp_values

psi_bilinear_matrix_rows

File : determinants/spindeterminants.irp.f

double precision, allocatable   :: psi_bilinear_matrix_values   (N_det,N_states)
integer, allocatable    :: psi_bilinear_matrix_rows     (N_det)
integer, allocatable    :: psi_bilinear_matrix_columns  (N_det)
integer, allocatable    :: psi_bilinear_matrix_order    (N_det)
Sparse coefficient matrix if the wave function is expressed in a bilinear form :

\(D_\alpha^\dagger.C.D_\beta\)

Rows are \(\alpha\) determinants and columns are \(\beta\) .

Order refers to psi_det

Needs:

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det_sorted_bit

  • psi_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

Needed by:

  • act_2_rdm_aa_mo

  • act_2_rdm_ab_mo

  • act_2_rdm_bb_mo

  • act_2_rdm_spin_trace_mo

  • det_alpha_norm

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix

  • psi_bilinear_matrix_columns_loc

  • psi_bilinear_matrix_order_reverse

  • psi_bilinear_matrix_transp_values

  • state_av_act_2_rdm_aa_mo

  • state_av_act_2_rdm_ab_mo

  • state_av_act_2_rdm_bb_mo

  • state_av_act_2_rdm_spin_trace_mo

psi_bilinear_matrix_transp_columns

File : determinants/spindeterminants.irp.f

double precision, allocatable   :: psi_bilinear_matrix_transp_values    (N_det,N_states)
integer, allocatable    :: psi_bilinear_matrix_transp_rows      (N_det)
integer, allocatable    :: psi_bilinear_matrix_transp_columns   (N_det)
integer, allocatable    :: psi_bilinear_matrix_transp_order     (N_det)

Transpose of psi_bilinear_matrix

\(D_\beta^\dagger.C^\dagger.D_\alpha\)

Rows are \(\alpha\) determinants and columns are \(\beta\) , but the matrix is stored in row major format.

Needs:

  • n_det

  • n_states

  • psi_bilinear_matrix_values

  • psi_det_sorted_bit

  • psi_det_alpha_unique

  • psi_det_beta_unique

Needed by:

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix_order_transp_reverse

  • psi_bilinear_matrix_transp_rows_loc

psi_bilinear_matrix_transp_order

File : determinants/spindeterminants.irp.f

double precision, allocatable   :: psi_bilinear_matrix_transp_values    (N_det,N_states)
integer, allocatable    :: psi_bilinear_matrix_transp_rows      (N_det)
integer, allocatable    :: psi_bilinear_matrix_transp_columns   (N_det)
integer, allocatable    :: psi_bilinear_matrix_transp_order     (N_det)

Transpose of psi_bilinear_matrix

\(D_\beta^\dagger.C^\dagger.D_\alpha\)

Rows are \(\alpha\) determinants and columns are \(\beta\) , but the matrix is stored in row major format.

Needs:

  • n_det

  • n_states

  • psi_bilinear_matrix_values

  • psi_det_sorted_bit

  • psi_det_alpha_unique

  • psi_det_beta_unique

Needed by:

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix_order_transp_reverse

  • psi_bilinear_matrix_transp_rows_loc

psi_bilinear_matrix_transp_rows

File : determinants/spindeterminants.irp.f

double precision, allocatable   :: psi_bilinear_matrix_transp_values    (N_det,N_states)
integer, allocatable    :: psi_bilinear_matrix_transp_rows      (N_det)
integer, allocatable    :: psi_bilinear_matrix_transp_columns   (N_det)
integer, allocatable    :: psi_bilinear_matrix_transp_order     (N_det)

Transpose of psi_bilinear_matrix

\(D_\beta^\dagger.C^\dagger.D_\alpha\)

Rows are \(\alpha\) determinants and columns are \(\beta\) , but the matrix is stored in row major format.

Needs:

  • n_det

  • n_states

  • psi_bilinear_matrix_values

  • psi_det_sorted_bit

  • psi_det_alpha_unique

  • psi_det_beta_unique

Needed by:

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix_order_transp_reverse

  • psi_bilinear_matrix_transp_rows_loc

psi_bilinear_matrix_transp_rows_loc

File : determinants/spindeterminants.irp.f

integer, allocatable    :: psi_bilinear_matrix_transp_rows_loc  (N_det_alpha_unique+1)

Location of the columns in the psi_bilinear_matrix

Needs:

  • n_det

  • psi_bilinear_matrix_transp_values

  • psi_det_alpha_unique

psi_bilinear_matrix_transp_values

File : determinants/spindeterminants.irp.f

double precision, allocatable   :: psi_bilinear_matrix_transp_values    (N_det,N_states)
integer, allocatable    :: psi_bilinear_matrix_transp_rows      (N_det)
integer, allocatable    :: psi_bilinear_matrix_transp_columns   (N_det)
integer, allocatable    :: psi_bilinear_matrix_transp_order     (N_det)

Transpose of psi_bilinear_matrix

\(D_\beta^\dagger.C^\dagger.D_\alpha\)

Rows are \(\alpha\) determinants and columns are \(\beta\) , but the matrix is stored in row major format.

Needs:

  • n_det

  • n_states

  • psi_bilinear_matrix_values

  • psi_det_sorted_bit

  • psi_det_alpha_unique

  • psi_det_beta_unique

Needed by:

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix_order_transp_reverse

  • psi_bilinear_matrix_transp_rows_loc

psi_bilinear_matrix_values

File : determinants/spindeterminants.irp.f

double precision, allocatable   :: psi_bilinear_matrix_values   (N_det,N_states)
integer, allocatable    :: psi_bilinear_matrix_rows     (N_det)
integer, allocatable    :: psi_bilinear_matrix_columns  (N_det)
integer, allocatable    :: psi_bilinear_matrix_order    (N_det)
Sparse coefficient matrix if the wave function is expressed in a bilinear form :

\(D_\alpha^\dagger.C.D_\beta\)

Rows are \(\alpha\) determinants and columns are \(\beta\) .

Order refers to psi_det

Needs:

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det_sorted_bit

  • psi_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

Needed by:

  • act_2_rdm_aa_mo

  • act_2_rdm_ab_mo

  • act_2_rdm_bb_mo

  • act_2_rdm_spin_trace_mo

  • det_alpha_norm

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix

  • psi_bilinear_matrix_columns_loc

  • psi_bilinear_matrix_order_reverse

  • psi_bilinear_matrix_transp_values

  • state_av_act_2_rdm_aa_mo

  • state_av_act_2_rdm_ab_mo

  • state_av_act_2_rdm_bb_mo

  • state_av_act_2_rdm_spin_trace_mo

psi_cas

File : determinants/psi_cas.irp.f

integer(bit_kind), allocatable  :: psi_cas      (N_int,2,psi_det_size)
double precision, allocatable   :: psi_cas_coef (psi_det_size,n_states)
integer, allocatable    :: idx_cas      (psi_det_size)
integer :: n_det_cas

CAS wave function, defined from the application of the CAS bitmask on the determinants. idx_cas gives the indice of the CAS determinant in psi_det.

Needs:

  • act_bitmask

  • hf_bitmask

  • mpi_master

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • h_matrix_cas

  • psi_cas_energy

  • psi_cas_sorted_bit

  • psi_coef_cas_diagonalized

  • psi_non_cas

  • psi_non_cas_sorted_bit

psi_cas_coef

File : determinants/psi_cas.irp.f

integer(bit_kind), allocatable  :: psi_cas      (N_int,2,psi_det_size)
double precision, allocatable   :: psi_cas_coef (psi_det_size,n_states)
integer, allocatable    :: idx_cas      (psi_det_size)
integer :: n_det_cas

CAS wave function, defined from the application of the CAS bitmask on the determinants. idx_cas gives the indice of the CAS determinant in psi_det.

Needs:

  • act_bitmask

  • hf_bitmask

  • mpi_master

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • h_matrix_cas

  • psi_cas_energy

  • psi_cas_sorted_bit

  • psi_coef_cas_diagonalized

  • psi_non_cas

  • psi_non_cas_sorted_bit

psi_cas_coef_sorted_bit

File : determinants/psi_cas.irp.f

integer(bit_kind), allocatable  :: psi_cas_sorted_bit   (N_int,2,psi_det_size)
double precision, allocatable   :: psi_cas_coef_sorted_bit      (psi_det_size,N_states)

CAS determinants sorted to accelerate the search of a random determinant in the wave function.

Needs:

  • n_int

  • n_states

  • psi_cas

  • psi_det_size

psi_cas_energy

File : determinants/psi_cas.irp.f

double precision, allocatable   :: psi_cas_energy       (N_states)

Variational energy of \(\Psi_{CAS}\) , where \(\Psi_{CAS} = \sum_{I \in CAS} \I \rangle \langle I | \Psi \rangle\) .

Needs:

  • h_matrix_cas

  • n_states

  • psi_cas

psi_cas_energy_diagonalized

File : determinants/psi_cas.irp.f

double precision, allocatable   :: psi_coef_cas_diagonalized    (N_det_cas,N_states)
double precision, allocatable   :: psi_cas_energy_diagonalized  (N_states)

Needs:

  • h_matrix_cas

  • n_states

  • psi_cas

psi_cas_sorted_bit

File : determinants/psi_cas.irp.f

integer(bit_kind), allocatable  :: psi_cas_sorted_bit   (N_int,2,psi_det_size)
double precision, allocatable   :: psi_cas_coef_sorted_bit      (psi_det_size,N_states)

CAS determinants sorted to accelerate the search of a random determinant in the wave function.

Needs:

  • n_int

  • n_states

  • psi_cas

  • psi_det_size

psi_coef

File : determinants/determinants.irp.f

double precision, allocatable   :: psi_coef     (psi_det_size,N_states)

The wave function coefficients. Initialized with Hartree-Fock if the EZFIO file is empty.

Needs:

  • ezfio_filename

  • mo_label

  • mpi_master

  • n_det

  • n_states

  • psi_det_size

  • read_wf

Needed by:

  • act_2_rdm_aa_mo

  • act_2_rdm_ab_mo

  • act_2_rdm_bb_mo

  • act_2_rdm_spin_trace_mo

  • barycentric_electronic_energy

  • c0_weight

  • ci_electronic_energy

  • dressed_column_idx

  • psi_average_norm_contrib

  • psi_bilinear_matrix_values

  • psi_cas

  • psi_coef_max

  • psi_det_sorted

  • psi_det_sorted_bit

  • psi_energy

  • psi_energy_two_e

  • psi_non_cas

  • s2_values

  • state_av_act_2_rdm_aa_mo

  • state_av_act_2_rdm_ab_mo

  • state_av_act_2_rdm_bb_mo

  • state_av_act_2_rdm_spin_trace_mo

  • weight_occ_pattern

  • weight_occ_pattern_average

psi_coef_cas_diagonalized

File : determinants/psi_cas.irp.f

double precision, allocatable   :: psi_coef_cas_diagonalized    (N_det_cas,N_states)
double precision, allocatable   :: psi_cas_energy_diagonalized  (N_states)

Needs:

  • h_matrix_cas

  • n_states

  • psi_cas

psi_coef_max

File : determinants/determinants.irp.f

double precision, allocatable   :: psi_coef_max (N_states)
double precision, allocatable   :: psi_coef_min (N_states)
double precision, allocatable   :: abs_psi_coef_max     (N_states)
double precision, allocatable   :: abs_psi_coef_min     (N_states)

Max and min values of the coefficients

Needs:

  • mpi_master

  • n_states

  • psi_coef

psi_coef_min

File : determinants/determinants.irp.f

double precision, allocatable   :: psi_coef_max (N_states)
double precision, allocatable   :: psi_coef_min (N_states)
double precision, allocatable   :: abs_psi_coef_max     (N_states)
double precision, allocatable   :: abs_psi_coef_min     (N_states)

Max and min values of the coefficients

Needs:

  • mpi_master

  • n_states

  • psi_coef

psi_coef_sorted

File : determinants/determinants.irp.f

integer(bit_kind), allocatable  :: psi_det_sorted       (N_int,2,psi_det_size)
double precision, allocatable   :: psi_coef_sorted      (psi_det_size,N_states)
double precision, allocatable   :: psi_average_norm_contrib_sorted      (psi_det_size)
integer, allocatable    :: psi_det_sorted_order (psi_det_size)

Wave function sorted by determinants contribution to the norm (state-averaged)

psi_det_sorted_order(i) -> k : index in psi_det

Needs:

  • n_det

  • n_int

  • n_states

  • psi_average_norm_contrib

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • pruned

psi_coef_sorted_bit

File : determinants/determinants.irp.f

integer(bit_kind), allocatable  :: psi_det_sorted_bit   (N_int,2,psi_det_size)
double precision, allocatable   :: psi_coef_sorted_bit  (psi_det_size,N_states)

Determinants on which we apply \(\langle i|H|psi \rangle\) for perturbation. They are sorted by determinants interpreted as integers. Useful to accelerate the search of a random determinant in the wave function.

Needs:

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_values

psi_det

File : determinants/determinants.irp.f

integer(bit_kind), allocatable  :: psi_det      (N_int,2,psi_det_size)

The determinants of the wave function. Initialized with Hartree-Fock if the EZFIO file is empty.

Needs:

  • ezfio_filename

  • hf_bitmask

  • mo_coef

  • mo_label

  • mpi_master

  • n_det

  • n_int

  • psi_det_size

  • read_wf

Needed by:

  • ci_electronic_energy

  • det_to_occ_pattern

  • diagonal_h_matrix_on_psi_det

  • h_matrix_all_dets

  • max_degree_exc

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix_values

  • psi_cas

  • psi_det_alpha

  • psi_det_beta

  • psi_det_hii

  • psi_det_sorted

  • psi_det_sorted_bit

  • psi_energy

  • psi_energy_two_e

  • psi_non_cas

  • psi_occ_pattern

  • s2_matrix_all_dets

  • s2_values

psi_det_alpha

File : determinants/spindeterminants.irp.f

integer(bit_kind), allocatable  :: psi_det_alpha        (N_int,psi_det_size)

List of \(\alpha\) determinants of psi_det

Needs:

  • n_det

  • n_int

  • psi_det

  • psi_det_size

Needed by:

  • psi_det_alpha_unique

psi_det_alpha_unique

File : determinants/spindeterminants.irp.f_template_144

integer(bit_kind), allocatable  :: psi_det_alpha_unique (N_int,psi_det_size)
integer :: n_det_alpha_unique

Unique \(\alpha\) determinants

Needs:

  • mpi_master

  • n_det

  • n_int

  • psi_det_alpha

  • psi_det_size

Needed by:

  • det_alpha_norm

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix

  • psi_bilinear_matrix_transp_rows_loc

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_values

  • singles_alpha_csc

  • singles_alpha_csc_idx

psi_det_beta

File : determinants/spindeterminants.irp.f

integer(bit_kind), allocatable  :: psi_det_beta (N_int,psi_det_size)

List of \(\beta\) determinants of psi_det

Needs:

  • n_det

  • n_int

  • psi_det

  • psi_det_size

Needed by:

  • psi_det_beta_unique

psi_det_beta_unique

File : determinants/spindeterminants.irp.f_template_144

integer(bit_kind), allocatable  :: psi_det_beta_unique  (N_int,psi_det_size)
integer :: n_det_beta_unique

Unique \(\beta\) determinants

Needs:

  • mpi_master

  • n_det

  • n_int

  • psi_det_beta

  • psi_det_size

Needed by:

  • det_alpha_norm

  • one_e_dm_mo_alpha

  • psi_bilinear_matrix

  • psi_bilinear_matrix_columns_loc

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_values

  • singles_beta_csc

  • singles_beta_csc_idx

psi_det_hii

File : determinants/determinants.irp.f

double precision, allocatable   :: psi_det_hii  (N_det)

\(\langle i|h|i \rangle\) for all determinants.

Needs:

  • elec_alpha_num

  • elec_beta_num

  • elec_num

  • n_det

  • n_int

  • psi_det

  • ref_bitmask

  • ref_bitmask_energy

Needed by:

  • psi_occ_pattern_hii

psi_det_size

File : determinants/determinants.irp.f

integer :: psi_det_size

Size of the psi_det and psi_coef arrays

Needs:

  • ezfio_filename

  • mpi_master

Needed by:

  • psi_average_norm_contrib

  • psi_cas

  • psi_cas_sorted_bit

  • psi_coef

  • psi_det

  • psi_det_alpha

  • psi_det_alpha_unique

  • psi_det_beta

  • psi_det_beta_unique

  • psi_det_sorted

  • psi_det_sorted_bit

  • psi_energy

  • psi_energy_two_e

  • psi_non_cas

  • psi_non_cas_sorted_bit

  • psi_occ_pattern

  • s2_values

psi_det_sorted

File : determinants/determinants.irp.f

integer(bit_kind), allocatable  :: psi_det_sorted       (N_int,2,psi_det_size)
double precision, allocatable   :: psi_coef_sorted      (psi_det_size,N_states)
double precision, allocatable   :: psi_average_norm_contrib_sorted      (psi_det_size)
integer, allocatable    :: psi_det_sorted_order (psi_det_size)

Wave function sorted by determinants contribution to the norm (state-averaged)

psi_det_sorted_order(i) -> k : index in psi_det

Needs:

  • n_det

  • n_int

  • n_states

  • psi_average_norm_contrib

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • pruned

psi_det_sorted_bit

File : determinants/determinants.irp.f

integer(bit_kind), allocatable  :: psi_det_sorted_bit   (N_int,2,psi_det_size)
double precision, allocatable   :: psi_coef_sorted_bit  (psi_det_size,N_states)

Determinants on which we apply \(\langle i|H|psi \rangle\) for perturbation. They are sorted by determinants interpreted as integers. Useful to accelerate the search of a random determinant in the wave function.

Needs:

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_values

psi_det_sorted_order

File : determinants/determinants.irp.f

integer(bit_kind), allocatable  :: psi_det_sorted       (N_int,2,psi_det_size)
double precision, allocatable   :: psi_coef_sorted      (psi_det_size,N_states)
double precision, allocatable   :: psi_average_norm_contrib_sorted      (psi_det_size)
integer, allocatable    :: psi_det_sorted_order (psi_det_size)

Wave function sorted by determinants contribution to the norm (state-averaged)

psi_det_sorted_order(i) -> k : index in psi_det

Needs:

  • n_det

  • n_int

  • n_states

  • psi_average_norm_contrib

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • pruned

psi_energy_h_core

File : determinants/psi_energy_mono_elec.irp.f

double precision, allocatable   :: psi_energy_h_core    (N_states)

psi_energy_h_core = \(\langle \Psi | h_{core} |\Psi \rangle\)

computed using the one_e_dm_mo_alpha + one_e_dm_mo_beta and mo_one_e_integrals

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_num

  • mo_one_e_integrals

  • n_states

  • one_e_dm_mo_alpha

psi_non_cas

File : determinants/psi_cas.irp.f

integer(bit_kind), allocatable  :: psi_non_cas  (N_int,2,psi_det_size)
double precision, allocatable   :: psi_non_cas_coef     (psi_det_size,n_states)
integer, allocatable    :: idx_non_cas  (psi_det_size)
integer :: n_det_non_cas

Set of determinants which are not part of the CAS, defined from the application of the CAS bitmask on the determinants. idx_non_cas gives the indice of the determinant in psi_det.

Needs:

  • n_det

  • n_int

  • n_states

  • psi_cas

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • psi_non_cas_sorted_bit

psi_non_cas_coef

File : determinants/psi_cas.irp.f

integer(bit_kind), allocatable  :: psi_non_cas  (N_int,2,psi_det_size)
double precision, allocatable   :: psi_non_cas_coef     (psi_det_size,n_states)
integer, allocatable    :: idx_non_cas  (psi_det_size)
integer :: n_det_non_cas

Set of determinants which are not part of the CAS, defined from the application of the CAS bitmask on the determinants. idx_non_cas gives the indice of the determinant in psi_det.

Needs:

  • n_det

  • n_int

  • n_states

  • psi_cas

  • psi_coef

  • psi_det

  • psi_det_size

Needed by:

  • psi_non_cas_sorted_bit

psi_non_cas_coef_sorted_bit

File : determinants/psi_cas.irp.f

integer(bit_kind), allocatable  :: psi_non_cas_sorted_bit       (N_int,2,psi_det_size)
double precision, allocatable   :: psi_non_cas_coef_sorted_bit  (psi_det_size,N_states)

CAS determinants sorted to accelerate the search of a random determinant in the wave function.

Needs:

  • n_int

  • n_states

  • psi_cas

  • psi_det_size

  • psi_non_cas

psi_non_cas_sorted_bit

File : determinants/psi_cas.irp.f

integer(bit_kind), allocatable  :: psi_non_cas_sorted_bit       (N_int,2,psi_det_size)
double precision, allocatable   :: psi_non_cas_coef_sorted_bit  (psi_det_size,N_states)

CAS determinants sorted to accelerate the search of a random determinant in the wave function.

Needs:

  • n_int

  • n_states

  • psi_cas

  • psi_det_size

  • psi_non_cas

psi_occ_pattern

File : determinants/occ_pattern.irp.f

integer(bit_kind), allocatable  :: psi_occ_pattern      (N_int,2,psi_det_size)
integer :: n_occ_pattern

Array of the occ_patterns present in the wave function.

psi_occ_pattern(:,1,j) = j-th occ_pattern of the wave function : represents all the single occupations

psi_occ_pattern(:,2,j) = j-th occ_pattern of the wave function : represents all the double occupations

The occ patterns are sorted by occ_pattern_search_key()

Needs:

  • elec_alpha_num

  • n_det

  • n_int

  • psi_det

  • psi_det_size

Needed by:

  • det_to_occ_pattern

  • pruned

  • psi_occ_pattern_hii

  • psi_occ_pattern_sorted

  • weight_occ_pattern

  • weight_occ_pattern_average

psi_occ_pattern_hii

File : determinants/occ_pattern.irp.f

double precision, allocatable   :: psi_occ_pattern_hii  (N_occ_pattern)

\(\langle I|H|I \rangle\) where \(|I\rangle\) is an occupation pattern. This is the minimum \(H_{ii}\) , where the \(|i\rangle\) are the determinants of \(|I\rangle\) .

Needs:

  • det_to_occ_pattern

  • n_det

  • psi_det_hii

  • psi_occ_pattern

psi_occ_pattern_sorted

File : determinants/occ_pattern.irp.f

integer(bit_kind), allocatable  :: psi_occ_pattern_sorted       (N_int,2,N_occ_pattern)
double precision, allocatable   :: weight_occ_pattern_average_sorted    (N_occ_pattern)
integer, allocatable    :: psi_occ_pattern_sorted_order (N_occ_pattern)
integer, allocatable    :: psi_occ_pattern_sorted_order_reverse (N_occ_pattern)

Occupation patterns sorted by weight

Needs:

  • n_int

  • psi_occ_pattern

  • weight_occ_pattern_average

Needed by:

  • pruned

psi_occ_pattern_sorted_order

File : determinants/occ_pattern.irp.f

integer(bit_kind), allocatable  :: psi_occ_pattern_sorted       (N_int,2,N_occ_pattern)
double precision, allocatable   :: weight_occ_pattern_average_sorted    (N_occ_pattern)
integer, allocatable    :: psi_occ_pattern_sorted_order (N_occ_pattern)
integer, allocatable    :: psi_occ_pattern_sorted_order_reverse (N_occ_pattern)

Occupation patterns sorted by weight

Needs:

  • n_int

  • psi_occ_pattern

  • weight_occ_pattern_average

Needed by:

  • pruned

psi_occ_pattern_sorted_order_reverse

File : determinants/occ_pattern.irp.f

integer(bit_kind), allocatable  :: psi_occ_pattern_sorted       (N_int,2,N_occ_pattern)
double precision, allocatable   :: weight_occ_pattern_average_sorted    (N_occ_pattern)
integer, allocatable    :: psi_occ_pattern_sorted_order (N_occ_pattern)
integer, allocatable    :: psi_occ_pattern_sorted_order_reverse (N_occ_pattern)

Occupation patterns sorted by weight

Needs:

  • n_int

  • psi_occ_pattern

  • weight_occ_pattern_average

Needed by:

  • pruned

ref_bitmask_energy

File : determinants/ref_bitmask.irp.f

double precision        :: ref_bitmask_energy
double precision        :: ref_bitmask_one_e_energy
double precision        :: ref_bitmask_kinetic_energy
double precision        :: ref_bitmask_n_e_energy
double precision        :: ref_bitmask_two_e_energy
double precision        :: ref_bitmask_energy_ab
double precision        :: ref_bitmask_energy_bb
double precision        :: ref_bitmask_energy_aa

Energy of the reference bitmask used in Slater rules

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_n_e

  • mo_kinetic_integrals

  • mo_one_e_integrals

  • mo_two_e_integrals_jj

  • n_int

  • ref_bitmask

Needed by:

  • diagonal_h_matrix_on_psi_det

  • psi_det_hii

ref_bitmask_energy_aa

File : determinants/ref_bitmask.irp.f

double precision        :: ref_bitmask_energy
double precision        :: ref_bitmask_one_e_energy
double precision        :: ref_bitmask_kinetic_energy
double precision        :: ref_bitmask_n_e_energy
double precision        :: ref_bitmask_two_e_energy
double precision        :: ref_bitmask_energy_ab
double precision        :: ref_bitmask_energy_bb
double precision        :: ref_bitmask_energy_aa

Energy of the reference bitmask used in Slater rules

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_n_e

  • mo_kinetic_integrals

  • mo_one_e_integrals

  • mo_two_e_integrals_jj

  • n_int

  • ref_bitmask

Needed by:

  • diagonal_h_matrix_on_psi_det

  • psi_det_hii

ref_bitmask_energy_ab

File : determinants/ref_bitmask.irp.f

double precision        :: ref_bitmask_energy
double precision        :: ref_bitmask_one_e_energy
double precision        :: ref_bitmask_kinetic_energy
double precision        :: ref_bitmask_n_e_energy
double precision        :: ref_bitmask_two_e_energy
double precision        :: ref_bitmask_energy_ab
double precision        :: ref_bitmask_energy_bb
double precision        :: ref_bitmask_energy_aa

Energy of the reference bitmask used in Slater rules

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_n_e

  • mo_kinetic_integrals

  • mo_one_e_integrals

  • mo_two_e_integrals_jj

  • n_int

  • ref_bitmask

Needed by:

  • diagonal_h_matrix_on_psi_det

  • psi_det_hii

ref_bitmask_energy_bb

File : determinants/ref_bitmask.irp.f

double precision        :: ref_bitmask_energy
double precision        :: ref_bitmask_one_e_energy
double precision        :: ref_bitmask_kinetic_energy
double precision        :: ref_bitmask_n_e_energy
double precision        :: ref_bitmask_two_e_energy
double precision        :: ref_bitmask_energy_ab
double precision        :: ref_bitmask_energy_bb
double precision        :: ref_bitmask_energy_aa

Energy of the reference bitmask used in Slater rules

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_n_e

  • mo_kinetic_integrals

  • mo_one_e_integrals

  • mo_two_e_integrals_jj

  • n_int

  • ref_bitmask

Needed by:

  • diagonal_h_matrix_on_psi_det

  • psi_det_hii

ref_bitmask_kinetic_energy

File : determinants/ref_bitmask.irp.f

double precision        :: ref_bitmask_energy
double precision        :: ref_bitmask_one_e_energy
double precision        :: ref_bitmask_kinetic_energy
double precision        :: ref_bitmask_n_e_energy
double precision        :: ref_bitmask_two_e_energy
double precision        :: ref_bitmask_energy_ab
double precision        :: ref_bitmask_energy_bb
double precision        :: ref_bitmask_energy_aa

Energy of the reference bitmask used in Slater rules

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_n_e

  • mo_kinetic_integrals

  • mo_one_e_integrals

  • mo_two_e_integrals_jj

  • n_int

  • ref_bitmask

Needed by:

  • diagonal_h_matrix_on_psi_det

  • psi_det_hii

ref_bitmask_n_e_energy

File : determinants/ref_bitmask.irp.f

double precision        :: ref_bitmask_energy
double precision        :: ref_bitmask_one_e_energy
double precision        :: ref_bitmask_kinetic_energy
double precision        :: ref_bitmask_n_e_energy
double precision        :: ref_bitmask_two_e_energy
double precision        :: ref_bitmask_energy_ab
double precision        :: ref_bitmask_energy_bb
double precision        :: ref_bitmask_energy_aa

Energy of the reference bitmask used in Slater rules

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_n_e

  • mo_kinetic_integrals

  • mo_one_e_integrals

  • mo_two_e_integrals_jj

  • n_int

  • ref_bitmask

Needed by:

  • diagonal_h_matrix_on_psi_det

  • psi_det_hii

ref_bitmask_one_e_energy

File : determinants/ref_bitmask.irp.f

double precision        :: ref_bitmask_energy
double precision        :: ref_bitmask_one_e_energy
double precision        :: ref_bitmask_kinetic_energy
double precision        :: ref_bitmask_n_e_energy
double precision        :: ref_bitmask_two_e_energy
double precision        :: ref_bitmask_energy_ab
double precision        :: ref_bitmask_energy_bb
double precision        :: ref_bitmask_energy_aa

Energy of the reference bitmask used in Slater rules

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_n_e

  • mo_kinetic_integrals

  • mo_one_e_integrals

  • mo_two_e_integrals_jj

  • n_int

  • ref_bitmask

Needed by:

  • diagonal_h_matrix_on_psi_det

  • psi_det_hii

ref_bitmask_two_e_energy

File : determinants/ref_bitmask.irp.f

double precision        :: ref_bitmask_energy
double precision        :: ref_bitmask_one_e_energy
double precision        :: ref_bitmask_kinetic_energy
double precision        :: ref_bitmask_n_e_energy
double precision        :: ref_bitmask_two_e_energy
double precision        :: ref_bitmask_energy_ab
double precision        :: ref_bitmask_energy_bb
double precision        :: ref_bitmask_energy_aa

Energy of the reference bitmask used in Slater rules

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_n_e

  • mo_kinetic_integrals

  • mo_one_e_integrals

  • mo_two_e_integrals_jj

  • n_int

  • ref_bitmask

Needed by:

  • diagonal_h_matrix_on_psi_det

  • psi_det_hii

ref_closed_shell_bitmask

File : determinants/single_excitations.irp.f

integer(bit_kind), allocatable  :: ref_closed_shell_bitmask     (N_int,2)

Needs:

  • elec_alpha_num

  • elec_beta_num

  • n_int

  • ref_bitmask

Needed by:

  • fock_operator_closed_shell_ref_bitmask

  • fock_wee_closed_shell

s2_matrix_all_dets

File : determinants/utils.irp.f

double precision, allocatable   :: s2_matrix_all_dets   (N_det,N_det)

\(\widehat{S^2}\) matrix on the basis of the Slater determinants defined by psi_det

Needs:

  • n_det

  • n_int

  • psi_det

Needed by:

  • ci_electronic_energy

  • psi_energy

s2_values

File : determinants/s2.irp.f

double precision, allocatable   :: s2_values    (N_states)

array of the averaged values of the S^2 operator on the various states

Needs:

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det

  • psi_det_size

s_z

File : determinants/s2.irp.f

double precision        :: s_z
double precision        :: s_z2_sz

z component of the Spin

Needs:

  • elec_alpha_num

  • elec_beta_num

s_z2_sz

File : determinants/s2.irp.f

double precision        :: s_z
double precision        :: s_z2_sz

z component of the Spin

Needs:

  • elec_alpha_num

  • elec_beta_num

single_exc_bitmask

File : determinants/determinants_bitmasks.irp.f

integer(bit_kind), allocatable  :: single_exc_bitmask   (N_int,2,N_single_exc_bitmasks)

single_exc_bitmask(:,1,i) is the bitmask for holes

single_exc_bitmask(:,2,i) is the bitmask for particles

for a given couple of hole/particle excitations i.

Needs:

  • hf_bitmask

  • n_int

  • n_single_exc_bitmasks

singles_alpha_csc

File : determinants/spindeterminants.irp.f

integer, allocatable    :: singles_alpha_csc    (singles_alpha_csc_size)

Indices of all single excitations

Needs:

  • n_int

  • psi_det_alpha_unique

  • singles_alpha_csc_idx

singles_alpha_csc_idx

File : determinants/spindeterminants.irp.f

integer*8, allocatable  :: singles_alpha_csc_idx        (N_det_alpha_unique+1)
integer*8       :: singles_alpha_csc_size

singles_alpha_csc_size : Dimension of the singles_alpha_csc array

singles_alpha_csc_idx : Index where the single excitations of determinant i start

Needs:

  • elec_alpha_num

  • mo_num

  • n_int

  • psi_det_alpha_unique

Needed by:

  • singles_alpha_csc

singles_alpha_csc_size

File : determinants/spindeterminants.irp.f

integer*8, allocatable  :: singles_alpha_csc_idx        (N_det_alpha_unique+1)
integer*8       :: singles_alpha_csc_size

singles_alpha_csc_size : Dimension of the singles_alpha_csc array

singles_alpha_csc_idx : Index where the single excitations of determinant i start

Needs:

  • elec_alpha_num

  • mo_num

  • n_int

  • psi_det_alpha_unique

Needed by:

  • singles_alpha_csc

singles_beta_csc

File : determinants/spindeterminants.irp.f

integer, allocatable    :: singles_beta_csc     (singles_beta_csc_size)

Indices of all single excitations

Needs:

  • n_int

  • psi_det_beta_unique

  • singles_beta_csc_idx

singles_beta_csc_idx

File : determinants/spindeterminants.irp.f

integer*8, allocatable  :: singles_beta_csc_idx (N_det_beta_unique+1)
integer*8       :: singles_beta_csc_size

singles_beta_csc_size : Dimension of the singles_beta_csc array

singles_beta_csc_idx : Index where the single excitations of determinant i start

Needs:

  • elec_beta_num

  • mo_num

  • n_int

  • psi_det_beta_unique

Needed by:

  • singles_beta_csc

singles_beta_csc_size

File : determinants/spindeterminants.irp.f

integer*8, allocatable  :: singles_beta_csc_idx (N_det_beta_unique+1)
integer*8       :: singles_beta_csc_size

singles_beta_csc_size : Dimension of the singles_beta_csc array

singles_beta_csc_idx : Index where the single excitations of determinant i start

Needs:

  • elec_beta_num

  • mo_num

  • n_int

  • psi_det_beta_unique

Needed by:

  • singles_beta_csc

state_average_weight

File : determinants/density_matrix.irp.f

double precision, allocatable   :: state_average_weight (N_states)

Weights in the state-average calculation of the density matrix

Needs:

  • c0_weight

  • n_states

  • weight_one_e_dm

Needed by:

  • det_alpha_norm

  • one_e_dm_average_alpha_mo_for_dft

  • one_e_dm_average_beta_mo_for_dft

  • one_e_dm_mo_alpha_average

  • psi_average_norm_contrib

  • state_av_act_2_rdm_aa_mo

  • state_av_act_2_rdm_ab_mo

  • state_av_act_2_rdm_bb_mo

  • state_av_act_2_rdm_spin_trace_mo

  • weight_occ_pattern_average

weight_occ_pattern

File : determinants/occ_pattern.irp.f

double precision, allocatable   :: weight_occ_pattern   (N_occ_pattern,N_states)

Weight of the occupation patterns in the wave function

Needs:

  • det_to_occ_pattern

  • n_det

  • n_states

  • psi_coef

  • psi_occ_pattern

weight_occ_pattern_average

File : determinants/occ_pattern.irp.f

double precision, allocatable   :: weight_occ_pattern_average   (N_occ_pattern)

State-average weight of the occupation patterns in the wave function

Needs:

  • det_to_occ_pattern

  • n_det

  • n_states

  • psi_coef

  • psi_occ_pattern

  • state_average_weight

Needed by:

  • psi_occ_pattern_sorted

weight_occ_pattern_average_sorted

File : determinants/occ_pattern.irp.f

integer(bit_kind), allocatable  :: psi_occ_pattern_sorted       (N_int,2,N_occ_pattern)
double precision, allocatable   :: weight_occ_pattern_average_sorted    (N_occ_pattern)
integer, allocatable    :: psi_occ_pattern_sorted_order (N_occ_pattern)
integer, allocatable    :: psi_occ_pattern_sorted_order_reverse (N_occ_pattern)

Occupation patterns sorted by weight

Needs:

  • n_int

  • psi_occ_pattern

  • weight_occ_pattern_average

Needed by:

  • pruned

Subroutines / functions

a_operator:

File : determinants/slater_rules.irp.f

subroutine a_operator(iorb,ispin,key,hjj,Nint,na,nb)

Needed for diag_H_mat_elem().

Needs:

  • mo_one_e_integrals

  • mo_two_e_integrals_jj

Called by:

  • diag_h_mat_elem()

Calls:

  • bitstring_to_list_ab()

a_operator_two_e:

File : determinants/slater_rules_wee_mono.irp.f

subroutine a_operator_two_e(iorb,ispin,key,hjj,Nint,na,nb)

Needed for diag_Wee_mat_elem().

Needs:

  • mo_two_e_integrals_jj

Called by:

  • diag_wee_mat_elem()

Calls:

  • bitstring_to_list_ab()

ac_operator:

File : determinants/slater_rules.irp.f

subroutine ac_operator(iorb,ispin,key,hjj,Nint,na,nb)

Needed for diag_H_mat_elem().

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_num

  • mo_one_e_integrals

  • mo_two_e_integrals_jj

Called by:

  • diag_h_mat_elem()

Calls:

  • bitstring_to_list_ab()

ac_operator_two_e:

File : determinants/slater_rules_wee_mono.irp.f

subroutine ac_operator_two_e(iorb,ispin,key,hjj,Nint,na,nb)

Needed for diag_Wee_mat_elem().

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_two_e_integrals_jj

Called by:

  • diag_wee_mat_elem()

Calls:

  • bitstring_to_list_ab()

apply_excitation:

File : determinants/determinants.irp.f

subroutine apply_excitation(det, exc, res, ok, Nint)
apply_hole:

File : determinants/determinants.irp.f

subroutine apply_hole(det, s1, h1, res, ok, Nint)
apply_holes:

File : determinants/determinants.irp.f

subroutine apply_holes(det, s1, h1, s2, h2, res, ok, Nint)
apply_particle:

File : determinants/determinants.irp.f

subroutine apply_particle(det, s1, p1, res, ok, Nint)
apply_particles:

File : determinants/determinants.irp.f

subroutine apply_particles(det, s1, p1, s2, p2, res, ok, Nint)
bitstring_to_list_ab:

File : determinants/slater_rules.irp.f

subroutine bitstring_to_list_ab( string, list, n_elements, Nint)

Gives the inidices(+1) of the bits set to 1 in the bit string For alpha/beta determinants.

Called by:

  • a_operator()

  • a_operator_two_e()

  • ac_operator()

  • ac_operator_two_e()

  • build_fock_tmp()

  • diag_h_mat_elem()

  • diag_h_mat_elem_one_e()

  • diag_wee_mat_elem()

  • example_determinants()

  • fock_operator_closed_shell_ref_bitmask

  • fock_wee_closed_shell

  • get_occupation_from_dets()

  • get_single_excitation_from_fock()

  • i_h_j()

  • i_h_j_s2()

  • i_h_j_two_e()

  • i_h_j_verbose()

  • one_e_dm_mo_alpha

  • orb_range_diag_to_all_2_rdm_dm_buffer()

  • orb_range_diag_to_all_states_2_rdm_dm_buffer()

  • orb_range_off_diag_single_to_2_rdm_aa_dm_buffer()

  • orb_range_off_diag_single_to_2_rdm_ab_dm_buffer()

  • orb_range_off_diag_single_to_2_rdm_bb_dm_buffer()

  • orb_range_off_diag_single_to_all_states_aa_dm_buffer()

  • orb_range_off_diag_single_to_all_states_ab_dm_buffer()

  • orb_range_off_diag_single_to_all_states_bb_dm_buffer()

  • ref_closed_shell_bitmask

  • single_excitation_wee()

build_fock_tmp:

File : determinants/fock_diag.irp.f

subroutine build_fock_tmp(fock_diag_tmp,det_ref,Nint)

Build the diagonal of the Fock matrix corresponding to a generator determinant. $F_{00}$ is $langle i|H|i rangle = E_0$.

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_num

  • mo_one_e_integrals

  • mo_two_e_integrals_jj

  • n_int

Calls:

  • bitstring_to_list_ab()

  • debug_det()

build_singly_excited_wavefunction:

File : determinants/create_excitations.irp.f

subroutine build_singly_excited_wavefunction(i_hole,i_particle,ispin,det_out,coef_out)

Applies the single excitation operator : a^{dager}_(i_particle) a_(i_hole) of spin = ispin to the current wave function (psi_det, psi_coef)

Needs:

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det

Calls:

  • do_single_excitation()

  • get_phase()

connected_to_hf:

File : determinants/slater_rules.irp.f

subroutine connected_to_hf(key_i,yes_no)

Needs:

  • mo_one_e_integrals

  • n_int

  • ref_bitmask

  • thresh_sym

Calls:

  • get_excitation_degree()

  • get_single_excitation()

  • i_h_j()

connected_to_ref:

File : determinants/connected_to_ref.irp.f

integer function connected_to_ref(key,keys,Nint,N_past_in,Ndet)

input : key : a given Slater determinant

: keys: a list of Slater determinants

: Ndet: the number of Slater determinants in keys

: N_past_in the number of Slater determinants for the connectivity research

output : 0 : key not connected to the N_past_in first Slater determinants in keys

i : key is connected to determinant i of keys

-i : key is the ith determinant of the reference wf keys

Needs:

  • n_int

connected_to_ref_by_single:

File : determinants/connected_to_ref.irp.f

integer function connected_to_ref_by_single(key,keys,Nint,N_past_in,Ndet)

Returns true is key is connected to the reference by a single excitation. input : key : a given Slater determinant

: keys: a list of Slater determinants

: Ndet: the number of Slater determinants in keys

: N_past_in the number of Slater determinants for the connectivity research

output : 0 : key not connected by a MONO EXCITATION to the N_past_in first Slater determinants in keys

i : key is connected by a MONO EXCITATION to determinant i of keys

-i : key is the ith determinant of the reference wf keys

Needs:

  • n_int

copy_h_apply_buffer_to_wf:

File : determinants/h_apply.irp.f

subroutine copy_H_apply_buffer_to_wf

Copies the H_apply buffer to psi_coef. After calling this subroutine, N_det, psi_det and psi_coef need to be touched

Needs:

  • elec_alpha_num

  • elec_beta_num

  • h_apply_buffer_allocated

  • n_det

  • n_int

  • n_states

  • nproc

  • pruned

  • psi_coef

  • psi_det

  • psi_det_size

Called by:

  • generate_all_alpha_beta_det_products()

  • make_s2_eigenfunction()

Calls:

  • normalize()

  • remove_duplicates_in_psi_det()

Touches:

  • n_det

  • c0_weight

  • psi_coef

  • psi_det_sorted_bit

  • psi_det

  • psi_det_size

  • psi_det_sorted_bit

copy_psi_bilinear_to_psi:

File : determinants/spindeterminants.irp.f

subroutine copy_psi_bilinear_to_psi(psi, isize)

Overwrites psi_det and psi_coef with the wave function in bilinear order

Needs:

  • n_int

  • psi_bilinear_matrix_values

  • psi_det_alpha_unique

  • psi_det_beta_unique

create_microlist:

File : determinants/filter_connected.irp.f

subroutine create_microlist(minilist, N_minilist, key_mask, microlist, idx_microlist, N_microlist, ptr_microlist, Nint)

Needs:

  • mo_num

Calls:

  • bitstring_to_list()

create_minilist:

File : determinants/slater_rules.irp.f

subroutine create_minilist(key_mask, fullList, miniList, idx_miniList, N_fullList, N_miniList, Nint)
create_minilist_find_previous:

File : determinants/slater_rules.irp.f

subroutine create_minilist_find_previous(key_mask, fullList, miniList, N_fullList, N_miniList, fullMatch, Nint)
create_wf_of_psi_bilinear_matrix:

File : determinants/spindeterminants.irp.f

subroutine create_wf_of_psi_bilinear_matrix(truncate)

Generates a wave function containing all possible products of $alpha$ and $beta$ determinants

Needs:

  • n_det

  • n_int

  • n_states

  • psi_bilinear_matrix

  • psi_coef

  • psi_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • psi_det_sorted

  • psi_det_sorted_bit

Calls:

  • generate_all_alpha_beta_det_products()

Touches:

  • n_det

  • c0_weight

  • psi_coef

  • psi_det_sorted_bit

  • psi_det

  • psi_det_size

  • psi_det_sorted_bit

decode_exc:

File : determinants/slater_rules.irp.f

subroutine decode_exc(exc,degree,h1,p1,h2,p2,s1,s2)

Decodes the exc arrays returned by get_excitation. h1,h2 : Holes p1,p2 : Particles s1,s2 : Spins (1:alpha, 2:beta) degree : Degree of excitation

Called by:

  • diag_h_mat_elem_fock()

  • example_determinants()

decode_exc_spin:

File : determinants/slater_rules.irp.f

subroutine decode_exc_spin(exc,h1,p1,h2,p2)

Decodes the exc arrays returned by get_excitation.

h1,h2 : Holes

p1,p2 : Particles

Called by:

  • one_e_dm_mo_alpha

det_inf:

File : determinants/sort_dets_ab.irp.f

logical function det_inf(key1, key2, Nint)

Ordering function for determinants.

det_search_key:

File : determinants/connected_to_ref.irp.f

integer*8 function det_search_key(det,Nint)

Return an integer*8 corresponding to a determinant index for searching

Needs:

  • elec_alpha_num

detcmp:

File : determinants/determinants.irp.f

integer function detCmp(a,b,Nint)
deteq:

File : determinants/determinants.irp.f

logical function detEq(a,b,Nint)
diag_h_mat_elem:

File : determinants/slater_rules.irp.f

double precision function diag_H_mat_elem(det_in,Nint)

Computes $langle i|H|i rangle$.

Needs:

  • elec_alpha_num

  • elec_beta_num

  • elec_num

  • ref_bitmask

  • ref_bitmask_energy

Calls:

  • a_operator()

  • ac_operator()

  • bitstring_to_list_ab()

diag_h_mat_elem_fock:

File : determinants/slater_rules.irp.f

double precision function diag_H_mat_elem_fock(det_ref,det_pert,fock_diag_tmp,Nint)

Computes $langle i|H|i rangle$ when $i$ is at most a double excitation from a reference.

Needs:

  • mo_num

  • mo_two_e_integrals_jj

Calls:

  • decode_exc()

  • get_double_excitation()

  • get_excitation_degree()

  • get_single_excitation()

diag_h_mat_elem_one_e:

File : determinants/slater_rules_wee_mono.irp.f

double precision function diag_H_mat_elem_one_e(det_in,Nint)

Computes $langle i|H|i rangle$.

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_one_e_integrals

Calls:

  • bitstring_to_list_ab()

diag_s_mat_elem:

File : determinants/s2.irp.f

double precision function diag_S_mat_elem(key_i,Nint)

Returns <i|S^2|i>

diag_wee_mat_elem:

File : determinants/slater_rules_wee_mono.irp.f

double precision function diag_wee_mat_elem(det_in,Nint)

Computes $langle i|H|i rangle$.

Needs:

  • elec_alpha_num

  • elec_beta_num

  • elec_num

  • ref_bitmask

  • ref_bitmask_energy

Calls:

  • a_operator_two_e()

  • ac_operator_two_e()

  • bitstring_to_list_ab()

do_single_excitation:

File : determinants/create_excitations.irp.f

subroutine do_single_excitation(key_in,i_hole,i_particle,ispin,i_ok)

Apply the single excitation operator : a^{dager}_(i_particle) a_(i_hole) of spin = ispin on key_in ispin = 1 == alpha ispin = 2 == beta i_ok = 1 == the excitation is possible i_ok = -1 == the excitation is not possible

Needs:

  • mo_num

  • n_int

Called by:

  • build_singly_excited_wavefunction()

  • example_determinants()

example_determinants:

File : determinants/example.irp.f

subroutine example_determinants

subroutine that illustrates the main features available in determinants

Needs:

  • elec_alpha_num

  • n_int

  • ref_bitmask

Calls:

  • bitstring_to_list_ab()

  • debug_det()

  • decode_exc()

  • do_single_excitation()

  • get_excitation()

  • get_excitation_degree()

  • i_h_j()

  • print_det()

example_determinants_psi_det:

File : determinants/example.irp.f

subroutine example_determinants_psi_det

subroutine that illustrates the main features available in determinants using the psi_det/psi_coef

Needs:

  • read_wf

Calls:

  • routine_example_psi_det()

Touches:

  • read_wf

fill_h_apply_buffer_no_selection:

File : determinants/h_apply.irp.f

subroutine fill_H_apply_buffer_no_selection(n_selected,det_buffer,Nint,iproc)

Fill the H_apply buffer with determiants for CISD

Needs:

  • elec_alpha_num

  • elec_beta_num

  • h_apply_buffer_allocated

  • n_det

  • n_int

  • n_states

Called by:

  • generate_all_alpha_beta_det_products()

  • make_s2_eigenfunction()

Calls:

  • omp_set_lock()

  • omp_unset_lock()

  • resize_h_apply_buffer()

filter_connected:

File : determinants/filter_connected.irp.f

subroutine filter_connected(key1,key2,Nint,sze,idx)

Filters out the determinants that are not connected by H

returns the array idx which contains the index of the

determinants in the array key1 that interact

via the H operator with key2.

idx(0) is the number of determinants that interact with key1

Called by:

  • get_uj_s2_ui()

filter_connected_i_h_psi0:

File : determinants/filter_connected.irp.f

subroutine filter_connected_i_H_psi0(key1,key2,Nint,sze,idx)

Returns the array idx which contains the index of the

determinants in the array key1 that interact

via the H operator with key2.

idx(0) is the number of determinants that interact with key1

Needs:

  • n_int

Called by:

  • i_h_psi()

  • i_h_psi_minilist()

  • i_s2_psi_minilist()

filter_not_connected:

File : determinants/filter_connected.irp.f

subroutine filter_not_connected(key1,key2,Nint,sze,idx)

Returns the array idx which contains the index of the

determinants in the array key1 that DO NOT interact

via the H operator with key2.

idx(0) is the number of determinants that DO NOT interact with key1

generate_all_alpha_beta_det_products:

File : determinants/spindeterminants.irp.f

subroutine generate_all_alpha_beta_det_products

Creates a wave function from all possible $alpha times beta$ determinants

Needs:

  • h_apply_buffer_allocated

  • n_det

  • n_int

  • psi_coef

  • psi_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

Called by:

  • create_wf_of_psi_bilinear_matrix()

Calls:

  • copy_h_apply_buffer_to_wf()

  • fill_h_apply_buffer_no_selection()

Touches:

  • n_det

  • c0_weight

  • psi_coef

  • psi_det_sorted_bit

  • psi_det

  • psi_det_size

  • psi_det_sorted_bit

get_all_spin_doubles:

File : determinants/spindeterminants.irp.f

subroutine get_all_spin_doubles(buffer, idx, spindet, Nint, size_buffer, doubles, n_doubles)

Returns the indices of all the double excitations in the list of unique $alpha$ determinants.

Needs:

  • n_int

Calls:

  • get_all_spin_doubles_1()

  • get_all_spin_doubles_2()

  • get_all_spin_doubles_3()

  • get_all_spin_doubles_4()

  • get_all_spin_doubles_n_int()

get_all_spin_doubles_1:

File : determinants/spindeterminants.irp.f

subroutine get_all_spin_doubles_1(buffer, idx, spindet, size_buffer, doubles, n_doubles)

Returns the indices of all the double excitations in the list of unique $alpha$ determinants.

Called by:

  • get_all_spin_doubles()

get_all_spin_doubles_2:

File : determinants/spindeterminants.irp.f_template_1316

subroutine get_all_spin_doubles_2(buffer, idx, spindet, size_buffer, doubles, n_doubles)

Returns the indices of all the double excitations in the list of unique $lpha$ determinants.

Called by:

  • get_all_spin_doubles()

get_all_spin_doubles_3:

File : determinants/spindeterminants.irp.f_template_1316

subroutine get_all_spin_doubles_3(buffer, idx, spindet, size_buffer, doubles, n_doubles)

Returns the indices of all the double excitations in the list of unique $lpha$ determinants.

Called by:

  • get_all_spin_doubles()

get_all_spin_doubles_4:

File : determinants/spindeterminants.irp.f_template_1316

subroutine get_all_spin_doubles_4(buffer, idx, spindet, size_buffer, doubles, n_doubles)

Returns the indices of all the double excitations in the list of unique $lpha$ determinants.

Called by:

  • get_all_spin_doubles()

get_all_spin_doubles_n_int:

File : determinants/spindeterminants.irp.f_template_1316

subroutine get_all_spin_doubles_N_int(buffer, idx, spindet, size_buffer, doubles, n_doubles)

Returns the indices of all the double excitations in the list of unique $lpha$ determinants.

Needs:

  • n_int

Called by:

  • get_all_spin_doubles()

get_all_spin_singles:

File : determinants/spindeterminants.irp.f

subroutine get_all_spin_singles(buffer, idx, spindet, Nint, size_buffer, singles, n_singles)

Returns the indices of all the single excitations in the list of unique $alpha$ determinants.

Needs:

  • n_int

Called by:

  • singles_alpha_csc

  • singles_alpha_csc_idx

  • singles_beta_csc

  • singles_beta_csc_idx

Calls:

  • get_all_spin_singles_1()

  • get_all_spin_singles_2()

  • get_all_spin_singles_3()

  • get_all_spin_singles_4()

  • get_all_spin_singles_n_int()

get_all_spin_singles_1:

File : determinants/spindeterminants.irp.f

subroutine get_all_spin_singles_1(buffer, idx, spindet, size_buffer, singles, n_singles)

Returns the indices of all the single excitations in the list of unique $alpha$ determinants.

Called by:

  • get_all_spin_singles()

  • h_s2_u_0_nstates_openmp_work_1()

  • h_s2_u_0_two_e_nstates_openmp_work_1()

  • orb_range_2_rdm_openmp_work_1()

  • orb_range_2_rdm_state_av_openmp_work_1()

get_all_spin_singles_2:

File : determinants/spindeterminants.irp.f_template_1316

subroutine get_all_spin_singles_2(buffer, idx, spindet, size_buffer, singles, n_singles)

Returns the indices of all the single excitations in the list of unique $lpha$ determinants.

Called by:

  • get_all_spin_singles()

  • h_s2_u_0_nstates_openmp_work_2()

  • h_s2_u_0_two_e_nstates_openmp_work_2()

  • orb_range_2_rdm_openmp_work_2()

  • orb_range_2_rdm_state_av_openmp_work_2()

get_all_spin_singles_3:

File : determinants/spindeterminants.irp.f_template_1316

subroutine get_all_spin_singles_3(buffer, idx, spindet, size_buffer, singles, n_singles)

Returns the indices of all the single excitations in the list of unique $lpha$ determinants.

Called by:

  • get_all_spin_singles()

  • h_s2_u_0_nstates_openmp_work_3()

  • h_s2_u_0_two_e_nstates_openmp_work_3()

  • orb_range_2_rdm_openmp_work_3()

  • orb_range_2_rdm_state_av_openmp_work_3()

get_all_spin_singles_4:

File : determinants/spindeterminants.irp.f_template_1316

subroutine get_all_spin_singles_4(buffer, idx, spindet, size_buffer, singles, n_singles)

Returns the indices of all the single excitations in the list of unique $lpha$ determinants.

Called by:

  • get_all_spin_singles()

  • h_s2_u_0_nstates_openmp_work_4()

  • h_s2_u_0_two_e_nstates_openmp_work_4()

  • orb_range_2_rdm_openmp_work_4()

  • orb_range_2_rdm_state_av_openmp_work_4()

get_all_spin_singles_and_doubles:

File : determinants/spindeterminants.irp.f

subroutine get_all_spin_singles_and_doubles(buffer, idx, spindet, Nint, size_buffer, singles, doubles, n_singles, n_doubles)

Returns the indices of all the single and double excitations in the list of unique $alpha$ determinants.

Warning: The buffer is transposed.

Calls:

  • get_all_spin_singles_and_doubles_1()

  • get_all_spin_singles_and_doubles_2()

  • get_all_spin_singles_and_doubles_3()

  • get_all_spin_singles_and_doubles_4()

  • get_all_spin_singles_and_doubles_n_int()

get_all_spin_singles_and_doubles_1:

File : determinants/spindeterminants.irp.f

subroutine get_all_spin_singles_and_doubles_1(buffer, idx, spindet, size_buffer, singles, doubles, n_singles, n_doubles)

Returns the indices of all the single and double excitations in the list of unique $alpha$ determinants.

/!: The buffer is transposed !

Called by:

  • get_all_spin_singles_and_doubles()

  • h_s2_u_0_nstates_openmp_work_1()

  • h_s2_u_0_two_e_nstates_openmp_work_1()

  • orb_range_2_rdm_openmp_work_1()

  • orb_range_2_rdm_state_av_openmp_work_1()

get_all_spin_singles_and_doubles_2:

File : determinants/spindeterminants.irp.f_template_1316

subroutine get_all_spin_singles_and_doubles_2(buffer, idx, spindet, size_buffer, singles, doubles, n_singles, n_doubles)

Returns the indices of all the single and double excitations in the list of unique $lpha$ determinants.

/!: The buffer is transposed !

Called by:

  • get_all_spin_singles_and_doubles()

  • h_s2_u_0_nstates_openmp_work_2()

  • h_s2_u_0_two_e_nstates_openmp_work_2()

  • orb_range_2_rdm_openmp_work_2()

  • orb_range_2_rdm_state_av_openmp_work_2()

get_all_spin_singles_and_doubles_3:

File : determinants/spindeterminants.irp.f_template_1316

subroutine get_all_spin_singles_and_doubles_3(buffer, idx, spindet, size_buffer, singles, doubles, n_singles, n_doubles)

Returns the indices of all the single and double excitations in the list of unique $lpha$ determinants.

/!: The buffer is transposed !

Called by:

  • get_all_spin_singles_and_doubles()

  • h_s2_u_0_nstates_openmp_work_3()

  • h_s2_u_0_two_e_nstates_openmp_work_3()

  • orb_range_2_rdm_openmp_work_3()

  • orb_range_2_rdm_state_av_openmp_work_3()

get_all_spin_singles_and_doubles_4:

File : determinants/spindeterminants.irp.f_template_1316

subroutine get_all_spin_singles_and_doubles_4(buffer, idx, spindet, size_buffer, singles, doubles, n_singles, n_doubles)

Returns the indices of all the single and double excitations in the list of unique $lpha$ determinants.

/!: The buffer is transposed !

Called by:

  • get_all_spin_singles_and_doubles()

  • h_s2_u_0_nstates_openmp_work_4()

  • h_s2_u_0_two_e_nstates_openmp_work_4()

  • orb_range_2_rdm_openmp_work_4()

  • orb_range_2_rdm_state_av_openmp_work_4()

get_all_spin_singles_and_doubles_n_int:

File : determinants/spindeterminants.irp.f_template_1316

subroutine get_all_spin_singles_and_doubles_N_int(buffer, idx, spindet, size_buffer, singles, doubles, n_singles, n_doubles)

Returns the indices of all the single and double excitations in the list of unique $lpha$ determinants.

/!: The buffer is transposed !

Needs:

  • n_int

Called by:

  • get_all_spin_singles_and_doubles()

  • h_s2_u_0_nstates_openmp_work_n_int()

  • h_s2_u_0_two_e_nstates_openmp_work_n_int()

  • orb_range_2_rdm_openmp_work_n_int()

  • orb_range_2_rdm_state_av_openmp_work_n_int()

get_all_spin_singles_n_int:

File : determinants/spindeterminants.irp.f_template_1316

subroutine get_all_spin_singles_N_int(buffer, idx, spindet, size_buffer, singles, n_singles)

Returns the indices of all the single excitations in the list of unique $lpha$ determinants.

Needs:

  • n_int

Called by:

  • get_all_spin_singles()

  • h_s2_u_0_nstates_openmp_work_n_int()

  • h_s2_u_0_two_e_nstates_openmp_work_n_int()

  • orb_range_2_rdm_openmp_work_n_int()

  • orb_range_2_rdm_state_av_openmp_work_n_int()

get_double_excitation:

File : determinants/slater_rules.irp.f

subroutine get_double_excitation(det1,det2,exc,phase,Nint)

Returns the two excitation operators between two doubly excited determinants and the phase.

Called by:

  • diag_h_mat_elem_fock()

  • get_excitation()

  • get_s2()

  • i_h_j()

  • i_h_j_s2()

  • i_h_j_two_e()

  • i_h_j_verbose()

  • orb_range_off_diag_double_to_2_rdm_ab_dm_buffer()

  • orb_range_off_diag_double_to_all_states_ab_dm_buffer()

get_double_excitation_spin:

File : determinants/slater_rules.irp.f

subroutine get_double_excitation_spin(det1,det2,exc,phase,Nint)

Returns the two excitation operators between two doubly excited spin-determinants and the phase.

Called by:

  • get_excitation_spin()

  • i_h_j_double_spin()

  • orb_range_off_diag_double_to_2_rdm_aa_dm_buffer()

  • orb_range_off_diag_double_to_2_rdm_bb_dm_buffer()

  • orb_range_off_diag_double_to_all_states_aa_dm_buffer()

  • orb_range_off_diag_double_to_all_states_bb_dm_buffer()

get_excitation:

File : determinants/slater_rules.irp.f

subroutine get_excitation(det1,det2,exc,degree,phase,Nint)

Returns the excitation operators between two determinants and the phase.

Called by:

  • example_determinants()

  • get_phase()

Calls:

  • get_double_excitation()

  • get_excitation_degree()

  • get_single_excitation()

get_excitation_degree:

File : determinants/slater_rules.irp.f

subroutine get_excitation_degree(key1,key2,degree,Nint)

Returns the excitation degree between two determinants.

Called by:

  • connected_to_hf()

  • diag_h_mat_elem_fock()

  • example_determinants()

  • get_excitation()

  • get_s2()

  • i_h_j()

  • i_h_j_one_e()

  • i_h_j_s2()

  • i_h_j_two_e()

  • i_h_j_verbose()

  • max_degree_exc

  • psi_non_cas

get_excitation_degree_spin:

File : determinants/slater_rules.irp.f

subroutine get_excitation_degree_spin(key1,key2,degree,Nint)

Returns the excitation degree between two determinants.

Called by:

  • get_excitation_spin()

  • one_e_dm_mo_alpha

get_excitation_degree_vector:

File : determinants/slater_rules.irp.f

subroutine get_excitation_degree_vector(key1,key2,degree,Nint,sze,idx)

Applies get_excitation_degree to an array of determinants.

Called by:

  • routine_example_psi_det()

get_excitation_degree_vector_double_alpha_beta:

File : determinants/slater_rules.irp.f

subroutine get_excitation_degree_vector_double_alpha_beta(key1,key2,degree,Nint,sze,idx)

Applies get_excitation_degree to an array of determinants and return only the single excitations and the connections through exchange integrals.

get_excitation_degree_vector_single:

File : determinants/slater_rules.irp.f

subroutine get_excitation_degree_vector_single(key1,key2,degree,Nint,sze,idx)

Applies get_excitation_degree to an array of determinants and returns only the single excitations.

get_excitation_degree_vector_single_or_exchange:

File : determinants/slater_rules.irp.f

subroutine get_excitation_degree_vector_single_or_exchange(key1,key2,degree,Nint,sze,idx)

Applies get_excitation_degree to an array of determinants and return only the single excitations and the connections through exchange integrals.

get_excitation_degree_vector_single_or_exchange_verbose:

File : determinants/slater_rules.irp.f

subroutine get_excitation_degree_vector_single_or_exchange_verbose(key1,key2,degree,Nint,sze,idx)

Applies get_excitation_degree to an array of determinants and return only the single excitations and the connections through exchange integrals.

Needs:

  • n_int

Calls:

  • debug_det()

get_excitation_spin:

File : determinants/slater_rules.irp.f

subroutine get_excitation_spin(det1,det2,exc,degree,phase,Nint)

Returns the excitation operators between two determinants and the phase.

Calls:

  • get_double_excitation_spin()

  • get_excitation_degree_spin()

  • get_single_excitation_spin()

get_index_in_psi_det_alpha_unique:

File : determinants/spindeterminants.irp.f

integer function get_index_in_psi_det_alpha_unique(key,Nint)

Returns the index of the determinant in the psi_det_alpha_unique array

Needs:

  • psi_det_alpha_unique

get_index_in_psi_det_beta_unique:

File : determinants/spindeterminants.irp.f

integer function get_index_in_psi_det_beta_unique(key,Nint)

Returns the index of the determinant in the psi_det_beta_unique array

Needs:

  • psi_det_beta_unique

get_index_in_psi_det_sorted_bit:

File : determinants/connected_to_ref.irp.f

integer function get_index_in_psi_det_sorted_bit(key,Nint)

Returns the index of the determinant in the psi_det_sorted_bit array

Needs:

  • n_det

  • psi_det_sorted_bit

get_occupation_from_dets:

File : determinants/density_matrix.irp.f

subroutine get_occupation_from_dets(istate,occupation)

Returns the average occupation of the MOs

Needs:

  • mo_num

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det

Calls:

  • bitstring_to_list_ab()

get_phase:

File : determinants/slater_rules.irp.f

subroutine get_phase(key1,key2,phase,Nint)

Returns the phase between key1 and key2.

Called by:

  • build_singly_excited_wavefunction()

Calls:

  • get_excitation()

get_phasemask_bit:

File : determinants/slater_rules.irp.f

subroutine get_phasemask_bit(det1, pm, Nint)
get_s2:

File : determinants/s2.irp.f

subroutine get_s2(key_i,key_j,Nint,s2)

Returns $langle S^2 rangle - S_z^2 S_z$

Called by:

  • get_uj_s2_ui()

  • h_s2_u_0_nstates_openmp_work_1()

  • h_s2_u_0_nstates_openmp_work_2()

  • h_s2_u_0_nstates_openmp_work_3()

  • h_s2_u_0_nstates_openmp_work_4()

  • h_s2_u_0_nstates_openmp_work_n_int()

  • h_s2_u_0_two_e_nstates_openmp_work_1()

  • h_s2_u_0_two_e_nstates_openmp_work_2()

  • h_s2_u_0_two_e_nstates_openmp_work_3()

  • h_s2_u_0_two_e_nstates_openmp_work_4()

  • h_s2_u_0_two_e_nstates_openmp_work_n_int()

  • i_s2_psi_minilist()

  • s2_matrix_all_dets

  • s2_u_0_nstates()

Calls:

  • get_double_excitation()

  • get_excitation_degree()

get_single_excitation:

File : determinants/slater_rules.irp.f

subroutine get_single_excitation(det1,det2,exc,phase,Nint)

Returns the excitation operator between two singly excited determinants and the phase.

Called by:

  • connected_to_hf()

  • diag_h_mat_elem_fock()

  • get_excitation()

  • i_h_j()

  • i_h_j_one_e()

  • i_h_j_s2()

  • i_h_j_two_e()

  • i_h_j_verbose()

  • orb_range_off_diag_single_to_2_rdm_aa_dm_buffer()

  • orb_range_off_diag_single_to_2_rdm_ab_dm_buffer()

  • orb_range_off_diag_single_to_2_rdm_bb_dm_buffer()

  • orb_range_off_diag_single_to_all_states_aa_dm_buffer()

  • orb_range_off_diag_single_to_all_states_ab_dm_buffer()

  • orb_range_off_diag_single_to_all_states_bb_dm_buffer()

get_single_excitation_from_fock:

File : determinants/single_excitations.irp.f

subroutine get_single_excitation_from_fock(det_1,det_2,h,p,spin,phase,hij)

Needs:

  • big_array_coulomb_integrals

  • fock_operator_closed_shell_ref_bitmask

  • mo_num

  • n_int

  • ref_closed_shell_bitmask

Called by:

  • i_h_j()

  • i_h_j_s2()

  • i_h_j_single_spin()

Calls:

  • bitstring_to_list_ab()

get_single_excitation_spin:

File : determinants/slater_rules.irp.f

subroutine get_single_excitation_spin(det1,det2,exc,phase,Nint)

Returns the excitation operator between two singly excited determinants and the phase.

Called by:

  • get_excitation_spin()

  • i_h_j_double_alpha_beta()

  • i_h_j_mono_spin_one_e()

  • i_h_j_single_spin()

  • i_wee_j_single()

  • one_e_dm_mo_alpha

get_uj_s2_ui:

File : determinants/s2.irp.f

subroutine get_uJ_s2_uI(psi_keys_tmp,psi_coefs_tmp,n,nmax_coefs,nmax_keys,s2,nstates)

returns the matrix elements of S^2 “s2(i,j)” between the “nstates” states psi_coefs_tmp(:,i) and psi_coefs_tmp(:,j)

Needs:

  • n_int

Calls:

  • filter_connected()

  • get_s2()

getmobiles:

File : determinants/filter_connected.irp.f

subroutine getMobiles(key,key_mask, mobiles,Nint)

Needs:

  • mo_num

Calls:

  • bitstring_to_list()

i_h_j:

File : determinants/slater_rules.irp.f

subroutine i_H_j(key_i,key_j,Nint,hij)

Returns $langle i|H|j rangle$ where $i$ and $j$ are determinants.

Needs:

  • big_array_coulomb_integrals

  • big_array_coulomb_integrals

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_map

  • mo_two_e_integrals_in_map

  • n_int

Called by:

  • connected_to_hf()

  • example_determinants()

  • h_matrix_all_dets

  • h_matrix_cas

  • i_h_psi()

  • i_h_psi_minilist()

  • routine_example_psi_det()

Calls:

  • bitstring_to_list_ab()

  • get_double_excitation()

  • get_excitation_degree()

  • get_single_excitation()

  • get_single_excitation_from_fock()

i_h_j_double_alpha_beta:

File : determinants/slater_rules.irp.f

subroutine i_H_j_double_alpha_beta(key_i,key_j,Nint,hij)

Returns $langle i|H|j rangle$ where $i$ and $j$ are determinants differing by an opposite-spin double excitation.

Needs:

  • big_array_coulomb_integrals

  • big_array_coulomb_integrals

  • mo_integrals_map

  • mo_two_e_integrals_in_map

Called by:

  • h_s2_u_0_nstates_openmp_work_1()

  • h_s2_u_0_nstates_openmp_work_2()

  • h_s2_u_0_nstates_openmp_work_3()

  • h_s2_u_0_nstates_openmp_work_4()

  • h_s2_u_0_nstates_openmp_work_n_int()

  • h_s2_u_0_two_e_nstates_openmp_work_1()

  • h_s2_u_0_two_e_nstates_openmp_work_2()

  • h_s2_u_0_two_e_nstates_openmp_work_3()

  • h_s2_u_0_two_e_nstates_openmp_work_4()

  • h_s2_u_0_two_e_nstates_openmp_work_n_int()

Calls:

  • get_single_excitation_spin()

i_h_j_double_spin:

File : determinants/slater_rules.irp.f

subroutine i_H_j_double_spin(key_i,key_j,Nint,hij)

Returns $langle i|H|j rangle$ where $i$ and $j$ are determinants differing by a same-spin double excitation.

Needs:

  • big_array_coulomb_integrals

  • mo_integrals_map

  • mo_two_e_integrals_in_map

Called by:

  • h_s2_u_0_nstates_openmp_work_1()

  • h_s2_u_0_nstates_openmp_work_2()

  • h_s2_u_0_nstates_openmp_work_3()

  • h_s2_u_0_nstates_openmp_work_4()

  • h_s2_u_0_nstates_openmp_work_n_int()

  • h_s2_u_0_two_e_nstates_openmp_work_1()

  • h_s2_u_0_two_e_nstates_openmp_work_2()

  • h_s2_u_0_two_e_nstates_openmp_work_3()

  • h_s2_u_0_two_e_nstates_openmp_work_4()

  • h_s2_u_0_two_e_nstates_openmp_work_n_int()

Calls:

  • get_double_excitation_spin()

i_h_j_mono_spin_one_e:

File : determinants/slater_rules_wee_mono.irp.f

subroutine i_H_j_mono_spin_one_e(key_i,key_j,Nint,spin,hij)

Returns $langle i|H|j rangle$ where $i$ and $j$ are determinants differing by a single excitation.

Needs:

  • mo_one_e_integrals

Calls:

  • get_single_excitation_spin()

i_h_j_one_e:

File : determinants/slater_rules_wee_mono.irp.f

subroutine i_H_j_one_e(key_i,key_j,Nint,hij)

Returns $langle i|H|j rangle$ where $i$ and $j$ are determinants.

Needs:

  • mo_one_e_integrals

  • n_int

Calls:

  • get_excitation_degree()

  • get_single_excitation()

i_h_j_s2:

File : determinants/slater_rules.irp.f

subroutine i_H_j_s2(key_i,key_j,Nint,hij,s2)

Returns $langle i|H|j rangle$ and $langle i|S^2|j rangle$ where $i$ and $j$ are determinants.

Needs:

  • big_array_coulomb_integrals

  • big_array_coulomb_integrals

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_map

  • mo_two_e_integrals_in_map

  • n_int

Calls:

  • bitstring_to_list_ab()

  • get_double_excitation()

  • get_excitation_degree()

  • get_single_excitation()

  • get_single_excitation_from_fock()

i_h_j_single_spin:

File : determinants/slater_rules.irp.f

subroutine i_H_j_single_spin(key_i,key_j,Nint,spin,hij)

Returns $langle i|H|j rangle$ where $i$ and $j$ are determinants differing by a single excitation.

Needs:

  • big_array_coulomb_integrals

  • mo_two_e_integrals_in_map

Called by:

  • h_s2_u_0_nstates_openmp_work_1()

  • h_s2_u_0_nstates_openmp_work_2()

  • h_s2_u_0_nstates_openmp_work_3()

  • h_s2_u_0_nstates_openmp_work_4()

  • h_s2_u_0_nstates_openmp_work_n_int()

Calls:

  • get_single_excitation_from_fock()

  • get_single_excitation_spin()

i_h_j_two_e:

File : determinants/slater_rules_wee_mono.irp.f

subroutine i_H_j_two_e(key_i,key_j,Nint,hij)

Returns $langle i|H|j rangle$ where $i$ and $j$ are determinants.

Needs:

  • big_array_coulomb_integrals

  • big_array_coulomb_integrals

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_map

  • mo_two_e_integrals_in_map

  • n_int

  • ref_bitmask_energy

Calls:

  • bitstring_to_list_ab()

  • get_double_excitation()

  • get_excitation_degree()

  • get_single_excitation()

  • single_excitation_wee()

i_h_j_verbose:

File : determinants/slater_rules.irp.f

subroutine i_H_j_verbose(key_i,key_j,Nint,hij,hmono,hdouble,phase)

Returns $langle i|H|j rangle$ where $i$ and $j$ are determinants.

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mo_integrals_map

  • mo_one_e_integrals

  • mo_two_e_integrals_in_map

  • n_int

Calls:

  • bitstring_to_list_ab()

  • get_double_excitation()

  • get_excitation_degree()

  • get_single_excitation()

i_h_psi:

File : determinants/slater_rules.irp.f

subroutine i_H_psi(key,keys,coef,Nint,Ndet,Ndet_max,Nstate,i_H_psi_array)

Computes $langle i|H|Psi rangle = sum_J c_J langle i | H | J rangle$.

Uses filter_connected_i_H_psi0 to get all the $|J rangle$ to which $|i rangle$ is connected. The i_H_psi_minilist is much faster but requires to build the minilists.

Needs:

  • n_int

Calls:

  • filter_connected_i_h_psi0()

  • i_h_j()

i_h_psi_minilist:

File : determinants/slater_rules.irp.f

subroutine i_H_psi_minilist(key,keys,idx_key,N_minilist,coef,Nint,Ndet,Ndet_max,Nstate,i_H_psi_array)

Computes $langle i|H|Psi rangle = sum_J c_J langle i|H|Jrangle$.

Uses filter_connected_i_H_psi0 to get all the $|J rangle$ to which $|i rangle$ is connected. The $|Jrangle$ are searched in short pre-computed lists.

Needs:

  • n_int

Calls:

  • filter_connected_i_h_psi0()

  • i_h_j()

i_s2_psi_minilist:

File : determinants/s2.irp.f

subroutine i_S2_psi_minilist(key,keys,idx_key,N_minilist,coef,Nint,Ndet,Ndet_max,Nstate,i_S2_psi_array)

Computes $langle i|S^2|Psi rangle = sum_J c_J langle i|S^2|J rangle$.

Uses filter_connected_i_H_psi0 to get all the $|Jrangle$ to which $|irangle$ is connected. The $|Jrangle$ are searched in short pre-computed lists.

Needs:

  • n_int

Calls:

  • filter_connected_i_h_psi0()

  • get_s2()

i_wee_j_single:

File : determinants/slater_rules_wee_mono.irp.f

subroutine i_Wee_j_single(key_i,key_j,Nint,spin,hij)

Returns $langle i|H|j rangle$ where $i$ and $j$ are determinants differing by a single excitation.

Needs:

  • big_array_coulomb_integrals

  • mo_two_e_integrals_in_map

Called by:

  • h_s2_u_0_two_e_nstates_openmp_work_1()

  • h_s2_u_0_two_e_nstates_openmp_work_2()

  • h_s2_u_0_two_e_nstates_openmp_work_3()

  • h_s2_u_0_two_e_nstates_openmp_work_4()

  • h_s2_u_0_two_e_nstates_openmp_work_n_int()

Calls:

  • get_single_excitation_spin()

  • single_excitation_wee()

is_connected_to:

File : determinants/connected_to_ref.irp.f

logical function is_connected_to(key,keys,Nint,Ndet)

Returns true if determinant key is connected to keys

Needs:

  • n_int

is_connected_to_by_single:

File : determinants/connected_to_ref.irp.f

logical function is_connected_to_by_single(key,keys,Nint,Ndet)

Returns true is key is connected to keys by a single excitation.

Needs:

  • n_int

is_in_wavefunction:

File : determinants/connected_to_ref.irp.f

logical function is_in_wavefunction(key,Nint)

true if the determinant det is in the wave function

is_spin_flip_possible:

File : determinants/create_excitations.irp.f

logical function is_spin_flip_possible(key_in,i_flip,ispin)

returns true if the spin-flip of spin ispin in the orbital i_flip is possible on key_in

Needs:

  • n_int

make_s2_eigenfunction:

File : determinants/occ_pattern.irp.f

subroutine make_s2_eigenfunction

Needs:

  • elec_alpha_num

  • n_det

  • n_int

  • psi_occ_pattern

  • psi_coef

  • psi_det

  • psi_occ_pattern

Calls:

  • copy_h_apply_buffer_to_wf()

  • fill_h_apply_buffer_no_selection()

  • occ_pattern_to_dets()

  • occ_pattern_to_dets_size()

  • write_int()

  • write_time()

Touches:

  • n_det

  • psi_occ_pattern

  • c0_weight

  • psi_coef

  • psi_det_sorted_bit

  • psi_det

  • psi_det_size

  • psi_det_sorted_bit

  • psi_occ_pattern

occ_pattern_of_det:

File : determinants/occ_pattern.irp.f

subroutine occ_pattern_of_det(d,o,Nint)

Transforms a determinant to an occupation pattern

occ(:,1) : Single occupations

occ(:,2) : Double occupations

occ_pattern_search_key:

File : determinants/connected_to_ref.irp.f

integer*8 function occ_pattern_search_key(det,Nint)

Return an integer*8 corresponding to a determinant index for searching

Needs:

  • elec_alpha_num

occ_pattern_to_dets:

File : determinants/occ_pattern.irp.f

subroutine occ_pattern_to_dets(o,d,sze,n_alpha,Nint)

Generate all possible determinants for a given occ_pattern

Input :

o : occupation pattern : (doubly occupied, singly occupied) sze : Number of produced determinants, computed by occ_pattern_to_dets_size n_alpha : Number of $alpha$ electrons Nint : N_int

Output:

d : determinants

Needs:

  • binom_int

Called by:

  • make_s2_eigenfunction()

occ_pattern_to_dets_size:

File : determinants/occ_pattern.irp.f

subroutine occ_pattern_to_dets_size(o,sze,n_alpha,Nint)

Number of possible determinants for a given occ_pattern

Needs:

  • binom_int

Called by:

  • make_s2_eigenfunction()

read_dets:

File : determinants/determinants.irp.f

subroutine read_dets(det,Nint,Ndet)

Reads the determinants from the EZFIO file

Called by:

  • psi_det

Calls:

  • ezfio_get_determinants_bit_kind()

  • ezfio_get_determinants_n_int()

  • ezfio_get_determinants_psi_det()

read_spindeterminants:

File : determinants/spindeterminants.irp.f

subroutine read_spindeterminants

Needs:

  • n_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • n_states

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_det_alpha_unique

  • psi_det_beta_unique

Calls:

  • ezfio_get_spindeterminants_n_det()

  • ezfio_get_spindeterminants_n_det_alpha()

  • ezfio_get_spindeterminants_n_det_beta()

  • ezfio_get_spindeterminants_n_states()

  • ezfio_get_spindeterminants_psi_coef_matrix_columns()

  • ezfio_get_spindeterminants_psi_coef_matrix_rows()

  • ezfio_get_spindeterminants_psi_coef_matrix_values()

  • ezfio_get_spindeterminants_psi_det_alpha()

  • ezfio_get_spindeterminants_psi_det_beta()

  • wf_of_psi_bilinear_matrix()

Touches:

  • n_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • n_states

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_coef

  • psi_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

remove_duplicates_in_psi_det:

File : determinants/h_apply.irp.f

subroutine remove_duplicates_in_psi_det(found_duplicates)

Removes duplicate determinants in the wave function.

Needs:

  • c0_weight

  • n_det

  • n_int

  • psi_coef

  • psi_det_sorted_bit

  • psi_det

  • psi_det_sorted

  • psi_det_sorted_bit

Called by:

  • copy_h_apply_buffer_to_wf()

Touches:

  • n_det

  • c0_weight

  • psi_coef

  • psi_det_sorted_bit

  • psi_det

  • psi_det_sorted_bit

resize_h_apply_buffer:

File : determinants/h_apply.irp.f

subroutine resize_H_apply_buffer(new_size,iproc)

Resizes the H_apply buffer of proc iproc. The buffer lock should be set before calling this function.

Needs:

  • elec_alpha_num

  • elec_beta_num

  • h_apply_buffer_allocated

  • n_det

  • n_int

  • n_states

  • nproc

Called by:

  • fill_h_apply_buffer_no_selection()

routine_example_psi_det:

File : determinants/example.irp.f

subroutine routine_example_psi_det

subroutine that illustrates the main features available in determinants using many determinants

Needs:

  • n_det

  • n_int

  • n_states

  • psi_coef

  • psi_det

Called by:

  • example_determinants_psi_det()

Calls:

  • debug_det()

  • get_excitation_degree_vector()

  • i_h_j()

s2_u_0:

File : determinants/s2.irp.f

subroutine S2_u_0(v_0,u_0,n,keys_tmp,Nint)

Computes v_0 = S^2|u_0>

n : number of determinants

Calls:

  • s2_u_0_nstates()

s2_u_0_nstates:

File : determinants/s2.irp.f

subroutine S2_u_0_nstates(v_0,u_0,n,keys_tmp,Nint,N_st,sze_8)

Computes v_0 = S^2|u_0>

n : number of determinants

Needs:

  • n_int

  • ref_bitmask_energy

Called by:

  • s2_u_0()

  • u_0_s2_u_0()

Calls:

  • get_s2()

  • sort_dets_ab_v()

  • sort_dets_ba_v()

save_natural_mos:

File : determinants/density_matrix.irp.f

subroutine save_natural_mos

Save natural orbitals, obtained by diagonalization of the one-body density matrix in the MO basis

Needs:

  • ao_num

  • mo_coef

  • mo_num

Calls:

  • nullify_small_elements()

  • orthonormalize_mos()

  • save_mos()

  • set_natural_mos()

Touches:

  • mo_coef

  • mo_occ

save_ref_determinant:

File : determinants/determinants.irp.f

subroutine save_ref_determinant

Needs:

  • n_states

  • ref_bitmask

Calls:

  • save_wavefunction_general()

save_wavefunction:

File : determinants/determinants.irp.f

subroutine save_wavefunction

Save the wave function into the EZFIO file

Needs:

  • mpi_master

  • n_det

  • n_states

  • psi_det_sorted

  • read_wf

Calls:

  • save_wavefunction_general()

save_wavefunction_general:

File : determinants/determinants.irp.f

subroutine save_wavefunction_general(ndet,nstates,psidet,dim_psicoef,psicoef)

Save the wave function into the EZFIO file

Needs:

  • mo_label

  • mpi_master

  • n_det_qp_edit

  • n_int

Called by:

  • save_ref_determinant()

  • save_wavefunction()

  • save_wavefunction_truncated()

  • save_wavefunction_unsorted()

Calls:

  • ezfio_set_determinants_bit_kind()

  • ezfio_set_determinants_mo_label()

  • ezfio_set_determinants_n_det()

  • ezfio_set_determinants_n_det_qp_edit()

  • ezfio_set_determinants_n_int()

  • ezfio_set_determinants_n_states()

  • ezfio_set_determinants_psi_coef()

  • ezfio_set_determinants_psi_coef_qp_edit()

  • ezfio_set_determinants_psi_det()

  • ezfio_set_determinants_psi_det_qp_edit()

  • normalize()

  • write_int()

save_wavefunction_specified:

File : determinants/determinants.irp.f

subroutine save_wavefunction_specified(ndet,nstates,psidet,psicoef,ndetsave,index_det_save)

Save the wave function into the EZFIO file

Needs:

  • mo_label

  • mpi_master

  • n_det_qp_edit

  • n_int

Calls:

  • ezfio_set_determinants_bit_kind()

  • ezfio_set_determinants_mo_label()

  • ezfio_set_determinants_n_det()

  • ezfio_set_determinants_n_det_qp_edit()

  • ezfio_set_determinants_n_int()

  • ezfio_set_determinants_n_states()

  • ezfio_set_determinants_psi_coef()

  • ezfio_set_determinants_psi_det()

  • ezfio_set_determinants_psi_det_qp_edit()

  • write_int()

save_wavefunction_truncated:

File : determinants/determinants.irp.f

subroutine save_wavefunction_truncated(thr)

Save the wave function into the EZFIO file

Needs:

  • mpi_master

  • n_det

  • n_states

  • psi_coef

  • psi_det_sorted

Calls:

  • nullify_small_elements()

  • save_wavefunction_general()

Touches:

  • psi_coef

save_wavefunction_unsorted:

File : determinants/determinants.irp.f

subroutine save_wavefunction_unsorted

Save the wave function into the EZFIO file

Needs:

  • mpi_master

  • n_det

  • n_states

  • psi_coef

  • psi_det

Calls:

  • save_wavefunction_general()

set_natural_mos:

File : determinants/density_matrix.irp.f

subroutine set_natural_mos

Set natural orbitals, obtained by diagonalization of the one-body density matrix in the MO basis

Needs:

  • list_core_inact_act

  • list_virt

  • mo_num

  • mo_occ

  • n_core_inact_act_orb

  • n_virt_orb

  • one_e_dm_mo

Called by:

  • save_natural_mos()

Calls:

  • mo_as_svd_vectors_of_mo_matrix_eig()

Touches:

  • mo_occ

single_excitation_wee:

File : determinants/single_excitation_two_e.irp.f

subroutine single_excitation_wee(det_1,det_2,h,p,spin,phase,hij)

Needs:

  • big_array_coulomb_integrals

  • fock_wee_closed_shell

  • n_int

  • ref_closed_shell_bitmask

Called by:

  • i_h_j_two_e()

  • i_wee_j_single()

Calls:

  • bitstring_to_list_ab()

sort_dets_ab:

File : determinants/sort_dets_ab.irp.f

subroutine sort_dets_ab(key, idx, shortcut, N_key, Nint)

Deprecated routine

Calls:

  • tamiser()

sort_dets_ab_v:

File : determinants/sort_dets_ab.irp.f

subroutine sort_dets_ab_v(key_in, key_out, idx, shortcut, version, N_key, Nint)

Deprecated routine

Called by:

  • s2_u_0_nstates()

  • sort_dets_ba_v()

Calls:

  • tamiser()

sort_dets_ba_v:

File : determinants/sort_dets_ab.irp.f

subroutine sort_dets_ba_v(key_in, key_out, idx, shortcut, version, N_key, Nint)

Deprecated routine

Called by:

  • s2_u_0_nstates()

Calls:

  • sort_dets_ab_v()

sort_dets_by_det_search_key:

File : determinants/determinants.irp.f

subroutine sort_dets_by_det_search_key(Ndet, det_in, coef_in, sze, det_out, coef_out, N_st)

Determinants are sorted according to their det_search_key(). Useful to accelerate the search of a random determinant in the wave function.

/!The first dimension of coef_out and coef_in need to be psi_det_size

Needs:

  • n_int

Called by:

  • psi_cas_sorted_bit

  • psi_det_sorted_bit

  • psi_non_cas_sorted_bit

Calls:

  • i8sort()

spin_det_search_key:

File : determinants/spindeterminants.irp.f

integer*8 function spin_det_search_key(det,Nint)

Returns an integer(8) corresponding to a determinant index for searching

tamiser:

File : determinants/sort_dets_ab.irp.f

subroutine tamiser(key, idx, no, n, Nint, N_key)

Called by:

  • sort_dets_ab()

  • sort_dets_ab_v()

u_0_s2_u_0:

File : determinants/s2.irp.f

subroutine u_0_S2_u_0(e_0,u_0,n,keys_tmp,Nint,N_st,sze_8)

Computes e_0 = <u_0|S2|u_0>/<u_0|u_0>

n : number of determinants

Called by:

  • ci_electronic_energy

  • s2_values

Calls:

  • s2_u_0_nstates()

wf_of_psi_bilinear_matrix:

File : determinants/spindeterminants.irp.f

subroutine wf_of_psi_bilinear_matrix(truncate)

Generate a wave function containing all possible products of $alpha$ and $beta$ determinants

Needs:

  • n_det

  • n_int

  • n_states

  • psi_bilinear_matrix_values

  • psi_coef

  • psi_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • psi_det_sorted

Called by:

  • read_spindeterminants()

Touches:

  • n_det

  • psi_coef

  • psi_det

write_spindeterminants:

File : determinants/spindeterminants.irp.f

subroutine write_spindeterminants

Needs:

  • n_det

  • n_int

  • n_states

  • psi_bilinear_matrix_values

  • psi_det_alpha_unique

  • psi_det_beta_unique

Calls:

  • ezfio_set_spindeterminants_bit_kind()

  • ezfio_set_spindeterminants_n_det()

  • ezfio_set_spindeterminants_n_det_alpha()

  • ezfio_set_spindeterminants_n_det_beta()

  • ezfio_set_spindeterminants_n_int()

  • ezfio_set_spindeterminants_n_states()

  • ezfio_set_spindeterminants_psi_coef_matrix_columns()

  • ezfio_set_spindeterminants_psi_coef_matrix_rows()

  • ezfio_set_spindeterminants_psi_coef_matrix_values()

  • ezfio_set_spindeterminants_psi_det_alpha()

  • ezfio_set_spindeterminants_psi_det_beta()

zmq_get_n_det:

File : determinants/zmq.irp.f_template_379

integer function zmq_get_N_det(zmq_to_qp_run_socket, worker_id)

Get N_det from the qp_run scheduler

Needs:

  • mpi_master

  • n_det

  • zmq_state

zmq_get_n_det_alpha_unique:

File : determinants/zmq.irp.f_template_379

integer function zmq_get_N_det_alpha_unique(zmq_to_qp_run_socket, worker_id)

Get N_det_alpha_unique from the qp_run scheduler

Needs:

  • mpi_master

  • psi_det_alpha_unique

  • zmq_state

zmq_get_n_det_beta_unique:

File : determinants/zmq.irp.f_template_379

integer function zmq_get_N_det_beta_unique(zmq_to_qp_run_socket, worker_id)

Get N_det_beta_unique from the qp_run scheduler

Needs:

  • mpi_master

  • psi_det_beta_unique

  • zmq_state

zmq_get_n_states:

File : determinants/zmq.irp.f_template_379

integer function zmq_get_N_states(zmq_to_qp_run_socket, worker_id)

Get N_states from the qp_run scheduler

Needs:

  • mpi_master

  • n_states

  • zmq_state

zmq_get_psi:

File : determinants/zmq.irp.f

integer function zmq_get_psi(zmq_to_qp_run_socket, worker_id)

Get the wave function from the qp_run scheduler

Needs:

  • n_det

  • n_states

  • psi_coef

  • psi_det

  • psi_det_size

Touches:

  • n_det

  • n_states

  • psi_coef

  • psi_det

  • psi_det_size

zmq_get_psi_bilinear:

File : determinants/zmq.irp.f

integer function zmq_get_psi_bilinear(zmq_to_qp_run_socket, worker_id)

Get the wave function from the qp_run scheduler

Needs:

  • n_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • n_states

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_coef

  • psi_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • psi_det_size

Touches:

  • n_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • n_states

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_coef

  • psi_det

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • psi_det_size

zmq_get_psi_bilinear_matrix_columns:

File : determinants/zmq.irp.f_template_500

integer*8 function zmq_get_psi_bilinear_matrix_columns(zmq_to_qp_run_socket,worker_id)

Get psi_bilinear_matrix_columns on the qp_run scheduler

Needs:

  • psi_bilinear_matrix_values

zmq_get_psi_bilinear_matrix_order:

File : determinants/zmq.irp.f_template_500

integer*8 function zmq_get_psi_bilinear_matrix_order(zmq_to_qp_run_socket,worker_id)

Get psi_bilinear_matrix_order on the qp_run scheduler

Needs:

  • psi_bilinear_matrix_values

zmq_get_psi_bilinear_matrix_rows:

File : determinants/zmq.irp.f_template_500

integer*8 function zmq_get_psi_bilinear_matrix_rows(zmq_to_qp_run_socket,worker_id)

Get psi_bilinear_matrix_rows on the qp_run scheduler

Needs:

  • psi_bilinear_matrix_values

zmq_get_psi_bilinear_matrix_values:

File : determinants/zmq.irp.f_template_564

integer*8 function zmq_get_psi_bilinear_matrix_values(zmq_to_qp_run_socket,worker_id)

get psi_bilinear_matrix_values on the qp_run scheduler

Needs:

  • psi_bilinear_matrix_values

zmq_get_psi_coef:

File : determinants/zmq.irp.f_template_564

integer*8 function zmq_get_psi_coef(zmq_to_qp_run_socket,worker_id)

get psi_coef on the qp_run scheduler

Needs:

  • psi_coef

zmq_get_psi_det:

File : determinants/zmq.irp.f_template_440

integer*8 function zmq_get_psi_det(zmq_to_qp_run_socket,worker_id)

Get psi_det on the qp_run scheduler

Needs:

  • psi_det

zmq_get_psi_det_alpha_unique:

File : determinants/zmq.irp.f_template_440

integer*8 function zmq_get_psi_det_alpha_unique(zmq_to_qp_run_socket,worker_id)

Get psi_det_alpha_unique on the qp_run scheduler

Needs:

  • psi_det_alpha_unique

zmq_get_psi_det_beta_unique:

File : determinants/zmq.irp.f_template_440

integer*8 function zmq_get_psi_det_beta_unique(zmq_to_qp_run_socket,worker_id)

Get psi_det_beta_unique on the qp_run scheduler

Needs:

  • psi_det_beta_unique

zmq_get_psi_det_size:

File : determinants/zmq.irp.f_template_379

integer function zmq_get_psi_det_size(zmq_to_qp_run_socket, worker_id)

Get psi_det_size from the qp_run scheduler

Needs:

  • mpi_master

  • psi_det_size

  • zmq_state

zmq_get_psi_notouch:

File : determinants/zmq.irp.f

integer function zmq_get_psi_notouch(zmq_to_qp_run_socket, worker_id)

Get the wave function from the qp_run scheduler

Needs:

  • n_int

  • n_states

  • psi_coef

  • psi_det

  • psi_det_size

zmq_put_n_det:

File : determinants/zmq.irp.f_template_379

integer function zmq_put_N_det(zmq_to_qp_run_socket,worker_id)

Put N_det on the qp_run scheduler

Needs:

  • n_det

  • zmq_state

zmq_put_n_det_alpha_unique:

File : determinants/zmq.irp.f_template_379

integer function zmq_put_N_det_alpha_unique(zmq_to_qp_run_socket,worker_id)

Put N_det_alpha_unique on the qp_run scheduler

Needs:

  • psi_det_alpha_unique

  • zmq_state

zmq_put_n_det_beta_unique:

File : determinants/zmq.irp.f_template_379

integer function zmq_put_N_det_beta_unique(zmq_to_qp_run_socket,worker_id)

Put N_det_beta_unique on the qp_run scheduler

Needs:

  • psi_det_beta_unique

  • zmq_state

zmq_put_n_states:

File : determinants/zmq.irp.f_template_379

integer function zmq_put_N_states(zmq_to_qp_run_socket,worker_id)

Put N_states on the qp_run scheduler

Needs:

  • n_states

  • zmq_state

zmq_put_psi:

File : determinants/zmq.irp.f

integer function zmq_put_psi(zmq_to_qp_run_socket,worker_id)

Put the wave function on the qp_run scheduler

zmq_put_psi_bilinear:

File : determinants/zmq.irp.f

integer function zmq_put_psi_bilinear(zmq_to_qp_run_socket,worker_id)

Put the wave function on the qp_run scheduler

zmq_put_psi_bilinear_matrix_columns:

File : determinants/zmq.irp.f_template_500

integer*8 function zmq_put_psi_bilinear_matrix_columns(zmq_to_qp_run_socket,worker_id)

Put psi_bilinear_matrix_columns on the qp_run scheduler

Needs:

  • psi_bilinear_matrix_values

zmq_put_psi_bilinear_matrix_order:

File : determinants/zmq.irp.f_template_500

integer*8 function zmq_put_psi_bilinear_matrix_order(zmq_to_qp_run_socket,worker_id)

Put psi_bilinear_matrix_order on the qp_run scheduler

Needs:

  • psi_bilinear_matrix_values

zmq_put_psi_bilinear_matrix_rows:

File : determinants/zmq.irp.f_template_500

integer*8 function zmq_put_psi_bilinear_matrix_rows(zmq_to_qp_run_socket,worker_id)

Put psi_bilinear_matrix_rows on the qp_run scheduler

Needs:

  • psi_bilinear_matrix_values

zmq_put_psi_bilinear_matrix_values:

File : determinants/zmq.irp.f_template_564

integer*8 function zmq_put_psi_bilinear_matrix_values(zmq_to_qp_run_socket,worker_id)

Put psi_bilinear_matrix_values on the qp_run scheduler

Needs:

  • psi_bilinear_matrix_values

zmq_put_psi_coef:

File : determinants/zmq.irp.f_template_564

integer*8 function zmq_put_psi_coef(zmq_to_qp_run_socket,worker_id)

Put psi_coef on the qp_run scheduler

Needs:

  • psi_coef

zmq_put_psi_det:

File : determinants/zmq.irp.f_template_440

integer*8 function zmq_put_psi_det(zmq_to_qp_run_socket,worker_id)

Put psi_det on the qp_run scheduler

Needs:

  • psi_det

zmq_put_psi_det_alpha_unique:

File : determinants/zmq.irp.f_template_440

integer*8 function zmq_put_psi_det_alpha_unique(zmq_to_qp_run_socket,worker_id)

Put psi_det_alpha_unique on the qp_run scheduler

Needs:

  • psi_det_alpha_unique

zmq_put_psi_det_beta_unique:

File : determinants/zmq.irp.f_template_440

integer*8 function zmq_put_psi_det_beta_unique(zmq_to_qp_run_socket,worker_id)

Put psi_det_beta_unique on the qp_run scheduler

Needs:

  • psi_det_beta_unique

zmq_put_psi_det_size:

File : determinants/zmq.irp.f_template_379

integer function zmq_put_psi_det_size(zmq_to_qp_run_socket,worker_id)

Put psi_det_size on the qp_run scheduler

Needs:

  • psi_det_size

  • zmq_state