LCPQ/quantum_package
10
0
Fork 0
mirror of https://github.com/LCPQ/quantum_package synced 2024-11-03 20:54:00 +01:00
Code Issues Releases Wiki Activity

Merge pull request #204 from scemama/master

Browse Source 
Fixes ddcc1a5
This commit is contained in:
Thomas Applencourt 2017-06-26 15:06:54 -05:00 committed by GitHub
commit 5b70d171fe
39 changed files with 149 additions and 1698 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

10
plugins/CASSCF/EZFIO.cfg
View File

@ -1,10 +0,0 @@
[energy]
type: double precision
doc: "Calculated CAS-SCF energy"
interface: ezfio
[energy_pt2]
type: double precision
doc: "Calculated selected CAS-SCF energy with PT2 correction"
interface: ezfio

39
plugins/CASSCF/H_apply.irp.f
View File

@ -1,39 +0,0 @@
use bitmasks
BEGIN_SHELL [ /usr/bin/env python ]
from generate_h_apply import *
s = H_apply("CAS_SD")
print s
s = H_apply("CAS_SD_selected_no_skip")
s.set_selection_pt2("epstein_nesbet_2x2")
s.unset_skip()
print s
s = H_apply("CAS_SD_selected")
s.set_selection_pt2("epstein_nesbet_2x2")
print s
s = H_apply("CAS_SD_PT2")
s.set_perturbation("epstein_nesbet_2x2")
print s
s = H_apply("CAS_S",do_double_exc=False)
print s
s = H_apply("CAS_S_selected_no_skip",do_double_exc=False)
s.set_selection_pt2("epstein_nesbet_2x2")
s.unset_skip()
print s
s = H_apply("CAS_S_selected",do_double_exc=False)
s.set_selection_pt2("epstein_nesbet_2x2")
print s
s = H_apply("CAS_S_PT2",do_double_exc=False)
s.set_perturbation("epstein_nesbet_2x2")
print s
END_SHELL

1
plugins/CASSCF/NEEDED_CHILDREN_MODULES
View File

@ -1 +0,0 @@
Generators_CAS Perturbation Selectors_full

20
plugins/CASSCF/README.rst
View File

@ -1,20 +0,0 @@
======
CASSCF
======
This module is not a "real" CAS-SCF. It is an orbital optimization step done by :
1) Doing the CAS+SD
2) Taking one-electron density matrix
3) Cancelling all active-active rotations
4) Finding the order which matches with the input MOs
Needed Modules
==============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
Documentation
=============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.

217
plugins/CASSCF/casscf.irp.f
View File

@ -1,217 +0,0 @@
program casscf
implicit none
BEGIN_DOC
! Optimize MOs and CI coefficients of the CAS
END_DOC
double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:)
integer(bit_kind), allocatable :: generators_bitmask_save(:,:,:,:)
integer :: degree, N_generators_bitmask_save, N_det_ci
double precision :: E_old, E_CI
double precision :: selection_criterion_save, selection_criterion_min_save
integer :: N_det_old
integer :: i, j, k, l
integer :: i_bit, j_bit, i_int, j_int
integer(bit_kind), allocatable :: bit_tmp(:), cas_bm(:)
character*(64) :: label
allocate( pt2(N_states), norm_pert(N_states),H_pert_diag(N_states) )
allocate( generators_bitmask_save(N_int,2,6,N_generators_bitmask) )
allocate( bit_tmp(N_int), cas_bm(N_int) )
PROVIDE N_det_cas
N_det_old = 0
pt2 = 1.d0
E_CI = 1.d0
E_old = 0.d0
diag_algorithm = "Lapack"
selection_criterion_save = selection_criterion
selection_criterion_min_save = selection_criterion_min
cas_bm = 0_bit_kind
do i=1,N_cas_bitmask
do j=1,N_int
cas_bm(j) = ior(cas_bm(j), cas_bitmask(j,1,i))
cas_bm(j) = ior(cas_bm(j), cas_bitmask(j,2,i))
enddo
enddo
! Save CAS-SD bitmask
generators_bitmask_save = generators_bitmask
N_generators_bitmask_save = N_generators_bitmask
! Set the CAS bitmask
do i=1,6
generators_bitmask(:,:,i,:) = cas_bitmask
enddo
N_generators_bitmask = N_cas_bitmask
SOFT_TOUCH generators_bitmask N_generators_bitmask
! If the number of dets already in the file is larger than the requested
! number of determinants, truncate the wf
if (N_det > N_det_max) then
call diagonalize_CI
call save_wavefunction
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
endif
! Start MCSCF iteration
! CAS-CI
! ------
E_old = E_CI
! Reset the selection criterion
selection_criterion = selection_criterion_save
selection_criterion_min = selection_criterion_min_save
SOFT_TOUCH selection_criterion_min selection_criterion selection_criterion_factor
! Set the CAS bitmask
do i=1,6
generators_bitmask(:,:,i,:) = cas_bitmask
enddo
N_generators_bitmask = N_cas_bitmask
SOFT_TOUCH generators_bitmask N_generators_bitmask
do while (N_det < N_det_max.and.maxval(abs(pt2(1:N_states))) > pt2_max)
N_det_old = N_det
call H_apply_CAS_SD_selected_no_skip(pt2, norm_pert, H_pert_diag, N_states)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
if (N_det > N_det_max) then
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
endif
call diagonalize_CI
call save_wavefunction
print *, '======'
print *, 'CAS-CI'
print *, '======'
print *, ''
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E(CAS) = ', CI_energy
print *, 'E(CAS)+PT2 = ', CI_energy+pt2
print *, '-----'
print *, ''
E_CI = sum(CI_energy(1:N_states)+pt2(1:N_states))/dble(N_states)
call ezfio_set_casscf_energy(CI_energy(1))
if (N_det == N_det_old) then
exit
endif
enddo
! Super-CI
! --------
selection_criterion_min = 1.d-12
selection_criterion = 1.d-12
! Set the CAS bitmask
generators_bitmask = generators_bitmask_save
N_generators_bitmask = N_generators_bitmask_save
SOFT_TOUCH generators_bitmask N_generators_bitmask selection_criterion selection_criterion_min selection_criterion_factor
N_det_ci = N_det
call H_apply_CAS_SD_selected(pt2, norm_pert, H_pert_diag, N_states)
do i=1,mo_tot_num
i_int = ishft(i-1,-bit_kind_shift)+1
i_bit = j-ishft(i_int-1,bit_kind_shift)-1
bit_tmp(:) = 0_bit_kind
bit_tmp(i_int) = ibset(0_bit_kind,i_bit)
if (iand(bit_tmp(i_int), cas_bm(i_int)) == 0_bit_kind) then
! Not a CAS MO
cycle
endif
do j=1,mo_tot_num
if (j == i) then
cycle
endif
j_int = ishft(j-1,-bit_kind_shift)+1
j_bit = j-ishft(j_int-1,bit_kind_shift)-1
bit_tmp(:) = 0_bit_kind
bit_tmp(j_int) = ibset(0_bit_kind,j_bit)
if (iand(bit_tmp(j_int), cas_bm(j_int)) == 0_bit_kind) then
! Not a CAS MO
cycle
endif
! Now, both i and j are MOs of the CAS. De-couple them in the DM
one_body_dm_mo(i,j) = 0.d0
enddo
enddo
SOFT_TOUCH one_body_dm_mo
double precision :: mx, ov
double precision, allocatable :: mo_coef_old(:,:)
integer, allocatable :: iorder(:)
logical, allocatable :: selected(:)
allocate( mo_coef_old(size(mo_coef,1), size(mo_coef,2)), iorder(mo_tot_num), selected(mo_tot_num) )
mo_coef_old = mo_coef
label = "Canonical"
call mo_as_eigvectors_of_mo_matrix(one_body_dm_mo,size(one_body_dm_mo,1),size(one_body_dm_mo,2),label,-1)
selected = .False.
do j=1,mo_tot_num
mx = -1.d0
iorder(j) = j
do i=1,mo_tot_num
if (selected(i)) then
cycle
endif
ov = 0.d0
do l=1,ao_num
do k=1,ao_num
ov = ov + mo_coef_old(k,j) * ao_overlap(k,l) * mo_coef(l,i)
enddo
enddo
ov= dabs(ov)
if (ov > mx) then
mx = ov
iorder(j) = i
endif
enddo
selected( iorder(j) ) = .True.
enddo
mo_coef_old = mo_coef
do i=1,mo_tot_num
mo_coef(:,i) = mo_coef_old(:,iorder(i))
enddo
call save_mos
call write_double(6,E_CI,"Energy(CAS)")
deallocate( mo_coef_old )
deallocate( pt2, norm_pert,H_pert_diag )
deallocate( generators_bitmask_save )
deallocate( bit_tmp, cas_bm, iorder )
end

