cipsi

CIPSI algorithm.

The run_stochastic_cipsi() and run_cipsi() subroutines start with a single determinant, or with the wave function in the EZFIO database if determinants read_wf is true.

The run_cipsi() subroutine iteratively:

  • Selects the most important determinants from the external space and adds them to the internal space

  • If determinants s2_eig is true, it adds all the necessary determinants to allow the eigenstates of \(\hat H\) to be eigenstates of \(\widehat{S^2}\)

  • Diagonalizes \(\hat H\) in the enlarged internal space

  • Computes the PT2 contribution to the energy stochastically [22] or deterministically, depending on perturbation do_pt2

  • Extrapolates the variational energy by fitting \(E=E_\text{FCI} - \alpha\, E_\text{PT2}\)

The difference between run_stochastic_cipsi() and run_cipsi() is that run_stochastic_cipsi() selects the determinants on the fly with the computation of the stochastic PT2 [6]. Hence, it is a semi-stochastic selection. It

  • Selects the most important determinants from the external space and adds them to the internal space, on the fly with the computation of the PT2 with the stochastic algorithm presented in [6].

  • If determinants s2_eig is true, it adds all the necessary determinants to allow the eigenstates of \(\hat H\) to be eigenstates of \(\widehat{S^2}\)

  • Extrapolates the variational energy by fitting \(E=E_\text{FCI} - \alpha\, E_\text{PT2}\)

  • Diagonalizes \(\hat H\) in the enlarged internal space

The number of selected determinants at each iteration will be such that the size of the wave function will double at every iteration. If determinants s2_eig is true, then the number of selected determinants will be 1.5x the current number, and then all the additional determinants will be added.

By default, the program will stop when more than one million determinants have been selected, or when the PT2 energy is below \(10^{-4}\).

The variational and PT2 energies of the iterations are stored in the EZFIO database, in the iterations module.

Computation of the PT2 energy

At each iteration, the PT2 energy is computed considering the Epstein-Nesbet zeroth-order Hamiltonian:

\[E_{\text{PT2}} = \sum_{ \alpha } \frac{|\langle \Psi_S | \hat{H} | \alpha \rangle|^2} {E - \langle \alpha | \hat{H} | \alpha \rangle}\]

where the \(|\alpha \rangle\) determinants are generated by applying all the single and double excitation operators to all the determinants of the wave function \(\Psi_G\).

When the hybrid-deterministic/stochastic algorithm is chosen (default), \(Psi_G = \Psi_S = \Psi\), the full wavefunction expanded in the internal space. When the deterministic algorithm is chosen (perturbation do_pt2 is set to false), \(Psi_G\) is a truncation of \(|\Psi \rangle\) using determinants threshold_generators, and \(Psi_S\) is a truncation of \(|\Psi \rangle\) using determinants threshold_selectors, and re-weighted by \(1/\langle \Psi_s | \Psi_s \rangle\).

At every iteration, while computing the PT2, the variance of the wave function is also computed:

\[\begin{split}\sigma^2 & = \langle \Psi | \hat{H}^2 | \Psi \rangle - \langle \Psi | \hat{H} | \Psi \rangle^2 \\ & = \sum_{i \in \text{FCI}} \langle \Psi | \hat{H} | i \rangle \langle i | \hat{H} | \Psi \rangle - \langle \Psi | \hat{H} | \Psi \rangle^2 \\ & = \sum_{ \alpha } \langle |\Psi | \hat{H} | \alpha \rangle|^2.\end{split}\]

The expression of the variance is the same as the expression of the PT2, with a denominator of 1. It measures how far the wave function is from the FCI solution. Note that the absence of denominator in the Heat-Bath selected CI method is selection method by minimization of the variance, whereas CIPSI is a selection method by minimization of the energy.

If perturbation do_pt2 is set to false, then the stochastic PT2 is not computed, and an approximate value is obtained from the CIPSI selection. The calculation is faster, but the extrapolated FCI value is less accurate. This way of running the code should be used when the only goal is to generate a wave function, as for using CIPSI wave functions as trial wave functions of QMC calculations for example.

The PT2 program reads the wave function of the EZFIO database and computes the energy and the PT2 contribution.

State-averaging

Extrapolated FCI energy

An estimate of the FCI energy is computed by extrapolating

\[E=E_\text{FCI} - \alpha\, E_\text{PT2}\]

This extrapolation is done for all the requested states, and excitation energies are printed as energy differences between the extrapolated energies of the excited states and the extrapolated energy of the ground state.

The extrapolations are given considering the 2 last points, the 3 last points, …, the 7 last points. The extrapolated value should be chosen such that the extrpolated value is stable with the number of points.

EZFIO parameters

pert_2rdm

If true, computes the one- and two-body rdms with perturbation theory

Default: False

Providers

global_selection_buffer

File : cipsi/run_pt2_slave.irp.f

type(selection_buffer)  :: global_selection_buffer

Global buffer for the OpenMP selection

Needs:

  • global_selection_buffer_lock

  • n_det_generators

  • n_int

global_selection_buffer_lock

File : cipsi/run_pt2_slave.irp.f

integer(omp_lock_kind)  :: global_selection_buffer_lock

Global buffer for the OpenMP selection

Needed by:

  • global_selection_buffer

initialize_pt2_e0_denominator

File : cipsi/energy.irp.f

logical :: initialize_pt2_e0_denominator

If true, initialize pt2_E0_denominator

Needed by:

  • pt2_e0_denominator

list_orb_pert_rdm

File : cipsi/pert_rdm_providers.irp.f

integer, allocatable    :: list_orb_pert_rdm    (n_orb_pert_rdm)

Needs:

  • list_act

  • n_orb_pert_rdm

list_orb_reverse_pert_rdm

File : cipsi/pert_rdm_providers.irp.f

integer, allocatable    :: list_orb_reverse_pert_rdm    (mo_num)

Needs:

  • list_act

  • mo_num

n_orb_pert_rdm

File : cipsi/pert_rdm_providers.irp.f

integer :: n_orb_pert_rdm

Needs:

  • n_act_orb

Needed by:

  • list_orb_pert_rdm

  • pert_2rdm_provider

nthreads_pt2

File : cipsi/environment.irp.f

integer :: nthreads_pt2

Number of threads for Davidson

Needs:

  • mpi_master

  • nproc

pert_2rdm_lock

File : cipsi/pert_rdm_providers.irp.f

integer(omp_lock_kind)  :: pert_2rdm_lock
pert_2rdm_provider

File : cipsi/pert_rdm_providers.irp.f

double precision, allocatable   :: pert_2rdm_provider   (n_orb_pert_rdm,n_orb_pert_rdm,n_orb_pert_rdm,n_orb_pert_rdm)

Needs:

  • n_orb_pert_rdm

pt2_cw

File : cipsi/pt2_stoch_routines.irp.f

double precision, allocatable   :: pt2_w        (N_det_generators)
double precision, allocatable   :: pt2_cw       (0:N_det_generators)
double precision        :: pt2_w_t
double precision        :: pt2_u_0
integer, allocatable    :: pt2_n_0      (pt2_N_teeth+1)

Needs:

  • n_det_generators

  • psi_det_sorted_gen

  • pt2_n_teeth

  • pt2_stoch_istate

  • qp_max_mem

Needed by:

  • pt2_f

  • pt2_j

pt2_e0_denominator

File : cipsi/energy.irp.f

double precision, allocatable   :: pt2_e0_denominator   (N_states)

E0 in the denominator of the PT2

Needs:

  • barycentric_electronic_energy

  • h0_type

  • initialize_pt2_e0_denominator

  • mpi_master

  • n_states

  • nuclear_repulsion

  • psi_coef

  • psi_det_hii

  • psi_energy

pt2_f

File : cipsi/pt2_stoch_routines.irp.f

integer, allocatable    :: pt2_f        (N_det_generators)
integer :: pt2_n_tasks_max

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mpi_master

  • n_core_orb

  • n_det_generators

  • pt2_min_parallel_tasks

  • pt2_n_teeth

  • pt2_w

pt2_j

File : cipsi/pt2_stoch_routines.irp.f

integer, allocatable    :: pt2_j        (N_det_generators)
integer, allocatable    :: pt2_r        (N_det_generators)

pt2_J contains the list of generators after ordering them according to the Monte Carlo sampling.

pt2_R(i) is the number of combs drawn when determinant i is computed.

Needs:

  • n_det_generators

  • pt2_n_tasks

  • pt2_n_teeth

  • pt2_u

  • pt2_w

  • qp_max_mem

pt2_match_weight

File : cipsi/selection.irp.f

double precision, allocatable   :: pt2_match_weight     (N_states)

Weights adjusted along the selection to make the PT2 contributions of each state coincide.

Needs:

  • n_states

Needed by:

  • selection_weight

pt2_mindetinfirstteeth

File : cipsi/pt2_stoch_routines.irp.f

integer :: pt2_n_teeth
integer :: pt2_mindetinfirstteeth

Needs:

  • mpi_master

  • n_det_generators

  • psi_det_sorted_gen

  • pt2_stoch_istate

Needed by:

  • pt2_f

  • pt2_j

  • pt2_w

pt2_n_0

File : cipsi/pt2_stoch_routines.irp.f

double precision, allocatable   :: pt2_w        (N_det_generators)
double precision, allocatable   :: pt2_cw       (0:N_det_generators)
double precision        :: pt2_w_t
double precision        :: pt2_u_0
integer, allocatable    :: pt2_n_0      (pt2_N_teeth+1)

Needs:

  • n_det_generators

  • psi_det_sorted_gen

  • pt2_n_teeth

  • pt2_stoch_istate

  • qp_max_mem

Needed by:

  • pt2_f

  • pt2_j

pt2_n_tasks

File : cipsi/pt2_stoch_routines.irp.f

integer :: pt2_n_tasks

Number of parallel tasks for the Monte Carlo

Needs:

  • n_det_generators

Needed by:

  • pt2_j

pt2_n_tasks_max

File : cipsi/pt2_stoch_routines.irp.f

integer, allocatable    :: pt2_f        (N_det_generators)
integer :: pt2_n_tasks_max

Needs:

  • elec_alpha_num

  • elec_beta_num

  • mpi_master

  • n_core_orb

  • n_det_generators

  • pt2_min_parallel_tasks

  • pt2_n_teeth

  • pt2_w

pt2_n_teeth

File : cipsi/pt2_stoch_routines.irp.f

integer :: pt2_n_teeth
integer :: pt2_mindetinfirstteeth

Needs:

  • mpi_master

  • n_det_generators

  • psi_det_sorted_gen

  • pt2_stoch_istate

Needed by:

  • pt2_f

  • pt2_j

  • pt2_w