10
plugins/CAS_SD/EZFIO.cfg
View File

@ -1,10 +0,0 @@
[energy]
type: double precision
doc: "Calculated CAS-SD energy"
interface: ezfio
[energy_pt2]
type: double precision
doc: "Calculated selected CAS-SD energy with PT2 correction"
interface: ezfio

37
plugins/CAS_SD/H_apply.irp.f
View File

@ -1,37 +0,0 @@
use bitmasks
BEGIN_SHELL [ /usr/bin/env python ]
from generate_h_apply import *
s = H_apply("CAS_SD")
s.unset_skip()
print s
s = H_apply("CAS_SD_selected_no_skip")
s.set_selection_pt2("epstein_nesbet_2x2")
s.unset_skip()
print s
s = H_apply("CAS_SD_selected")
s.set_selection_pt2("epstein_nesbet_2x2")
s.unset_skip()
print s
s = H_apply("CAS_SD_PT2")
s.set_perturbation("epstein_nesbet_2x2")
print s
s = H_apply("CAS_S",do_double_exc=False)
print s
s = H_apply("CAS_S_selected",do_double_exc=False)
s.set_selection_pt2("epstein_nesbet_2x2")
s.unset_skip()
print s
s = H_apply("CAS_S_PT2",do_double_exc=False)
s.set_perturbation("epstein_nesbet_2x2")
print s
END_SHELL

1
plugins/CAS_SD/NEEDED_CHILDREN_MODULES
View File

@ -1 +0,0 @@
Perturbation Selectors_full Generators_CAS Davidson

295
plugins/CAS_SD/README.rst
View File

@ -1,295 +0,0 @@
======================
CAS_SD_selected Module
======================
Selected CAS + SD module.
1) Set the different MO classes using the ``qp_set_mo_class`` command
2) Run the selected CAS+SD program
Documentation
=============
.. Do not edit this section. It was auto-generated from the
.. by the `update_README.py` script.
`full_ci <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/cas_sd_selected.irp.f#L1>`_
Undocumented
`h_apply_cas_sd <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/H_apply.irp.f_shell_22#L414>`_
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cas_sd_diexc <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/H_apply.irp.f_shell_22#L1>`_
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_monoexc <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/H_apply.irp.f_shell_22#L269>`_
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_pt2 <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/H_apply.irp.f_shell_22#L2610>`_
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cas_sd_pt2_diexc <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/H_apply.irp.f_shell_22#L2118>`_
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_pt2_monoexc <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/H_apply.irp.f_shell_22#L2427>`_
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_selected <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/H_apply.irp.f_shell_22#L1872>`_
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cas_sd_selected_diexc <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/H_apply.irp.f_shell_22#L1346>`_
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_selected_monoexc <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/H_apply.irp.f_shell_22#L1675>`_
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_selected_no_skip <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/H_apply.irp.f_shell_22#L1128>`_
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_cas_sd_selected_no_skip_diexc <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/H_apply.irp.f_shell_22#L602>`_
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_cas_sd_selected_no_skip_monoexc <http://github.com/LCPQ/quantum_package/tree/master/src/CAS_SD/H_apply.irp.f_shell_22#L931>`_
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
Needed Modules
==============
.. Do not edit this section. It was auto-generated from the
.. by the `update_README.py` script.
.. image:: tree_dependency.png
* `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/src/Perturbation>`_
* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
* `Generators_CAS <http://github.com/LCPQ/quantum_package/tree/master/src/Generators_CAS>`_
Needed Modules
==============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
.. image:: tree_dependency.png
* `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation>`_
* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/plugins/Selectors_full>`_
* `Generators_CAS <http://github.com/LCPQ/quantum_package/tree/master/plugins/Generators_CAS>`_
* `Davidson <http://github.com/LCPQ/quantum_package/tree/master/src/Davidson>`_
Documentation
=============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
`full_ci <http://github.com/LCPQ/quantum_package/tree/master/plugins/CAS_SD/cas_sd_selected.irp.f#L1>`_
Undocumented
h_apply_cas_s
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_cas_s_diexc
Undocumented
h_apply_cas_s_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_s_diexcp
Undocumented
h_apply_cas_s_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_s_pt2
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_cas_s_pt2_diexc
Undocumented
h_apply_cas_s_pt2_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_s_pt2_diexcp
Undocumented
h_apply_cas_s_pt2_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_s_selected
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_cas_s_selected_diexc
Undocumented
h_apply_cas_s_selected_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_s_selected_diexcp
Undocumented
h_apply_cas_s_selected_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_sd
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_cas_sd_diexc
Undocumented
h_apply_cas_sd_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_sd_diexcp
Undocumented
h_apply_cas_sd_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_sd_pt2
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_cas_sd_pt2_diexc
Undocumented
h_apply_cas_sd_pt2_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_sd_pt2_diexcp
Undocumented
h_apply_cas_sd_pt2_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_sd_selected
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_cas_sd_selected_diexc
Undocumented
h_apply_cas_sd_selected_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_sd_selected_diexcp
Undocumented
h_apply_cas_sd_selected_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_sd_selected_no_skip
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_cas_sd_selected_no_skip_diexc
Undocumented
h_apply_cas_sd_selected_no_skip_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_cas_sd_selected_no_skip_diexcp
Undocumented
h_apply_cas_sd_selected_no_skip_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.

92
plugins/CAS_SD/cas_s.irp.f
View File

@ -1,92 +0,0 @@
program full_ci
implicit none
integer :: i,k
integer :: N_det_old
double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:)
integer :: N_st, degree
N_st = N_states
allocate (pt2(N_st), norm_pert(N_st),H_pert_diag(N_st))
character*(64) :: perturbation
PROVIDE N_det_cas
N_det_old = 0
pt2 = 1.d0
diag_algorithm = "Lapack"
if (N_det > N_det_max) then
call diagonalize_CI
call save_wavefunction
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
endif
do while (N_det < N_det_max.and.maxval(abs(pt2(1:N_st))) > pt2_max)
N_det_old = N_det
call H_apply_CAS_S(pt2, norm_pert, H_pert_diag, N_st)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
if (N_det > N_det_max) then
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
endif
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
call ezfio_set_cas_sd_energy(CI_energy(1))
if (N_det == N_det_old) then
exit
endif
enddo
call diagonalize_CI
if(do_pt2_end)then
print*,'Last iteration only to compute the PT2'
threshold_selectors = 1.d0
threshold_generators = 0.999d0
call H_apply_CAS_S_PT2(pt2, norm_pert, H_pert_diag, N_st)
print *, 'Final step'
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
call ezfio_set_cas_sd_energy_pt2(CI_energy(1)+pt2(1))
endif
integer :: exc_max, degree_min
exc_max = 0
print *, 'CAS determinants : ', N_det_cas
do i=1,min(N_det_cas,10)
do k=i,N_det_cas
call get_excitation_degree(psi_cas(1,1,k),psi_cas(1,1,i),degree,N_int)
exc_max = max(exc_max,degree)
enddo
call debug_det(psi_cas(1,1,i),N_int)
print *, ''
enddo
print *, 'Max excitation degree in the CAS :', exc_max
end

123
plugins/CAS_SD/cas_s_selected.irp.f
View File

@ -1,123 +0,0 @@
program full_ci
implicit none
integer :: i,k
double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:)
integer :: N_st, degree
N_st = N_states
allocate (pt2(N_st), norm_pert(N_st),H_pert_diag(N_st))
character*(64) :: perturbation
PROVIDE N_det_cas
pt2 = 1.d0
diag_algorithm = "Lapack"
if (N_det > N_det_max) then
call diagonalize_CI
call save_wavefunction
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
endif
double precision :: i_H_psi_array(N_states),diag_H_mat_elem,h,i_O1_psi_array(N_states)
double precision :: E_CI_before(N_states)
if(read_wf)then
call i_H_psi(psi_det(1,1,N_det),psi_det,psi_coef,N_int,N_det,psi_det_size,N_states,i_H_psi_array)
h = diag_H_mat_elem(psi_det(1,1,N_det),N_int)
selection_criterion = dabs(psi_coef(N_det,1) * (i_H_psi_array(1) - h * psi_coef(N_det,1))) * 0.1d0
soft_touch selection_criterion
endif
integer :: n_det_before
print*,'Beginning the selection ...'
E_CI_before(1:N_states) = CI_energy(1:N_states)
do while (N_det < N_det_max.and.maxval(abs(pt2(1:N_st))) > pt2_max)
n_det_before = N_det
call H_apply_CAS_SD_selected(pt2, norm_pert, H_pert_diag, N_st)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
call diagonalize_CI
if (N_det > N_det_max) then
N_det = N_det_max
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
touch N_det psi_det psi_coef psi_det_sorted psi_coef_sorted psi_average_norm_contrib_sorted
endif
call save_wavefunction
if(n_det_before == N_det)then
selection_criterion = selection_criterion * 0.5d0
endif
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
do k = 1, N_states
print*,'State ',k
print *, 'PT2 = ', pt2(k)
print *, 'E = ', CI_energy(k)
print *, 'E(before)+PT2 = ', E_CI_before(k)+pt2(k)
enddo
print *, '-----'
if(N_states.gt.1)then
print*,'Variational Energy difference'
do i = 2, N_states
print*,'Delta E = ',CI_energy(i) - CI_energy(1)
enddo
endif
if(N_states.gt.1)then
print*,'Variational + perturbative Energy difference'
do i = 2, N_states
print*,'Delta E = ',E_CI_before(i)+ pt2(i) - (E_CI_before(1) + pt2(1))
enddo
endif
E_CI_before(1:N_states) = CI_energy(1:N_states)
call ezfio_set_cas_sd_energy(CI_energy(1))
enddo
N_det = min(N_det_max,N_det)
touch N_det psi_det psi_coef
call diagonalize_CI
if(do_pt2_end)then
print*,'Last iteration only to compute the PT2'
threshold_selectors = 1.d0
threshold_generators = 0.999d0
call H_apply_CAS_SD_PT2(pt2, norm_pert, H_pert_diag, N_st)
print *, 'Final step'
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy(1:N_states)
print *, 'E+PT2 = ', CI_energy(1:N_states)+pt2(1:N_states)
print *, '-----'
call ezfio_set_cas_sd_energy_pt2(CI_energy(1)+pt2(1))
endif
integer :: exc_max, degree_min
exc_max = 0
print *, 'CAS determinants : ', N_det_cas
do i=1,min(N_det_cas,10)
do k=i,N_det_cas
call get_excitation_degree(psi_cas(1,1,k),psi_cas(1,1,i),degree,N_int)
exc_max = max(exc_max,degree)
enddo
print *, psi_coef_cas_diagonalized(i,:)
call debug_det(psi_cas(1,1,i),N_int)
print *, ''
enddo
print *, 'Max excitation degree in the CAS :', exc_max
end

123
plugins/CAS_SD/cas_sd.irp.f
View File

@ -1,123 +0,0 @@
program full_ci
implicit none
integer :: i,k
double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:)
integer :: N_st, degree
N_st = N_states
allocate (pt2(N_st), norm_pert(N_st),H_pert_diag(N_st))
character*(64) :: perturbation
PROVIDE N_det_cas
pt2 = 1.d0
diag_algorithm = "Lapack"
if (N_det > N_det_max) then
call diagonalize_CI
call save_wavefunction
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
endif
double precision :: i_H_psi_array(N_states),diag_H_mat_elem,h,i_O1_psi_array(N_states)
double precision :: E_CI_before(N_states)
if(read_wf)then
call i_H_psi(psi_det(1,1,N_det),psi_det,psi_coef,N_int,N_det,psi_det_size,N_states,i_H_psi_array)
h = diag_H_mat_elem(psi_det(1,1,N_det),N_int)
selection_criterion = dabs(psi_coef(N_det,1) * (i_H_psi_array(1) - h * psi_coef(N_det,1))) * 0.1d0
soft_touch selection_criterion
endif
integer :: n_det_before
print*,'Beginning the selection ...'
E_CI_before(1:N_states) = CI_energy(1:N_states)
do while (N_det < N_det_max.and.maxval(abs(pt2(1:N_st))) > pt2_max)
n_det_before = N_det
call H_apply_CAS_SD(pt2, norm_pert, H_pert_diag, N_st)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
call diagonalize_CI
if (N_det > N_det_max) then
N_det = N_det_max
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
touch N_det psi_det psi_coef psi_det_sorted psi_coef_sorted psi_average_norm_contrib_sorted
endif
call save_wavefunction
if(n_det_before == N_det)then
selection_criterion = selection_criterion * 0.5d0
endif
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
do k = 1, N_states
print*,'State ',k
print *, 'PT2 = ', pt2(k)
print *, 'E = ', CI_energy(k)
print *, 'E(before)+PT2 = ', E_CI_before(k)+pt2(k)
enddo
print *, '-----'
if(N_states.gt.1)then
print*,'Variational Energy difference'
do i = 2, N_states
print*,'Delta E = ',CI_energy(i) - CI_energy(1)
enddo
endif
if(N_states.gt.1)then
print*,'Variational + perturbative Energy difference'
do i = 2, N_states
print*,'Delta E = ',E_CI_before(i)+ pt2(i) - (E_CI_before(1) + pt2(1))
enddo
endif
E_CI_before(1:N_states) = CI_energy(1:N_states)
call ezfio_set_cas_sd_energy(CI_energy(1))
enddo
N_det = min(N_det_max,N_det)
touch N_det psi_det psi_coef
call diagonalize_CI
if(do_pt2_end)then
print*,'Last iteration only to compute the PT2'
threshold_selectors = 1.d0
threshold_generators = 0.999d0
call H_apply_CAS_SD_PT2(pt2, norm_pert, H_pert_diag, N_st)
print *, 'Final step'
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy(1:N_states)
print *, 'E+PT2 = ', CI_energy(1:N_states)+pt2(1:N_states)
print *, '-----'
call ezfio_set_cas_sd_energy_pt2(CI_energy(1)+pt2(1))
endif
integer :: exc_max, degree_min
exc_max = 0
print *, 'CAS determinants : ', N_det_cas
do i=1,min(N_det_cas,10)
do k=i,N_det_cas
call get_excitation_degree(psi_cas(1,1,k),psi_cas(1,1,i),degree,N_int)
exc_max = max(exc_max,degree)
enddo
print *, psi_coef_cas_diagonalized(i,:)
call debug_det(psi_cas(1,1,i),N_int)
print *, ''
enddo
print *, 'Max excitation degree in the CAS :', exc_max
end