pt2_overlap

File : cipsi/energy.irp.f

double precision, allocatable   :: pt2_overlap  (N_states,N_states)

Overlap between the perturbed wave functions

Needs:

  • n_states

pt2_r

File : cipsi/pt2_stoch_routines.irp.f

integer, allocatable    :: pt2_j        (N_det_generators)
integer, allocatable    :: pt2_r        (N_det_generators)

pt2_J contains the list of generators after ordering them according to the Monte Carlo sampling.

pt2_R(i) is the number of combs drawn when determinant i is computed.

Needs:

  • n_det_generators

  • pt2_n_tasks

  • pt2_n_teeth

  • pt2_u

  • pt2_w

  • qp_max_mem

pt2_stoch_istate

File : cipsi/pt2_stoch_routines.irp.f

integer :: pt2_stoch_istate

State for stochatsic PT2

Needed by:

  • pt2_n_teeth

  • pt2_w

pt2_u

File : cipsi/pt2_stoch_routines.irp.f

double precision, allocatable   :: pt2_u        (N_det_generators)

Needs:

  • n_det_generators

Needed by:

  • pt2_j

pt2_u_0

File : cipsi/pt2_stoch_routines.irp.f

double precision, allocatable   :: pt2_w        (N_det_generators)
double precision, allocatable   :: pt2_cw       (0:N_det_generators)
double precision        :: pt2_w_t
double precision        :: pt2_u_0
integer, allocatable    :: pt2_n_0      (pt2_N_teeth+1)

Needs:

  • n_det_generators

  • psi_det_sorted_gen

  • pt2_n_teeth

  • pt2_stoch_istate

  • qp_max_mem

Needed by:

  • pt2_f

  • pt2_j

pt2_w

File : cipsi/pt2_stoch_routines.irp.f

double precision, allocatable   :: pt2_w        (N_det_generators)
double precision, allocatable   :: pt2_cw       (0:N_det_generators)
double precision        :: pt2_w_t
double precision        :: pt2_u_0
integer, allocatable    :: pt2_n_0      (pt2_N_teeth+1)

Needs:

  • n_det_generators

  • psi_det_sorted_gen

  • pt2_n_teeth

  • pt2_stoch_istate

  • qp_max_mem

Needed by:

  • pt2_f

  • pt2_j

pt2_w_t

File : cipsi/pt2_stoch_routines.irp.f

double precision, allocatable   :: pt2_w        (N_det_generators)
double precision, allocatable   :: pt2_cw       (0:N_det_generators)
double precision        :: pt2_w_t
double precision        :: pt2_u_0
integer, allocatable    :: pt2_n_0      (pt2_N_teeth+1)

Needs:

  • n_det_generators

  • psi_det_sorted_gen

  • pt2_n_teeth

  • pt2_stoch_istate

  • qp_max_mem

Needed by:

  • pt2_f

  • pt2_j

selection_weight

File : cipsi/selection.irp.f

double precision, allocatable   :: selection_weight     (N_states)

Weights used in the selection criterion

Needs:

  • c0_weight

  • n_states

  • pt2_match_weight

  • state_average_weight

  • variance_match_weight

  • weight_selection

variance_match_weight

File : cipsi/selection.irp.f

double precision, allocatable   :: variance_match_weight        (N_states)

Weights adjusted along the selection to make the variances of each state coincide.

Needs:

  • n_states

Needed by:

  • selection_weight

Subroutines / functions

add_to_selection_buffer:

File : cipsi/selection_buffer.irp.f

subroutine add_to_selection_buffer(b, det, val)

Needs:

  • n_int

Called by:

  • fill_buffer_double()

  • fill_buffer_double_rdm()

  • pt2_collector()

  • selection_collector()

Calls:

  • sort_selection_buffer()

bitstring_to_list_in_selection:

File : cipsi/selection.irp.f

subroutine bitstring_to_list_in_selection( string, list, n_elements, Nint)

Gives the inidices(+1) of the bits set to 1 in the bit string

Called by:

  • splash_pq()

  • spot_isinwf()

create_selection_buffer:

File : cipsi/selection_buffer.irp.f

subroutine create_selection_buffer(N, size_in, res)

Allocates the memory for a selection buffer. The arrays have dimension size_in and the maximum number of elements is N

Needs:

  • n_int

Called by:

  • global_selection_buffer

  • pt2_collector()

  • run_pt2_slave_large()

  • run_pt2_slave_small()

  • run_selection_slave()

  • selection_collector()

  • zmq_pt2()

  • zmq_selection()

Calls:

  • check_mem()

delete_selection_buffer:

File : cipsi/selection_buffer.irp.f

subroutine delete_selection_buffer(b)

Called by:

  • global_selection_buffer

  • pt2_collector()

  • run_pt2_slave_large()

  • run_pt2_slave_small()

  • run_selection_slave()

  • selection_collector()

  • zmq_pt2()

  • zmq_selection()

fill_buffer_double:

File : cipsi/selection.irp.f

subroutine fill_buffer_double(i_generator, sp, h1, h2, bannedOrb, banned, fock_diag_tmp, E0, pt2_data, mat, buf)

Needs:

  • det_to_occ_pattern

  • do_only_1h1p

  • h0_type

  • mo_num

  • n_int

  • n_states

  • pseudo_sym

  • psi_det_generators

  • psi_det_hii

  • psi_occ_pattern_hii

  • selection_weight

  • thresh_sym

  • weight_selection