123
plugins/CAS_SD/cas_sd_selected.irp.f
View File

@ -1,123 +0,0 @@
program full_ci
implicit none
integer :: i,k
double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:)
integer :: N_st, degree
N_st = N_states
allocate (pt2(N_st), norm_pert(N_st),H_pert_diag(N_st))
character*(64) :: perturbation
PROVIDE N_det_cas
pt2 = 1.d0
diag_algorithm = "Lapack"
if (N_det > N_det_max) then
call diagonalize_CI
call save_wavefunction
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
endif
double precision :: i_H_psi_array(N_states),diag_H_mat_elem,h,i_O1_psi_array(N_states)
double precision :: E_CI_before(N_states)
if(read_wf)then
call i_H_psi(psi_det(1,1,N_det),psi_det,psi_coef,N_int,N_det,psi_det_size,N_states,i_H_psi_array)
h = diag_H_mat_elem(psi_det(1,1,N_det),N_int)
selection_criterion = dabs(psi_coef(N_det,1) * (i_H_psi_array(1) - h * psi_coef(N_det,1))) * 0.1d0
soft_touch selection_criterion
endif
integer :: n_det_before
print*,'Beginning the selection ...'
E_CI_before(1:N_states) = CI_energy(1:N_states)
do while (N_det < N_det_max.and.maxval(abs(pt2(1:N_st))) > pt2_max)
n_det_before = N_det
call H_apply_CAS_SD_selected(pt2, norm_pert, H_pert_diag, N_st)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
call diagonalize_CI
if (N_det > N_det_max) then
N_det = N_det_max
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
touch N_det psi_det psi_coef psi_det_sorted psi_coef_sorted psi_average_norm_contrib_sorted
endif
call save_wavefunction
if(n_det_before == N_det)then
selection_criterion = selection_criterion * 0.5d0
endif
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
do k = 1, N_states
print*,'State ',k
print *, 'PT2 = ', pt2(k)
print *, 'E = ', CI_energy(k)
print *, 'E(before)+PT2 = ', E_CI_before(k)+pt2(k)
enddo
print *, '-----'
if(N_states.gt.1)then
print*,'Variational Energy difference'
do i = 2, N_states
print*,'Delta E = ',CI_energy(i) - CI_energy(1)
enddo
endif
if(N_states.gt.1)then
print*,'Variational + perturbative Energy difference'
do i = 2, N_states
print*,'Delta E = ',E_CI_before(i)+ pt2(i) - (E_CI_before(1) + pt2(1))
enddo
endif
E_CI_before(1:N_states) = CI_energy(1:N_states)
call ezfio_set_cas_sd_energy(CI_energy(1))
enddo
N_det = min(N_det_max,N_det)
touch N_det psi_det psi_coef
call diagonalize_CI
if(do_pt2_end)then
print*,'Last iteration only to compute the PT2'
threshold_selectors = max(threshold_selectors,threshold_selectors_pt2)
threshold_generators = max(threshold_generators,threshold_generators_pt2)
call H_apply_CAS_SD_PT2(pt2, norm_pert, H_pert_diag, N_st)
print *, 'Final step'
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy(1:N_states)
print *, 'E+PT2 = ', CI_energy(1:N_states)+pt2(1:N_states)
print *, '-----'
call ezfio_set_cas_sd_energy_pt2(CI_energy(1)+pt2(1))
endif
integer :: exc_max, degree_min
exc_max = 0
print *, 'CAS determinants : ', N_det_cas
do i=1,min(N_det_cas,10)
do k=i,N_det_cas
call get_excitation_degree(psi_cas(1,1,k),psi_cas(1,1,i),degree,N_int)
exc_max = max(exc_max,degree)
enddo
print *, psi_cas_coef(i,:)
call debug_det(psi_cas(1,1,i),N_int)
print *, ''
enddo
print *, 'Max excitation degree in the CAS :', exc_max
end

BIN
plugins/CAS_SD/tree_dependency.png
View File

Binary file not shown.
Side by Side

Before

Width:  |  Height:  |  Size: 109 KiB

2
plugins/CAS_SD_ZMQ/cassd_zmq.irp.f
View File

@ -132,7 +132,7 @@ program cassd_zmq
enddo
print *, 'Max excitation degree in the CAS :', exc_max
if(do_pt2_end)then
if(do_pt2)then
print*,'Last iteration only to compute the PT2'
threshold_selectors = max(threshold_selectors,threshold_selectors_pt2)
threshold_generators = max(threshold_generators,threshold_generators_pt2)

4
plugins/DDCI_selected/EZFIO.cfg
View File

@ -1,4 +0,0 @@
[energy]
type: double precision
doc: "Calculated CAS-SD energy"
interface: ezfio

16
plugins/DDCI_selected/H_apply.irp.f
View File

@ -1,16 +0,0 @@
use bitmasks
BEGIN_SHELL [ /usr/bin/env python ]
from generate_h_apply import *
s = H_apply("DDCI_selection")
s.set_selection_pt2("epstein_nesbet_2x2")
s.set_filter_2h_2p()
print s
s = H_apply("DDCI_PT2")
s.set_perturbation("epstein_nesbet_2x2")
s.set_filter_2h_2p()
print s
END_SHELL

1
plugins/DDCI_selected/NEEDED_CHILDREN_MODULES
View File

@ -1 +0,0 @@
Perturbation Selectors_full Generators_CAS Davidson

130
plugins/DDCI_selected/README.rst
View File

@ -1,130 +0,0 @@
====================
DDCI_selected Module
====================
Documentation
=============
.. Do not edit this section. It was auto-generated from the
.. by the `update_README.py` script.
`ddci <http://github.com/LCPQ/quantum_package/tree/master/src/DDCI_selected/ddci.irp.f#L1>`_
Undocumented
`h_apply_ddci_pt2 <http://github.com/LCPQ/quantum_package/tree/master/src/DDCI_selected/H_apply.irp.f_shell_15#L1271>`_
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_ddci_pt2_diexc <http://github.com/LCPQ/quantum_package/tree/master/src/DDCI_selected/H_apply.irp.f_shell_15#L776>`_
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_ddci_pt2_monoexc <http://github.com/LCPQ/quantum_package/tree/master/src/DDCI_selected/H_apply.irp.f_shell_15#L1086>`_
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_ddci_selection <http://github.com/LCPQ/quantum_package/tree/master/src/DDCI_selected/H_apply.irp.f_shell_15#L530>`_
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
`h_apply_ddci_selection_diexc <http://github.com/LCPQ/quantum_package/tree/master/src/DDCI_selected/H_apply.irp.f_shell_15#L1>`_
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
`h_apply_ddci_selection_monoexc <http://github.com/LCPQ/quantum_package/tree/master/src/DDCI_selected/H_apply.irp.f_shell_15#L331>`_
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
Needed Modules
==============
.. Do not edit this section. It was auto-generated from the
.. by the `update_README.py` script.
.. image:: tree_dependency.png
* `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/src/Perturbation>`_
* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/src/Selectors_full>`_
* `Generators_CAS <http://github.com/LCPQ/quantum_package/tree/master/src/Generators_CAS>`_
Needed Modules
==============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
.. image:: tree_dependency.png
* `Perturbation <http://github.com/LCPQ/quantum_package/tree/master/plugins/Perturbation>`_
* `Selectors_full <http://github.com/LCPQ/quantum_package/tree/master/plugins/Selectors_full>`_
* `Generators_CAS <http://github.com/LCPQ/quantum_package/tree/master/plugins/Generators_CAS>`_
Documentation
=============
.. Do not edit this section It was auto-generated
.. by the `update_README.py` script.
`ddci <http://github.com/LCPQ/quantum_package/tree/master/plugins/DDCI_selected/ddci.irp.f#L1>`_
Undocumented
h_apply_ddci_pt2
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_ddci_pt2_diexc
Undocumented
h_apply_ddci_pt2_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_ddci_pt2_diexcp
Undocumented
h_apply_ddci_pt2_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_ddci_selection
Calls H_apply on the HF determinant and selects all connected single and double
excitations (of the same symmetry). Auto-generated by the ``generate_h_apply`` script.
h_apply_ddci_selection_diexc
Undocumented
h_apply_ddci_selection_diexcorg
Generate all double excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.
h_apply_ddci_selection_diexcp
Undocumented
h_apply_ddci_selection_monoexc
Generate all single excitations of key_in using the bit masks of holes and
particles.
Assume N_int is already provided.

77
plugins/DDCI_selected/ddci.irp.f
View File

@ -1,77 +0,0 @@
program ddci
implicit none
integer :: i,k
double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:),E_before(:)
integer :: N_st, degree
N_st = N_states_diag
allocate (pt2(N_st), norm_pert(N_st),H_pert_diag(N_st),E_before(N_st))
character*(64) :: perturbation
pt2 = 1.d0
diag_algorithm = "Lapack"
if (N_det > N_det_max) then
call diagonalize_CI
call save_wavefunction
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
endif
call set_bitmask_particl_as_input(reunion_of_bitmask)
call set_bitmask_hole_as_input(reunion_of_bitmask)
do while (N_det < N_det_max.and.maxval(abs(pt2(1:N_st))) > pt2_max)
call H_apply_DDCI_selection(pt2, norm_pert, H_pert_diag, N_st)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
if (N_det > N_det_max) then
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
endif
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', E_before+pt2
print *, '-----'
if(N_states_diag.gt.1)then
print*,'Variational Energy difference'
do i = 2, N_st
print*,'Delta E = ',CI_energy(i) - CI_energy(1)
enddo
endif
if(N_states.gt.1)then
print*,'Variational + perturbative Energy difference'
do i = 2, N_st
print*,'Delta E = ',E_before(i)+ pt2(i) - (E_before(1) + pt2(1))
enddo
endif
E_before = CI_energy
call ezfio_set_ddci_selected_energy(CI_energy)
enddo
if(do_pt2_end)then
call H_apply_DDCI_pt2(pt2, norm_pert, H_pert_diag, N_st)
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
endif
end

BIN
plugins/DDCI_selected/tree_dependency.png
View File

Binary file not shown.
Side by Side

Before

Width:  |  Height:  |  Size: 112 KiB

16
plugins/Full_CI/H_apply.irp.f
View File

@ -4,24 +4,13 @@ from generate_h_apply import *
s = H_apply("FCI")
s.set_selection_pt2("epstein_nesbet_2x2")
#s.unset_openmp()
s.unset_skip()
print s
s = H_apply("FCI_PT2")
s.set_perturbation("epstein_nesbet_2x2")
s.unset_openmp()
print s
s = H_apply("FCI_PT2_new")
s.set_perturbation("decontracted")
s.unset_openmp()
print s
s = H_apply("FCI_no_skip")
s.set_selection_pt2("epstein_nesbet_2x2")
s.unset_skip()
#s.unset_openmp()
s.unset_openmp()
print s
s = H_apply("FCI_no_selection")
@ -31,6 +20,7 @@ print s
s = H_apply("FCI_mono")
s.set_selection_pt2("epstein_nesbet_2x2")
s.unset_skip()
s.unset_double_excitations()
s.unset_openmp()
print s

2
plugins/Full_CI/full_ci.irp.f
View File

@ -95,7 +95,7 @@ program full_ci
N_det = min(N_det_max,N_det)
touch N_det psi_det psi_coef
call diagonalize_CI
if(do_pt2_end)then
if(do_pt2)then
print*,'Last iteration only to compute the PT2'
threshold_generators = threshold_generators_pt2
threshold_selectors = threshold_selectors_pt2

94
plugins/Full_CI/full_ci_no_skip.irp.f
View File

@ -1,94 +0,0 @@
program full_ci
implicit none
integer :: i,k
double precision, allocatable :: pt2(:), norm_pert(:), H_pert_diag(:)
integer :: N_st, degree
N_st = N_states
allocate (pt2(N_st), norm_pert(N_st),H_pert_diag(N_st))
character*(64) :: perturbation
pt2 = 1.d0
diag_algorithm = "Lapack"
if (N_det > N_det_max) then
call diagonalize_CI
call save_wavefunction
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
endif
double precision :: i_H_psi_array(N_states),diag_H_mat_elem,h,i_O1_psi_array(N_states)
double precision :: E_CI_before(N_states)
if(read_wf)then
call i_H_psi(psi_det(1,1,N_det),psi_det,psi_coef,N_int,N_det,psi_det_size,N_states,i_H_psi_array)
h = diag_H_mat_elem(psi_det(1,1,N_det),N_int)
selection_criterion = dabs(psi_coef(N_det,1) * (i_H_psi_array(1) - h * psi_coef(N_det,1))) * 0.1d0
soft_touch selection_criterion
endif
integer :: n_det_before
print*,'Beginning the selection ...'
E_CI_before = CI_energy
do while (N_det < N_det_max.and.maxval(abs(pt2(1:N_st))) > pt2_max)
n_det_before = N_det
call H_apply_FCI_no_skip(pt2, norm_pert, H_pert_diag, N_st)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
if (N_det > N_det_max) then
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
endif
call diagonalize_CI
call save_wavefunction
if(n_det_before == N_det)then
selection_criterion = selection_criterion * 0.5d0
endif
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', E_CI_before+pt2
print *, '-----'
E_CI_before = CI_energy
call ezfio_set_full_ci_energy(CI_energy)
enddo
N_det = min(N_det_max,N_det)
touch N_det psi_det psi_coef
call diagonalize_CI
if(do_pt2_end)then
print*,'Last iteration only to compute the PT2'
threshold_generators = threshold_generators_pt2
threshold_selectors = threshold_selectors_pt2
SOFT_TOUCH threshold_generators threshold_selectors
! print*,'The thres'
call H_apply_FCI_PT2(pt2, norm_pert, H_pert_diag, N_st)
print *, 'Final step'
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print *, 'PT2 = ', pt2
print *, 'E = ', CI_energy
print *, 'E+PT2 = ', CI_energy+pt2
print *, '-----'
call ezfio_set_full_ci_energy_pt2(CI_energy+pt2)
endif
call save_wavefunction
deallocate(pt2,norm_pert)
end

71
plugins/Full_CI_ZMQ/fci_zmq.irp.f
View File