Called by:

  • select_singles_and_doubles()

Calls:

  • add_to_selection_buffer()

  • apply_holes()

  • apply_particles()

  • dsyevd()

fill_buffer_double_rdm:

File : cipsi/pert_rdm_providers.irp.f

subroutine fill_buffer_double_rdm(i_generator, sp, h1, h2, bannedOrb, banned, fock_diag_tmp, E0, pt2_data, mat, buf, psi_det_connection, psi_coef_connection_reverse, n_det_connection)

Needs:

  • det_to_occ_pattern

  • do_only_1h1p

  • h0_type

  • hf_bitmask

  • mo_num

  • n_int

  • n_orb_pert_rdm

  • n_states

  • pert_2rdm_lock

  • pert_2rdm_provider

  • psi_det_generators

  • psi_det_hii

  • psi_occ_pattern_hii

  • selection_weight

  • weight_selection

Called by:

  • select_singles_and_doubles()

Calls:

  • add_to_selection_buffer()

  • apply_holes()

  • apply_particles()

  • get_excitation_degree()

  • give_2rdm_pert_contrib()

  • update_keys_values()

get_d0:

File : cipsi/selection.irp.f

subroutine get_d0(gen, phasemask, bannedOrb, banned, mat, mask, h, p, sp, coefs)

Needs:

  • mo_integrals_map

  • mo_num

  • n_int

  • n_states

Called by:

  • splash_pq()

Calls:

  • apply_particles()

  • get_mo_two_e_integrals()

  • i_h_j()

get_d0_reference:

File : cipsi/selection.irp.f

subroutine get_d0_reference(gen, phasemask, bannedOrb, banned, mat, mask, h, p, sp, coefs)

Needs:

  • mo_num

  • n_int

  • n_states

Calls:

  • apply_particles()

  • i_h_j()

get_d1:

File : cipsi/selection.irp.f

subroutine get_d1(gen, phasemask, bannedOrb, banned, mat, mask, h, p, sp, coefs)

Needs:

  • mo_integrals_map

  • mo_num

  • n_int

  • n_states

Called by:

  • splash_pq()

Calls:

  • apply_particles()

  • get_mo_two_e_integrals()

  • i_h_j()

get_d1_reference:

File : cipsi/selection.irp.f

subroutine get_d1_reference(gen, phasemask, bannedOrb, banned, mat, mask, h, p, sp, coefs)

Needs:

  • mo_num

  • n_int

  • n_states

Calls:

  • apply_particles()

  • i_h_j()

get_d2:

File : cipsi/selection.irp.f

subroutine get_d2(gen, phasemask, bannedOrb, banned, mat, mask, h, p, sp, coefs)

Needs:

  • mo_num

  • n_int

  • n_states

Called by:

  • splash_pq()

get_d2_reference:

File : cipsi/selection.irp.f

subroutine get_d2_reference(gen, phasemask, bannedOrb, banned, mat, mask, h, p, sp, coefs)

Needs:

  • mo_num

  • n_int

  • n_states

get_mask_phase:

File : cipsi/selection.irp.f

subroutine get_mask_phase(det1, pm, Nint)

Called by:

  • splash_pq()

get_phase_bi:

File : cipsi/selection.irp.f

double precision function get_phase_bi(phasemask, s1, s2, h1, p1, h2, p2, Nint)
give_2rdm_pert_contrib:

File : cipsi/update_2rdm.irp.f

subroutine give_2rdm_pert_contrib(det,coef,psi_det_connection,psi_coef_connection_reverse,n_det_connection,nkeys,keys,values,sze_buff)

Needs:

  • elec_alpha_num

  • n_int

  • n_orb_pert_rdm

  • n_states

  • pert_2rdm_lock

  • pert_2rdm_provider

  • state_average_weight

Called by:

  • fill_buffer_double_rdm()

Calls:

  • get_excitation()

  • update_buffer_single_exc_rdm()

  • update_keys_values()

make_selection_buffer_s2:

File : cipsi/selection_buffer.irp.f

subroutine make_selection_buffer_s2(b)

Needs:

  • elec_alpha_num

  • n_int

Called by:

  • zmq_pt2()

  • zmq_selection()

Calls:

  • check_mem()

  • dsort()

  • i8sort()

  • occ_pattern_to_dets()

  • occ_pattern_to_dets_size()

merge_selection_buffers:

File : cipsi/selection_buffer.irp.f

subroutine merge_selection_buffers(b1, b2)

Merges the selection buffers b1 and b2 into b2

Needs:

  • n_int

Called by:

  • run_pt2_slave_large()

Calls:

  • check_mem()

past_d1:

File : cipsi/selection.irp.f

subroutine past_d1(bannedOrb, p)

Needs:

  • mo_num

Called by:

  • splash_pq()

past_d2:

File : cipsi/selection.irp.f

subroutine past_d2(banned, p, sp)

Needs:

  • mo_num

Called by:

  • splash_pq()

provide_everything:

File : cipsi/slave_cipsi.irp.f

subroutine provide_everything

Needs:

  • ci_energy

  • generators_bitmask

  • h_apply_buffer_allocated

  • mo_num

  • mo_two_e_integrals_in_map

  • mpi_master

  • n_det

  • n_det_generators

  • n_det_selectors

  • n_int

  • n_states

  • n_states_diag

  • pseudo_sym

  • psi_coef

  • psi_det_generators

  • psi_det

  • psi_det_generators

  • psi_det_sorted_bit

  • psi_selectors

  • pt2_e0_denominator

  • pt2_stoch_istate

  • selection_weight

  • state_average_weight

  • threshold_generators

  • zmq_context

  • zmq_state