@ -10,6 +10,9 @@ program fci_zmq
double precision :: hf_energy_ref
logical :: has
double precision :: relative_error
relative_error=1.d-3
pt2 = -huge(1.d0)
threshold_davidson_in = threshold_davidson
threshold_davidson = threshold_davidson_in * 100.d0
@ -42,16 +45,17 @@ program fci_zmq
print *, '-----'
enddo
endif
double precision :: E_CI_before(N_states)
print*,'Beginning the selection ...'
if (.True.) then ! Avoid pre-calculation of CI_energy
E_CI_before(1:N_states) = CI_energy(1:N_states)
endif
n_det_before = 0
character*(8) :: pt2_string
double precision :: correlation_energy_ratio
double precision :: threshold_selectors_save, threshold_generators_save
threshold_selectors_save = threshold_selectors
threshold_generators_save = threshold_generators
correlation_energy_ratio = 0.d0
if (.True.) then ! Avoid pre-calculation of CI_energy
@ -61,8 +65,31 @@ program fci_zmq
(correlation_energy_ratio <= correlation_energy_ratio_max) &
)
if (do_pt2) then
pt2_string = ' '
pt2 = 0.d0
if (N_states == 1) then
threshold_selectors = 1.d0
threshold_generators = 1d0
SOFT_TOUCH threshold_selectors threshold_generators
call ZMQ_pt2(CI_energy, pt2,relative_error) ! Stochastic PT2
threshold_selectors = threshold_selectors_save
threshold_generators = threshold_generators_save
SOFT_TOUCH threshold_selectors threshold_generators
else
threshold_selectors = max(threshold_selectors,threshold_selectors_pt2)
threshold_generators = max(threshold_generators,threshold_generators_pt2)
SOFT_TOUCH threshold_selectors threshold_generators
call ZMQ_selection(0, pt2) ! Deterministic PT2
endif
else
pt2_string = '(approx)'
endif
correlation_energy_ratio = (CI_energy(1) - hf_energy_ref) / &
(E_CI_before(1) + pt2(1) - hf_energy_ref)
(CI_energy(1) + pt2(1) - hf_energy_ref)
correlation_energy_ratio = min(1.d0,correlation_energy_ratio)
@ -74,7 +101,7 @@ program fci_zmq
print*,'State ',k
print *, 'PT2 = ', pt2(k)
print *, 'E = ', CI_energy(k)
print *, 'E(before)+PT2 = ', E_CI_before(k)+pt2(k)
print *, 'E+PT2'//pt2_string//' = ', CI_energy(k)+pt2(k)
enddo
print *, '-----'
@ -87,10 +114,9 @@ program fci_zmq
if(N_states.gt.1)then
print*,'Variational + perturbative Energy difference'
do i = 2, N_states
print*,'Delta E = ',E_CI_before(i)+ pt2(i) - (E_CI_before(1) + pt2(1))
print*,'Delta E = ',CI_energy(i)+ pt2(i) - (CI_energy(1) + pt2(1))
enddo
endif
E_CI_before(1:N_states) = CI_energy(1:N_states)
n_det_before = N_det
to_select = N_det
@ -108,6 +134,8 @@ program fci_zmq
call diagonalize_CI
call save_wavefunction
call ezfio_set_full_ci_zmq_energy(CI_energy(1))
call ezfio_set_full_ci_zmq_energy_pt2(CI_energy(1)+pt2(1))
enddo
endif
@ -118,37 +146,24 @@ program fci_zmq
call ezfio_set_full_ci_zmq_energy(CI_energy(1))
endif
if(do_pt2_end)then
print*,'Last iteration only to compute the PT2'
E_CI_before(1:N_states) = CI_energy(1:N_states)
double precision :: relative_error
relative_error=1.d-3
if (do_pt2) then
pt2 = 0.d0
if (N_states == 1) then
threshold_selectors = 1.d0
threshold_generators = 1d0
SOFT_TOUCH threshold_selectors threshold_generators
print *, 'Stochastic PT2'
call ZMQ_pt2(E_CI_before(1), pt2,relative_error) ! Stochastic PT2
call ZMQ_pt2(CI_energy, pt2, relative_error) ! Stochastic PT2
threshold_selectors = threshold_selectors_save
threshold_generators = threshold_generators_save
SOFT_TOUCH threshold_selectors threshold_generators
else
threshold_selectors = max(threshold_selectors,threshold_selectors_pt2)
threshold_generators = max(threshold_generators,threshold_generators_pt2)
SOFT_TOUCH threshold_selectors threshold_generators
print *, 'Deterministic PT2'
call ZMQ_selection(0, pt2) ! Deterministic PT2
endif
print *, 'Final step'
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
do k=1,N_states
print *, 'State', k
print *, 'PT2 = ', pt2(k)
print *, 'E = ', E_CI_before(k)
print *, 'E+PT2 = ', E_CI_before(k)+pt2(k)
print *, '-----'
enddo
call ezfio_set_full_ci_zmq_energy(E_CI_before(1))
call ezfio_set_full_ci_zmq_energy_pt2(E_CI_before(1)+pt2(1))
call ezfio_set_full_ci_zmq_energy_pt2(CI_energy(1)+pt2(1))
endif
end

146
plugins/Full_CI_ZMQ/fci_zmq_pt2.irp.f
View File

@ -1,146 +0,0 @@
program fci_zmq
implicit none
integer :: i,j,k
double precision, allocatable :: pt2(:)
integer :: degree
integer :: n_det_before, to_select
double precision :: threshold_davidson_in
allocate (pt2(N_states))
double precision :: hf_energy_ref
logical :: has
double precision :: relative_error
relative_error=1.d-3
pt2 = -huge(1.d0)
threshold_davidson_in = threshold_davidson
threshold_davidson = threshold_davidson_in * 100.d0
SOFT_TOUCH threshold_davidson
call diagonalize_CI
call save_wavefunction
call ezfio_has_hartree_fock_energy(has)
if (has) then
call ezfio_get_hartree_fock_energy(hf_energy_ref)
else
hf_energy_ref = ref_bitmask_energy
endif
if (N_det > N_det_max) then
psi_det = psi_det_sorted
psi_coef = psi_coef_sorted
N_det = N_det_max
soft_touch N_det psi_det psi_coef
call diagonalize_CI
call save_wavefunction
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
do k=1,N_states
print*,'State ',k
print *, 'PT2 = ', pt2(k)
print *, 'E = ', CI_energy(k)
print *, 'E+PT2 = ', CI_energy(k) + pt2(k)
print *, '-----'
enddo
endif
print*,'Beginning the selection ...'
n_det_before = 0
double precision :: correlation_energy_ratio
double precision :: threshold_selectors_save, threshold_generators_save
threshold_selectors_save = threshold_selectors
threshold_generators_save = threshold_generators
correlation_energy_ratio = 0.d0
if (.True.) then ! Avoid pre-calculation of CI_energy
do while ( &
(N_det < N_det_max) .and. &
(maxval(abs(pt2(1:N_states))) > pt2_max) .and. &
(correlation_energy_ratio <= correlation_energy_ratio_max) &
)
pt2 = 0.d0
if (N_states == 1) then
threshold_selectors = 1.d0
threshold_generators = 1d0
SOFT_TOUCH threshold_selectors threshold_generators
call ZMQ_pt2(CI_energy, pt2,relative_error) ! Stochastic PT2
threshold_selectors = threshold_selectors_save
threshold_generators = threshold_generators_save
SOFT_TOUCH threshold_selectors threshold_generators
else
threshold_selectors = max(threshold_selectors,threshold_selectors_pt2)
threshold_generators = max(threshold_generators,threshold_generators_pt2)
SOFT_TOUCH threshold_selectors threshold_generators
call ZMQ_selection(0, pt2) ! Deterministic PT2
endif
correlation_energy_ratio = (CI_energy(1) - hf_energy_ref) / &
(CI_energy(1) + pt2(1) - hf_energy_ref)
correlation_energy_ratio = min(1.d0,correlation_energy_ratio)
print *, 'N_det = ', N_det
print *, 'N_states = ', N_states
print*, 'correlation_ratio = ', correlation_energy_ratio
do k=1, N_states
print*,'State ',k
print *, 'PT2 = ', pt2(k)
print *, 'E = ', CI_energy(k)
print *, 'E+PT2 = ', CI_energy(k)+pt2(k)
enddo
print *, '-----'
if(N_states.gt.1)then
print*,'Variational Energy difference'
do i = 2, N_states
print*,'Delta E = ',CI_energy(i) - CI_energy(1)
enddo
endif
if(N_states.gt.1)then
print*,'Variational + perturbative Energy difference'
do i = 2, N_states
print*,'Delta E = ',CI_energy(i)+ pt2(i) - (CI_energy(1) + pt2(1))
enddo
endif
n_det_before = N_det
to_select = N_det
to_select = max(N_det, to_select)
to_select = min(to_select, N_det_max-n_det_before)
call ZMQ_selection(to_select, pt2)
PROVIDE psi_coef
PROVIDE psi_det
PROVIDE psi_det_sorted
if (N_det >= N_det_max) then
threshold_davidson = threshold_davidson_in
end if
call diagonalize_CI
call save_wavefunction
call ezfio_set_full_ci_zmq_energy(CI_energy(1))
call ezfio_set_full_ci_zmq_energy_pt2(CI_energy(1)+pt2(1))
enddo
endif
if (N_det < N_det_max) then
threshold_davidson = threshold_davidson_in
call diagonalize_CI
call save_wavefunction
call ezfio_set_full_ci_zmq_energy(CI_energy(1))
call ezfio_set_full_ci_zmq_energy_pt2(CI_energy(1)+pt2(1))
endif
end