Called by:

  • run_slave_cipsi()

pt2_add:

File : cipsi/pt2_type.irp.f

subroutine pt2_add(p1, w, p2)

p1 =! p1 +( w * p2)

Called by:

  • pt2_collector()

  • selection_collector()

pt2_add2:

File : cipsi/pt2_type.irp.f

subroutine pt2_add2(p1, w, p2)

p1 =! p1 +( w * p2**2)

Called by:

  • pt2_collector()

pt2_addr:

File : cipsi/pt2_type.irp.f

subroutine pt2_addr(p1, a, b, p2)

p1 =! p1 +( a / b * p2)

pt2_alloc:

File : cipsi/pt2_type.irp.f

subroutine pt2_alloc(pt2_data,N)

Called by:

  • pt2_collector()

  • run_cipsi()

  • run_pt2_slave_large()

  • run_pt2_slave_small()

  • run_selection_slave()

  • run_stochastic_cipsi()

  • selection_collector()

pt2_collector:

File : cipsi/pt2_stoch_routines.irp.f

subroutine pt2_collector(zmq_socket_pull, E, relative_error, pt2_data, pt2_data_err, b, N_)

Needs:

  • n_det_generators

  • n_states

  • pt2_f

  • pt2_j

  • pt2_n_teeth

  • pt2_stoch_istate

  • pt2_u

  • pt2_w

Called by:

  • zmq_pt2()

Calls:

  • add_to_selection_buffer()

  • check_mem()

  • create_selection_buffer()

  • delete_selection_buffer()

  • end_zmq_to_qp_run_socket()

  • pt2_add()

  • pt2_add2()

  • pt2_alloc()

  • pt2_dealloc()

  • pull_pt2_results()

  • sleep()

  • sort_selection_buffer()

  • wall_time()

pt2_dealloc:

File : cipsi/pt2_type.irp.f

subroutine pt2_dealloc(pt2_data)

Called by:

  • pt2_collector()

  • run_cipsi()

  • run_pt2_slave_large()

  • run_pt2_slave_small()

  • run_selection_slave()

  • run_stochastic_cipsi()

  • selection_collector()

pt2_deserialize:

File : cipsi/pt2_type.irp.f

subroutine pt2_deserialize(pt2_data, n, x)

Called by:

  • pull_pt2_results()

  • pull_selection_results()

pt2_find_sample:

File : cipsi/pt2_stoch_routines.irp.f

integer function pt2_find_sample(v, w)

Needs:

  • n_det_generators

pt2_find_sample_lr:

File : cipsi/pt2_stoch_routines.irp.f

integer function pt2_find_sample_lr(v, w, l_in, r_in)

Needs:

  • n_det_generators

pt2_serialize:

File : cipsi/pt2_type.irp.f

subroutine pt2_serialize(pt2_data, n, x)

Called by:

  • push_pt2_results_async_send()

  • push_selection_results()

pt2_slave_inproc:

File : cipsi/pt2_stoch_routines.irp.f

subroutine pt2_slave_inproc(i)

Needs:

  • global_selection_buffer

  • pt2_e0_denominator

Called by:

  • zmq_pt2()

Calls:

  • run_pt2_slave()

pull_pt2_results:

File : cipsi/run_pt2_slave.irp.f

subroutine pull_pt2_results(zmq_socket_pull, index, pt2_data, task_id, n_tasks, b)

Needs:

  • n_int

  • n_states

Called by:

  • pt2_collector()

Calls:

  • pt2_deserialize()

pull_selection_results:

File : cipsi/run_selection_slave.irp.f

subroutine pull_selection_results(zmq_socket_pull, pt2_data, val, det, N, task_id, ntasks)

Needs:

  • n_int

  • n_states

Called by:

  • selection_collector()

Calls:

  • pt2_deserialize()

push_pt2_results:

File : cipsi/run_pt2_slave.irp.f

subroutine push_pt2_results(zmq_socket_push, index, pt2_data, b, task_id, n_tasks)

Called by:

  • run_pt2_slave_small()

Calls:

  • push_pt2_results_async_recv()

  • push_pt2_results_async_send()

push_pt2_results_async_recv:

File : cipsi/run_pt2_slave.irp.f

subroutine push_pt2_results_async_recv(zmq_socket_push,mini,sending)

Called by:

  • push_pt2_results()

  • run_pt2_slave_large()

push_pt2_results_async_send:

File : cipsi/run_pt2_slave.irp.f

subroutine push_pt2_results_async_send(zmq_socket_push, index, pt2_data, b, task_id, n_tasks, sending)

Needs:

  • n_int

  • n_states

Called by:

  • push_pt2_results()

  • run_pt2_slave_large()

Calls:

  • pt2_serialize()

push_selection_results:

File : cipsi/run_selection_slave.irp.f

subroutine push_selection_results(zmq_socket_push, pt2_data, b, task_id, ntasks)

Needs:

  • n_int

  • n_states

Called by:

  • run_selection_slave()

Calls:

  • pt2_serialize()

remove_duplicates_in_selection_buffer:

File : cipsi/selection_buffer.irp.f

subroutine remove_duplicates_in_selection_buffer(b)