149
plugins/Full_CI_ZMQ/pt2_stoch_routines.irp.f
View File

@ -9,96 +9,102 @@ subroutine ZMQ_pt2(E, pt2,relative_error)
implicit none
character(len=64000) :: task
character(len=64000) :: task
integer(ZMQ_PTR) :: zmq_to_qp_run_socket, zmq_to_qp_run_socket2
type(selection_buffer) :: b
integer, external :: omp_get_thread_num
double precision, intent(in) :: relative_error, E
double precision, intent(out) :: pt2(N_states)
double precision, allocatable :: pt2_detail(:,:), comb(:)
logical, allocatable :: computed(:)
integer, allocatable :: tbc(:)
integer :: i, j, k, Ncomb, generator_per_task, i_generator_end
integer, external :: pt2_find
double precision :: sumabove(comb_teeth), sum2above(comb_teeth), Nabove(comb_teeth)
double precision, external :: omp_get_wtime
double precision :: time
allocate(pt2_detail(N_states, N_det_generators+1), comb(N_det_generators), computed(N_det_generators), tbc(0:size_tbc))
sumabove = 0d0
sum2above = 0d0
Nabove = 0d0
provide nproc fragment_first fragment_count mo_bielec_integrals_in_map mo_mono_elec_integral pt2_weight psi_selectors
computed = .false.
tbc(0) = first_det_of_comb - 1
do i=1, tbc(0)
tbc(i) = i
computed(i) = .true.
end do
double precision, allocatable :: pt2_detail(:,:), comb(:)
logical, allocatable :: computed(:)
integer, allocatable :: tbc(:)
integer :: i, j, k, Ncomb, generator_per_task, i_generator_end
integer, external :: pt2_find
pt2_detail = 0d0
generator_per_task = 1
print *, '========== ================= ================= ================='
print *, ' Samples Energy Stat. Error Seconds '
print *, '========== ================= ================= ================='
double precision :: sumabove(comb_teeth), sum2above(comb_teeth), Nabove(comb_teeth)
double precision, external :: omp_get_wtime
double precision :: time
if (N_det < 10) then
call ZMQ_selection(0, pt2)
return
else
call new_parallel_job(zmq_to_qp_run_socket,'pt2')
call zmq_put_psi(zmq_to_qp_run_socket,1,pt2_e0_denominator,size(pt2_e0_denominator))
call create_selection_buffer(1, 1*2, b)
Ncomb=size(comb)
call get_carlo_workbatch(computed, comb, Ncomb, tbc)
integer(ZMQ_PTR), external :: new_zmq_to_qp_run_socket
integer :: ipos
ipos=1
do i=1,tbc(0)
if(tbc(i) > fragment_first) then
write(task(ipos:ipos+20),'(I9,1X,I9,''|'')') 0, tbc(i)
ipos += 20
if (ipos > 63980) then
call add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos)))
ipos=1
endif
else
do j=1,fragment_count
write(task(ipos:ipos+20),'(I9,1X,I9,''|'')') j, tbc(i)
allocate(pt2_detail(N_states, N_det_generators+1), comb(N_det_generators), computed(N_det_generators), tbc(0:size_tbc))
sumabove = 0d0
sum2above = 0d0
Nabove = 0d0
provide nproc fragment_first fragment_count mo_bielec_integrals_in_map mo_mono_elec_integral pt2_weight psi_selectors
computed = .false.
tbc(0) = first_det_of_comb - 1
do i=1, tbc(0)
tbc(i) = i
computed(i) = .true.
end do
pt2_detail = 0d0
generator_per_task = 1
print *, '========== ================= ================= ================='
print *, ' Samples Energy Stat. Error Seconds '
print *, '========== ================= ================= ================='
call new_parallel_job(zmq_to_qp_run_socket,'pt2')
call zmq_put_psi(zmq_to_qp_run_socket,1,pt2_e0_denominator,size(pt2_e0_denominator))
call create_selection_buffer(1, 1*2, b)
Ncomb=size(comb)
call get_carlo_workbatch(computed, comb, Ncomb, tbc)
integer(ZMQ_PTR), external :: new_zmq_to_qp_run_socket
integer :: ipos
ipos=1
do i=1,tbc(0)
if(tbc(i) > fragment_first) then
write(task(ipos:ipos+20),'(I9,1X,I9,''|'')') 0, tbc(i)
ipos += 20
if (ipos > 63980) then
call add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos)))
ipos=1
endif
end do
end if
end do
if (ipos > 1) then
call add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos)))
endif
call zmq_set_running(zmq_to_qp_run_socket)
!$OMP PARALLEL DEFAULT(shared) NUM_THREADS(nproc+1) &
!$OMP PRIVATE(i)
else
do j=1,fragment_count
write(task(ipos:ipos+20),'(I9,1X,I9,''|'')') j, tbc(i)
ipos += 20
if (ipos > 63980) then
call add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos)))
ipos=1
endif
end do
end if
end do
if (ipos > 1) then
call add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos)))
endif
call zmq_set_running(zmq_to_qp_run_socket)
!$OMP PARALLEL DEFAULT(shared) NUM_THREADS(nproc+1) &
!$OMP PRIVATE(i)
i = omp_get_thread_num()
if (i==0) then
call pt2_collector(E, b, tbc, comb, Ncomb, computed, pt2_detail, sumabove, sum2above, Nabove, relative_error, pt2)
else
call pt2_slave_inproc(i)
endif
!$OMP END PARALLEL
call delete_selection_buffer(b)
call end_parallel_job(zmq_to_qp_run_socket, 'pt2')
print *, '========== ================= ================= ================='
deallocate(pt2_detail, comb, computed, tbc)
!$OMP END PARALLEL
call delete_selection_buffer(b)
call end_parallel_job(zmq_to_qp_run_socket, 'pt2')
print *, '========== ================= ================= ================='
deallocate(pt2_detail, comb, computed, tbc)
endif
end subroutine
@ -279,7 +285,6 @@ subroutine pt2_collector(E, b, tbc, comb, Ncomb, computed, pt2_detail, sumabove,
prop = prop * pt2_weight_inv(first_det_of_teeth(tooth))
E0 += pt2_detail(1,first_det_of_teeth(tooth)) * prop
pt2(1) = E0 + (sumabove(tooth) / Nabove(tooth))
eqt = sqrt(1d0 / (Nabove(tooth)-1) * abs(sum2above(tooth) / Nabove(tooth) - (sumabove(tooth)/Nabove(tooth))**2))
call end_zmq_to_qp_run_socket(zmq_to_qp_run_socket)
call end_zmq_pull_socket(zmq_socket_pull)
@ -471,7 +476,7 @@ end subroutine
end if
norm_left -= pt2_weight(i)
end do
first_det_of_comb = max(1,first_det_of_comb)
first_det_of_comb = max(2,first_det_of_comb)
call write_int(6, first_det_of_comb-1, 'Size of deterministic set')
comb_step = (1d0 - pt2_cweight(first_det_of_comb-1)) * comb_step

4
plugins/Perturbation/EZFIO.cfg
View File

@ -1,6 +1,6 @@
[do_pt2_end]
[do_pt2]
type: logical
doc: If true, compute the PT2 at the end of the selection
doc: If true, compute the PT2
interface: ezfio,provider,ocaml
default: True

7
plugins/mrcepa0/dressing.irp.f
View File

@ -410,13 +410,6 @@ end
END_PROVIDER
BEGIN_PROVIDER [ integer, HP, (2,N_det_non_ref) ]
integer :: i
do i=1,N_det_non_ref
call getHP(psi_non_ref(1,1,i), HP(1,i), HP(2,i), N_int)
end do
END_PROVIDER
BEGIN_PROVIDER [ integer, cepa0_shortcut, (0:N_det_non_ref+1) ]
&BEGIN_PROVIDER [ integer, det_cepa0_idx, (N_det_non_ref) ]
&BEGIN_PROVIDER [ integer(bit_kind), det_cepa0_active, (N_int,2,N_det_non_ref) ]

8
plugins/mrcepa0/dressing_slave.irp.f
View File

@ -39,6 +39,7 @@ subroutine mrsc2_dressing_slave(thread,iproc)
double precision, allocatable :: delta(:,:,:), delta_s2(:,:,:)
integer, allocatable :: hp(:,:)
integer :: i_state, i, i_I, J, k, k2, k1, kk, ll, degree, degree2, m, l, deg, ni, m2
@ -66,6 +67,12 @@ subroutine mrsc2_dressing_slave(thread,iproc)
allocate (dleat_s2(N_states, N_det_non_ref, 2), delta_s2(N_states,0:N_det_non_ref, 2))
allocate(komon(0:N_det_non_ref))
allocate(hp(2,N_det_non_ref))
do i=1,N_det_non_ref
call getHP(psi_non_ref(1,1,i), HP(1,i), HP(2,i), N_int)
end do
do
call get_task_from_taskserver(zmq_to_qp_run_socket,worker_id, task_id, task)
if (task_id == 0) exit
@ -171,7 +178,6 @@ subroutine mrsc2_dressing_slave(thread,iproc)
call apply_excitation(psi_non_ref(1,1,i),exc_Ik,det_tmp,ok,N_int)
if(.not. ok) cycle
if(HP(1,i) + HP(1,k) <= 2 .and. HP(2,i) + HP(2,k) <= 2) then
! if(is_in_wavefunction(det_tmp, N_int)) cycle
cycle
end if

2
plugins/mrcepa0/mrcc.irp.f
View File

@ -19,7 +19,7 @@ program mrsc2sub
SOFT_TOUCH psi_coef
endif
call run(N_states,energy)
if(do_pt2_end)then
if(do_pt2)then
call run_pt2(N_states,energy)
endif
deallocate(energy)

2
plugins/mrcepa0/mrcepa0.irp.f
View File

@ -21,7 +21,7 @@ program mrcepa0
call print_cas_coefs
call run(N_states,energy)
if(do_pt2_end)then
if(do_pt2)then
call run_pt2(N_states,energy)
endif
deallocate(energy)

2
plugins/mrcepa0/mrsc2.irp.f
View File

@ -18,7 +18,7 @@ program mrsc2
TOUCH psi_coef
endif
call run(N_states,energy)
if(do_pt2_end)then
if(do_pt2)then
call run_pt2(N_states,energy)
endif
deallocate(energy)

6
src/MO_Basis/ao_ortho_canonical.irp.f
View File

@ -106,9 +106,9 @@ END_PROVIDER
ao_ortho_canonical_coef(i,i) = 1.d0
enddo
call ortho_lowdin(ao_overlap,size(ao_overlap,1),ao_num,ao_ortho_canonical_coef,size(ao_ortho_canonical_coef,1),ao_num)
ao_ortho_canonical_num=ao_num
return
!call ortho_lowdin(ao_overlap,size(ao_overlap,1),ao_num,ao_ortho_canonical_coef,size(ao_ortho_canonical_coef,1),ao_num)
!ao_ortho_canonical_num=ao_num
!return
if (ao_cartesian) then

2
src/Utils/LinearAlgebra.irp.f
View File

@ -361,7 +361,7 @@ subroutine find_rotation(A,LDA,B,m,C,n)
double precision, allocatable :: A_inv(:,:)
allocate(A_inv(LDA,n))
call get_pseudo_inverse(A,m,n,A_inv,LDA)
call get_pseudo_inverse(A,LDA,m,n,A_inv,LDA)
integer :: i,j,k
call dgemm('N','N',n,n,m,1.d0,A_inv,n,B,LDA,0.d0,C,n)

5
src/ZMQ/utils.irp.f
View File

@ -844,7 +844,7 @@ subroutine get_task_from_taskserver(zmq_to_qp_run_socket,worker_id,task_id,task)
integer, intent(out) :: task_id
character*(512), intent(out) :: task
character*(512) :: message
character*(1024) :: message
character*(64) :: reply
integer :: rc, sze
@ -858,7 +858,8 @@ subroutine get_task_from_taskserver(zmq_to_qp_run_socket,worker_id,task_id,task)
endif
message = repeat(' ',512)
rc = f77_zmq_recv(zmq_to_qp_run_socket, message, 510, 0)
rc = f77_zmq_recv(zmq_to_qp_run_socket, message, 1024, 0)
rc = min(1024,rc)
read(message(1:rc),*) reply
if (trim(reply) == 'get_task_reply') then
read(message(1:rc),*) reply, task_id

4
tests/bats/cassd.bats
View File

@ -8,7 +8,7 @@ source $QP_ROOT/tests/bats/common.bats.sh
rm -rf work/h2o.ezfio/determinants/
qp_edit -c $INPUT
ezfio set_file $INPUT
ezfio set perturbation do_pt2_end True
ezfio set perturbation do_pt2 True
ezfio set determinants n_det_max 16384
qp_set_mo_class $INPUT -core "[1]" -inact "[2,5]" -act "[3,4,6,7]" -virt "[8-24]"
qp_run cassd_zmq $INPUT
@ -17,7 +17,7 @@ source $QP_ROOT/tests/bats/common.bats.sh
ezfio set determinants n_det_max 1024
ezfio set determinants read_wf True
ezfio set perturbation do_pt2_end True
ezfio set perturbation do_pt2 True
qp_run cassd_zmq $INPUT
ezfio set determinants read_wf False
energy="$(ezfio get cas_sd_zmq energy)"

4
tests/bats/fci.bats
View File

@ -7,7 +7,7 @@ function run_FCI() {
test_exe full_ci || skip
qp_edit -c $1
ezfio set_file $1
ezfio set perturbation do_pt2_end True
ezfio set perturbation do_pt2 True
ezfio set determinants n_det_max $2
ezfio set davidson threshold_davidson 1.e-10
@ -23,7 +23,7 @@ function run_FCI_ZMQ() {
test_exe fci_zmq || skip
qp_edit -c $1
ezfio set_file $1
ezfio set perturbation do_pt2_end True
ezfio set perturbation do_pt2 True
ezfio set determinants n_det_max $2
ezfio set davidson threshold_davidson 1.e-10

2
tests/bats/pseudo.bats
View File

@ -26,7 +26,7 @@ function run_FCI_ZMQ() {
test_exe fci_zmq|| skip
qp_edit -c $1
ezfio set_file $1
ezfio set perturbation do_pt2_end True
ezfio set perturbation do_pt2 True
ezfio set determinants n_det_max $2
ezfio set davidson threshold_davidson 1.e-10