Needs:

  • n_int

Called by:

  • zmq_pt2()

Calls:

  • check_mem()

  • i8sort()

run_cipsi:

File : cipsi/cipsi.irp.f

subroutine run_cipsi

Selected Full Configuration Interaction with deterministic selection and stochastic PT2.

Needs:

  • correlation_energy_ratio_max

  • do_pt2

  • h_apply_buffer_allocated

  • n_det

  • n_det_max

  • n_iter

  • n_states

  • n_states_diag

  • psi_coef

  • psi_det

  • psi_det_sorted

  • psi_energy

  • psi_energy_with_nucl_rep

  • psi_occ_pattern

  • pt2_max

  • pt2_relative_error

  • ref_bitmask_energy

  • s2_eig

  • selection_factor

  • threshold_generators

  • variance_max

Called by:

  • fci()

Calls:

  • check_mem()

  • copy_h_apply_buffer_to_wf()

  • diagonalize_ci()

  • ezfio_get_hartree_fock_energy()

  • ezfio_has_hartree_fock_energy()

  • make_s2_eigenfunction()

  • print_extrapolated_energy()

  • print_summary()

  • pt2_alloc()

  • pt2_dealloc()

  • save_energy()

  • save_iterations()

  • save_wavefunction()

  • write_double()

  • zmq_pt2()

  • zmq_selection()

Touches:

  • ci_electronic_energy

  • ci_electronic_energy

  • ci_energy

  • ci_electronic_energy

  • n_det

  • n_iter

  • psi_occ_pattern

  • c0_weight

  • psi_coef

  • psi_det_sorted_bit

  • psi_det

  • psi_det_size

  • psi_det_sorted_bit

  • psi_energy

  • psi_occ_pattern

  • psi_energy

  • pt2_match_weight

  • pt2_overlap

  • pt2_stoch_istate

  • selection_weight

  • state_average_weight

  • threshold_davidson_pt2

  • threshold_generators

  • variance_match_weight

run_pt2_slave:

File : cipsi/run_pt2_slave.irp.f

subroutine run_pt2_slave(thread,iproc,energy)

Needs:

  • n_states_diag

Called by:

  • pt2_slave_inproc()

  • run_slave_main()

Calls:

  • run_pt2_slave_large()

run_pt2_slave_large:

File : cipsi/run_pt2_slave.irp.f

subroutine run_pt2_slave_large(thread,iproc,energy)

Needs:

  • elec_alpha_num

  • global_selection_buffer

  • global_selection_buffer_lock

  • mo_num

  • n_states

  • n_states_diag

  • pt2_f

Called by:

  • run_pt2_slave()

Calls:

  • create_selection_buffer()

  • delete_selection_buffer()

  • end_zmq_push_socket()

  • end_zmq_to_qp_run_socket()

  • merge_selection_buffers()

  • omp_set_lock()

  • omp_unset_lock()

  • pt2_alloc()

  • pt2_dealloc()

  • push_pt2_results_async_recv()

  • push_pt2_results_async_send()

  • select_connected()

  • sleep()

  • sort_selection_buffer()

  • wall_time()

run_pt2_slave_small:

File : cipsi/run_pt2_slave.irp.f

subroutine run_pt2_slave_small(thread,iproc,energy)

Needs:

  • elec_alpha_num

  • mo_num

  • n_states

  • n_states_diag

  • nproc

  • pt2_f

Calls:

  • create_selection_buffer()

  • delete_selection_buffer()

  • end_zmq_push_socket()

  • end_zmq_to_qp_run_socket()

  • pt2_alloc()

  • pt2_dealloc()

  • push_pt2_results()

  • select_connected()

  • sort_selection_buffer()

  • usleep()

  • wall_time()

run_selection_slave:

File : cipsi/run_selection_slave.irp.f

subroutine run_selection_slave(thread,iproc,energy)

Needs:

  • n_int

  • n_states

  • pseudo_sym

  • psi_bilinear_matrix_columns_loc

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_transp_rows_loc

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • psi_det_sorted

  • psi_det_sorted

  • psi_selectors_coef_transp

  • pt2_f

  • weight_selection

Called by:

  • run_slave_main()

  • selection_slave_inproc()

Calls:

  • create_selection_buffer()

  • delete_selection_buffer()

  • end_zmq_push_socket()

  • end_zmq_to_qp_run_socket()

  • pt2_alloc()

  • pt2_dealloc()

  • push_selection_results()

  • select_connected()

  • sort_selection_buffer()

  • usleep()

run_slave_cipsi:

File : cipsi/slave_cipsi.irp.f

subroutine run_slave_cipsi

Helper program for distributed parallelism

Needs:

  • distributed_davidson

  • read_wf

Called by:

  • fci()

  • pt2()

Calls:

  • omp_set_nested()

  • provide_everything()

  • run_slave_main()

  • switch_qp_run_to_master()

Touches:

  • distributed_davidson

  • pt2_e0_denominator

  • pt2_stoch_istate

  • read_wf

  • selection_weight

  • state_average_weight

  • threshold_generators

run_slave_main:

File : cipsi/slave_cipsi.irp.f

subroutine run_slave_main

Needs:

  • det_to_occ_pattern

  • elec_alpha_num

  • global_selection_buffer

  • h0_type

  • mo_num

  • mpi_master

  • mpi_rank

  • n_det

  • n_det_generators

  • n_det_selectors

  • n_int

  • n_states

  • n_states_diag

  • nthreads_pt2

  • psi_coef

  • psi_det

  • pt2_e0_denominator

  • pt2_f

  • pt2_stoch_istate

  • qp_max_mem

  • selection_weight

  • state_average_weight

  • threshold_generators

  • zmq_context

  • zmq_state

Called by:

  • run_slave_cipsi()

Calls:

  • check_mem()

  • davidson_slave_tcp()

  • mpi_print()

  • omp_set_nested()

  • resident_memory()

  • run_pt2_slave()

  • run_selection_slave()

  • usleep()

  • wait_for_states()

  • wall_time()

  • write_double()

Touches:

  • pt2_e0_denominator

  • pt2_stoch_istate

  • selection_weight

  • state_average_weight

  • threshold_generators

run_stochastic_cipsi:

File : cipsi/stochastic_cipsi.irp.f

subroutine run_stochastic_cipsi

Selected Full Configuration Interaction with Stochastic selection and PT2.

Needs:

  • correlation_energy_ratio_max

  • h_apply_buffer_allocated

  • n_det

  • n_det_max

  • n_iter

  • n_states

  • n_states_diag

  • psi_coef

  • psi_det

  • psi_det_sorted

  • psi_energy

  • psi_energy_with_nucl_rep

  • psi_occ_pattern

  • pt2_max

  • pt2_relative_error

  • ref_bitmask_energy

  • s2_eig

  • selection_factor

  • threshold_generators

  • variance_max

Called by:

  • fci()

Calls:

  • check_mem()

  • copy_h_apply_buffer_to_wf()

  • diagonalize_ci()

  • ezfio_get_hartree_fock_energy()

  • ezfio_has_hartree_fock_energy()

  • make_s2_eigenfunction()

  • print_extrapolated_energy()

  • print_summary()

  • pt2_alloc()

  • pt2_dealloc()

  • save_energy()

  • save_iterations()

  • save_wavefunction()

  • write_double()

  • zmq_pt2()

Touches:

  • ci_electronic_energy

  • ci_electronic_energy

  • ci_energy

  • ci_electronic_energy

  • n_det

  • n_iter

  • psi_occ_pattern

  • c0_weight

  • psi_coef

  • psi_det_sorted_bit

  • psi_det

  • psi_det_size

  • psi_det_sorted_bit

  • psi_energy

  • psi_occ_pattern

  • psi_energy

  • pt2_match_weight

  • pt2_overlap

  • pt2_stoch_istate

  • selection_weight

  • state_average_weight

  • threshold_davidson_pt2

  • threshold_generators

  • variance_match_weight

select_connected:

File : cipsi/selection.irp.f

subroutine select_connected(i_generator,E0,pt2_data,b,subset,csubset)

Needs:

  • generators_bitmask

  • mo_num

  • n_int

  • n_states

  • psi_det_generators

Called by:

  • run_pt2_slave_large()

  • run_pt2_slave_small()

  • run_selection_slave()

Calls:

  • build_fock_tmp()

  • select_singles_and_doubles()

select_singles_and_doubles:

File : cipsi/selection.irp.f

subroutine select_singles_and_doubles(i_generator,hole_mask,particle_mask,fock_diag_tmp,E0,pt2_data,buf,subset,csubset)

WARNING /!: It is assumed that the generators and selectors are psi_det_sorted

Needs:

  • banned_excitation

  • mo_num

  • n_det

  • n_det_selectors

  • n_int

  • n_states

  • pert_2rdm

  • psi_bilinear_matrix_columns_loc

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_transp_rows_loc

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_values

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • psi_det_generators

  • psi_det_sorted

  • psi_det_sorted

  • psi_selectors_coef_transp

Called by:

  • select_connected()

Calls:

  • apply_hole()

  • bitstring_to_list_ab()

  • fill_buffer_double()

  • fill_buffer_double_rdm()

  • get_excitation_degree_spin()

  • isort()

  • splash_pq()

  • spot_isinwf()

selection_collector:

File : cipsi/zmq_selection.irp.f

subroutine selection_collector(zmq_socket_pull, b, N, pt2_data)

Needs:

  • n_det_generators

  • n_states

Called by:

  • zmq_selection()

Calls:

  • add_to_selection_buffer()

  • check_mem()

  • create_selection_buffer()

  • delete_selection_buffer()

  • end_zmq_to_qp_run_socket()

  • pt2_add()

  • pt2_alloc()

  • pt2_dealloc()

  • pull_selection_results()

  • sort_selection_buffer()

selection_slave_inproc:

File : cipsi/zmq_selection.irp.f

subroutine selection_slave_inproc(i)

Needs:

  • pt2_e0_denominator

Called by:

  • zmq_selection()

Calls:

  • run_selection_slave()

sort_selection_buffer:

File : cipsi/selection_buffer.irp.f

subroutine sort_selection_buffer(b)

Needs:

  • n_int

Called by:

  • add_to_selection_buffer()

  • pt2_collector()

  • run_pt2_slave_large()

  • run_pt2_slave_small()

  • run_selection_slave()

  • selection_collector()

Calls:

  • check_mem()

  • dsort()

splash_pq:

File : cipsi/selection.irp.f

subroutine splash_pq(mask, sp, det, i_gen, N_sel, bannedOrb, banned, mat, interesting)

Computes the contributions A(r,s) by comparing the external determinant to all the internal determinants det(i). an applying two particles (r,s) to the mask.

Needs:

  • mo_num

  • n_int

  • n_states

  • psi_det_sorted

  • psi_selectors_coef_transp

Called by:

  • select_singles_and_doubles()

Calls:

  • bitstring_to_list_in_selection()

  • get_d0()

  • get_d1()

  • get_d2()

  • get_mask_phase()

  • past_d1()

  • past_d2()

spot_isinwf:

File : cipsi/selection.irp.f

subroutine spot_isinwf(mask, det, i_gen, N, banned, fullMatch, interesting)

Identify the determinants in det which are in the internal space. These are the determinants that can be produced by creating two particles on the mask.

Needs:

  • mo_num

  • n_int

Called by:

  • select_singles_and_doubles()

Calls:

  • bitstring_to_list_in_selection()

testteethbuilding:

File : cipsi/pt2_stoch_routines.irp.f

logical function testTeethBuilding(minF, N)

Needs:

  • n_det_generators

  • psi_det_sorted_gen

  • pt2_stoch_istate

Calls:

  • check_mem()

update_buffer_double_exc_rdm:

File : cipsi/update_2rdm.irp.f

subroutine update_buffer_double_exc_rdm(exc,phase,contrib,nkeys,keys,values,sze_buff)

Needs:

  • list_orb_reverse_pert_rdm

update_buffer_single_exc_rdm:

File : cipsi/update_2rdm.irp.f

subroutine update_buffer_single_exc_rdm(det1,det2,exc,phase,contrib,nkeys,keys,values,sze_buff)

Needs:

  • list_orb_reverse_pert_rdm

  • n_int

Called by:

  • give_2rdm_pert_contrib()

Calls:

  • bitstring_to_list_ab()

update_pt2_and_variance_weights:

File : cipsi/selection.irp.f

subroutine update_pt2_and_variance_weights(pt2_data, N_st)

Updates the PT2- and Variance- matching weights.

Needs:

  • n_det

  • n_states

  • pt2_match_weight

  • pt2_relative_error

  • threshold_davidson

  • threshold_davidson_pt2

  • variance_match_weight

Called by:

  • zmq_pt2()

  • zmq_selection()

Touches:

  • pt2_match_weight

  • threshold_davidson_pt2

  • variance_match_weight

zmq_pt2:

File : cipsi/pt2_stoch_routines.irp.f

subroutine ZMQ_pt2(E, pt2_data, pt2_data_err, relative_error, N_in)

Needs:

  • det_to_occ_pattern

  • elec_alpha_num

  • global_selection_buffer

  • h0_type

  • mo_num

  • mo_one_e_integrals

  • mo_two_e_integrals_in_map

  • n_det

  • n_det_generators

  • n_int

  • n_states

  • nproc

  • nthreads_pt2

  • pseudo_sym

  • psi_bilinear_matrix_columns_loc

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_transp_rows_loc

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • psi_det_hii

  • psi_det_sorted

  • psi_det_sorted

  • psi_occ_pattern_hii

  • psi_selectors

  • psi_selectors_coef_transp

  • pt2_e0_denominator

  • pt2_f

  • pt2_j

  • pt2_n_teeth

  • pt2_overlap

  • pt2_j

  • pt2_stoch_istate

  • pt2_u

  • pt2_w

  • qp_max_mem

  • s2_eig

  • selection_weight

  • state_average_weight

  • threshold_generators

Called by:

  • run_cipsi()

  • run_stochastic_cipsi()

Calls:

  • check_mem()

  • create_selection_buffer()

  • delete_selection_buffer()

  • end_parallel_job()

  • fill_h_apply_buffer_no_selection()

  • make_selection_buffer_s2()

  • new_parallel_job()

  • omp_set_nested()

  • pt2_collector()

  • pt2_slave_inproc()

  • remove_duplicates_in_selection_buffer()

  • resident_memory()

  • update_pt2_and_variance_weights()

  • write_double()

  • write_int()

  • zmq_selection()

Touches:

  • pt2_match_weight

  • pt2_overlap

  • pt2_stoch_istate

  • selection_weight

  • state_average_weight

  • threshold_davidson_pt2

  • variance_match_weight

zmq_selection:

File : cipsi/zmq_selection.irp.f

subroutine ZMQ_selection(N_in, pt2_data)

Needs:

  • do_pt2

  • elec_alpha_num

  • mo_num

  • n_det

  • n_det_generators

  • n_det_selectors

  • n_int

  • n_states

  • nproc

  • pseudo_sym

  • psi_bilinear_matrix_columns_loc

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_values

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_transp_values

  • psi_bilinear_matrix_transp_rows_loc

  • psi_det_alpha_unique

  • psi_det_beta_unique

  • psi_det_sorted

  • psi_selectors

  • pt2_e0_denominator

  • pt2_f

  • pt2_overlap

  • qp_max_mem

  • s2_eig

  • selection_weight

  • state_average_weight

  • threshold_generators

Called by:

  • run_cipsi()

  • zmq_pt2()

Calls:

  • create_selection_buffer()

  • delete_selection_buffer()

  • end_parallel_job()

  • fill_h_apply_buffer_no_selection()

  • make_selection_buffer_s2()

  • new_parallel_job()

  • selection_collector()

  • selection_slave_inproc()

  • update_pt2_and_variance_weights()

  • write_double()

Touches:

  • pt2_match_weight

  • pt2_overlap

  • threshold_davidson_pt2

  • variance_match_weight