mirror of
https://github.com/LCPQ/quantum_package
synced 2024-06-02 11:25:26 +02:00
Renaming in lowercase
This commit is contained in:
parent
67e71e5931
commit
ba839928db
|
|
@ -4,42 +4,41 @@ irpf90.make
|
|||
irpf90_entities
|
||||
tags
|
||||
Makefile
|
||||
AO_Basis
|
||||
AO_one_e_integrals
|
||||
Bitmask
|
||||
CIS
|
||||
CISD
|
||||
Davidson
|
||||
DavidsonDressed
|
||||
DavidsonUndressed
|
||||
Determinants
|
||||
Dressing
|
||||
ao_basis
|
||||
ao_one_e_integrals
|
||||
bitmask
|
||||
cis
|
||||
cisd
|
||||
davidson
|
||||
davidson_dressed
|
||||
davidson_undressed
|
||||
determinants
|
||||
dressing
|
||||
dummy
|
||||
Electrons
|
||||
Ezfio_files
|
||||
FCI
|
||||
FourIdx
|
||||
Generators_CAS
|
||||
Generators_full
|
||||
gitignore
|
||||
Hartree_Fock
|
||||
Integrals_Bielec
|
||||
Iterations
|
||||
MO_Basis
|
||||
MOGuess
|
||||
MO_one_e_integrals
|
||||
MPI
|
||||
MRPT_Utils
|
||||
Nuclei
|
||||
Perturbation
|
||||
Pseudo
|
||||
Psiref_CAS
|
||||
Psiref_Utils
|
||||
Selectors_CASSD
|
||||
Selectors_full
|
||||
Selectors_Utils
|
||||
SingleRefMethod
|
||||
Slave
|
||||
Tools
|
||||
Utils
|
||||
ZMQ
|
||||
electrons
|
||||
ezfio_files
|
||||
fci
|
||||
four_idx
|
||||
generators_cas
|
||||
generators_full
|
||||
hartree_fock
|
||||
integrals_bielec
|
||||
iterations
|
||||
mo_basis
|
||||
mo_guess
|
||||
mo_one_e_integrals
|
||||
mpi
|
||||
mrpt_utils
|
||||
nuclei
|
||||
perturbation
|
||||
pseudo
|
||||
psiref_cas
|
||||
psiref_utils
|
||||
selectors_cassd
|
||||
selectors_full
|
||||
selectors_utils
|
||||
single_ref_method
|
||||
slave
|
||||
tools
|
||||
utils
|
||||
zmq
|
||||
|
|
|
|||
|
|
@ -1 +0,0 @@
|
|||
Nuclei
|
||||
|
|
@ -1,2 +0,0 @@
|
|||
AO_Basis
|
||||
Pseudo
|
||||
|
|
@ -1 +0,0 @@
|
|||
MO_Basis
|
||||
|
|
@ -1,3 +0,0 @@
|
|||
Selectors_full
|
||||
SingleRefMethod
|
||||
DavidsonUndressed
|
||||
|
|
@ -1,3 +0,0 @@
|
|||
Selectors_full
|
||||
SingleRefMethod
|
||||
DavidsonUndressed
|
||||
|
|
@ -1 +0,0 @@
|
|||
Determinants
|
||||
|
|
@ -1 +0,0 @@
|
|||
Davidson
|
||||
|
|
@ -1 +0,0 @@
|
|||
Davidson
|
||||
|
|
@ -1,3 +0,0 @@
|
|||
MO_Basis
|
||||
MO_one_e_integrals
|
||||
Integrals_Bielec
|
||||
|
|
@ -1 +0,0 @@
|
|||
ZMQ
|
||||
|
|
@ -1 +0,0 @@
|
|||
Ezfio_files
|
||||
|
|
@ -1 +0,0 @@
|
|||
MPI
|
||||
|
|
@ -1,8 +0,0 @@
|
|||
===========
|
||||
Ezfio_files
|
||||
===========
|
||||
|
||||
This modules essentially contains the name of the |EZFIO| directory in the
|
||||
:c:data:`ezfio_filename` variable. This is read as the first argument of the
|
||||
command-line, or as the :envvar:`QP_INPUT` environment variable.
|
||||
|
||||
|
|
@ -1,53 +0,0 @@
|
|||
BEGIN_PROVIDER [ character*(128), ezfio_filename ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Name of EZFIO file. It is obtained from the QPACKAGE_INPUT environment
|
||||
! variable if it is set, or as the 1st argument of the command line.
|
||||
END_DOC
|
||||
|
||||
PROVIDE mpi_initialized
|
||||
|
||||
! Get the QPACKAGE_INPUT environment variable
|
||||
call getenv('QPACKAGE_INPUT',ezfio_filename)
|
||||
if (ezfio_filename == '') then
|
||||
! Get from the command line
|
||||
integer :: iargc
|
||||
call getarg(0,ezfio_filename)
|
||||
if (iargc() /= 1) then
|
||||
print *, 'Missing EZFIO file name in the command line:'
|
||||
print *, trim(ezfio_filename)//' <ezfio_file>'
|
||||
stop 1
|
||||
endif
|
||||
call getarg(1,ezfio_filename)
|
||||
endif
|
||||
|
||||
! Check that file exists
|
||||
logical :: exists
|
||||
inquire(file=trim(ezfio_filename)//'/ezfio/creation',exist=exists)
|
||||
if (.not.exists) then
|
||||
print *, 'Error: file '//trim(ezfio_filename)//' does not exist'
|
||||
stop 1
|
||||
endif
|
||||
|
||||
call ezfio_set_file(ezfio_filename)
|
||||
|
||||
! Adjust out-of-memory killer flag such that the current process will be
|
||||
! killed first by the OOM killer, allowing compute nodes to survive
|
||||
integer :: getpid
|
||||
character*(64) :: command, pidc
|
||||
write(pidc,*) getpid()
|
||||
write(command,*) 'echo 15 > /proc//'//trim(adjustl(pidc))//'/oom_adj'
|
||||
call system(command)
|
||||
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ character*(128), ezfio_work_dir ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! EZFIO/work/
|
||||
END_DOC
|
||||
call ezfio_set_work_empty(.False.)
|
||||
ezfio_work_dir = trim(ezfio_filename)//'/work/'
|
||||
END_PROVIDER
|
||||
|
||||
|
|
@ -1,58 +0,0 @@
|
|||
integer function getUnitAndOpen(f,mode)
|
||||
implicit none
|
||||
|
||||
BEGIN_DOC
|
||||
! :f:
|
||||
! file name
|
||||
!
|
||||
! :mode:
|
||||
! 'R' : READ, UNFORMATTED
|
||||
! 'W' : WRITE, UNFORMATTED
|
||||
! 'r' : READ, FORMATTED
|
||||
! 'w' : WRITE, FORMATTED
|
||||
! 'a' : APPEND, FORMATTED
|
||||
! 'x' : READ/WRITE, FORMATTED
|
||||
!
|
||||
END_DOC
|
||||
|
||||
character*(*) :: f
|
||||
character*(128) :: new_f
|
||||
integer :: iunit
|
||||
logical :: is_open, exists
|
||||
character :: mode
|
||||
|
||||
is_open = .True.
|
||||
iunit = 10
|
||||
new_f = f
|
||||
do while (is_open)
|
||||
inquire(unit=iunit,opened=is_open)
|
||||
if (.not.is_open) then
|
||||
getUnitAndOpen = iunit
|
||||
endif
|
||||
iunit = iunit+1
|
||||
enddo
|
||||
if (mode.eq.'r') then
|
||||
inquire(file=f,exist=exists)
|
||||
if (.not.exists) then
|
||||
open(unit=getUnitAndOpen,file=f,status='NEW',action='WRITE',form='FORMATTED')
|
||||
close(unit=getUnitAndOpen)
|
||||
endif
|
||||
open(unit=getUnitAndOpen,file=f,status='OLD',action='READ',form='FORMATTED')
|
||||
else if (mode.eq.'R') then
|
||||
inquire(file=f,exist=exists)
|
||||
if (.not.exists) then
|
||||
open(unit=getUnitAndOpen,file=f,status='NEW',action='WRITE',form='UNFORMATTED')
|
||||
close(unit=getUnitAndOpen)
|
||||
endif
|
||||
open(unit=getUnitAndOpen,file=f,status='OLD',action='READ',form='UNFORMATTED')
|
||||
else if (mode.eq.'W') then
|
||||
open(unit=getUnitAndOpen,file=new_f,status='UNKNOWN',action='WRITE',form='UNFORMATTED')
|
||||
else if (mode.eq.'w') then
|
||||
open(unit=getUnitAndOpen,file=new_f,status='UNKNOWN',action='WRITE',form='FORMATTED')
|
||||
else if (mode.eq.'a') then
|
||||
open(unit=getUnitAndOpen,file=new_f,status='UNKNOWN',action='WRITE',position='APPEND',form='FORMATTED')
|
||||
else if (mode.eq.'x') then
|
||||
open(unit=getUnitAndOpen,file=new_f,form='FORMATTED')
|
||||
endif
|
||||
end function getUnitAndOpen
|
||||
|
||||
|
|
@ -1,85 +0,0 @@
|
|||
BEGIN_PROVIDER [ double precision, output_wall_time_0 ]
|
||||
&BEGIN_PROVIDER [ double precision, output_cpu_time_0 ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Initial CPU and wall times when printing in the output files
|
||||
END_DOC
|
||||
call cpu_time(output_wall_time_0)
|
||||
call wall_time(output_wall_time_0)
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
subroutine write_time(iunit)
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Write a time stamp in the output for chronological reconstruction
|
||||
END_DOC
|
||||
integer, intent(in) :: iunit
|
||||
double precision :: wt, ct
|
||||
if (.not.mpi_master) then
|
||||
return
|
||||
endif
|
||||
write(6,*)
|
||||
call print_memory_usage()
|
||||
call cpu_time(ct)
|
||||
call wall_time(wt)
|
||||
write(6,'(A,F15.6,A,F15.6,A)') &
|
||||
'.. >>>>> [ WALL TIME: ', wt-output_wall_time_0, &
|
||||
' s ] [ CPU TIME: ', ct-output_cpu_time_0, ' s ] <<<<< ..'
|
||||
write(6,*)
|
||||
end
|
||||
|
||||
subroutine write_double(iunit,value,label)
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Write a double precision value in output
|
||||
END_DOC
|
||||
if (.not.mpi_master) then
|
||||
return
|
||||
endif
|
||||
integer, intent(in) :: iunit
|
||||
double precision :: value
|
||||
character*(*) :: label
|
||||
character*(64), parameter :: f = '(A50,G24.16)'
|
||||
character*(50) :: newlabel
|
||||
write(newlabel,'(A,A)') '* ',trim(label)
|
||||
write(6,f) newlabel, value
|
||||
end
|
||||
|
||||
|
||||
subroutine write_int(iunit,value,label)
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Write an integer value in output
|
||||
END_DOC
|
||||
if (.not.mpi_master) then
|
||||
return
|
||||
endif
|
||||
integer, intent(in) :: iunit
|
||||
integer :: value
|
||||
character*(*) :: label
|
||||
character*(64), parameter :: f = '(A50,I16)'
|
||||
character*(50) :: newlabel
|
||||
write(newlabel,'(A,A)') '* ',trim(label)
|
||||
write(6,f) newlabel, value
|
||||
end
|
||||
|
||||
|
||||
subroutine write_bool(iunit,value,label)
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Write an logical value in output
|
||||
END_DOC
|
||||
if (.not.mpi_master) then
|
||||
return
|
||||
endif
|
||||
integer, intent(in) :: iunit
|
||||
logical :: value
|
||||
character*(*) :: label
|
||||
character*(64), parameter :: f = '(A50,L1)'
|
||||
character*(50) :: newlabel
|
||||
write(newlabel,'(A,A)') '* ',trim(label)
|
||||
write(6,f) newlabel, value
|
||||
end
|
||||
|
||||
|
||||
|
|
@ -1,6 +0,0 @@
|
|||
work
|
||||
empty logical
|
||||
|
||||
save
|
||||
empty logical
|
||||
|
||||
|
|
@ -1,13 +0,0 @@
|
|||
[energy]
|
||||
type: double precision
|
||||
doc: Calculated Selected |FCI| energy
|
||||
interface: ezfio
|
||||
size: (determinants.n_states)
|
||||
|
||||
[energy_pt2]
|
||||
type: double precision
|
||||
doc: Calculated |FCI| energy + |PT2|
|
||||
interface: ezfio
|
||||
size: (determinants.n_states)
|
||||
|
||||
|
||||
|
|
@ -1,147 +0,0 @@
|
|||
program fci_zmq
|
||||
implicit none
|
||||
integer :: i,j,k
|
||||
double precision, allocatable :: pt2(:), variance(:), norm(:), rpt2(:)
|
||||
integer :: n_det_before, to_select
|
||||
|
||||
double precision :: rss
|
||||
double precision, external :: memory_of_double
|
||||
rss = memory_of_double(N_states)*4.d0
|
||||
call check_mem(rss,irp_here)
|
||||
|
||||
allocate (pt2(N_states), rpt2(N_states), norm(N_states), variance(N_states))
|
||||
|
||||
double precision :: hf_energy_ref
|
||||
logical :: has
|
||||
double precision :: relative_error
|
||||
|
||||
PROVIDE H_apply_buffer_allocated
|
||||
|
||||
relative_error=PT2_relative_error
|
||||
|
||||
pt2 = -huge(1.e0)
|
||||
rpt2 = -huge(1.e0)
|
||||
norm = 0.d0
|
||||
variance = 0.d0
|
||||
|
||||
if (s2_eig) then
|
||||
call make_s2_eigenfunction
|
||||
endif
|
||||
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
|
||||
endif
|
||||
|
||||
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
|
||||
double precision :: error(N_states)
|
||||
|
||||
correlation_energy_ratio = 0.d0
|
||||
|
||||
do while ( &
|
||||
(N_det < N_det_max) .and. &
|
||||
(maxval(abs(pt2(1:N_states))) > pt2_max) .and. &
|
||||
(correlation_energy_ratio <= correlation_energy_ratio_max) &
|
||||
)
|
||||
write(*,'(A)') '--------------------------------------------------------------------------------'
|
||||
|
||||
|
||||
if (do_pt2) then
|
||||
pt2 = 0.d0
|
||||
variance = 0.d0
|
||||
norm = 0.d0
|
||||
threshold_selectors = 1.d0
|
||||
threshold_generators = 1.d0
|
||||
SOFT_TOUCH threshold_selectors threshold_generators
|
||||
call ZMQ_pt2(psi_energy_with_nucl_rep,pt2,relative_error,error, variance, norm) ! Stochastic PT2
|
||||
threshold_selectors = threshold_selectors_save
|
||||
threshold_generators = threshold_generators_save
|
||||
SOFT_TOUCH threshold_selectors threshold_generators
|
||||
endif
|
||||
|
||||
|
||||
correlation_energy_ratio = (psi_energy_with_nucl_rep(1) - hf_energy_ref) / &
|
||||
(psi_energy_with_nucl_rep(1) + pt2(1) - hf_energy_ref)
|
||||
correlation_energy_ratio = min(1.d0,correlation_energy_ratio)
|
||||
|
||||
call ezfio_set_fci_energy_pt2(psi_energy_with_nucl_rep+pt2)
|
||||
call write_double(6,correlation_energy_ratio, 'Correlation ratio')
|
||||
call print_summary(psi_energy_with_nucl_rep(1:N_states),pt2,error,variance,norm)
|
||||
|
||||
do k=1,N_states
|
||||
rpt2(:) = pt2(:)/(1.d0 + norm(k))
|
||||
enddo
|
||||
|
||||
call save_iterations(psi_energy_with_nucl_rep(1:N_states),rpt2,N_det)
|
||||
call print_extrapolated_energy(psi_energy_with_nucl_rep(1:N_states),rpt2)
|
||||
N_iter += 1
|
||||
|
||||
n_det_before = N_det
|
||||
to_select = N_det
|
||||
to_select = max(N_states_diag, to_select)
|
||||
to_select = min(to_select, N_det_max-n_det_before)
|
||||
call ZMQ_selection(to_select, pt2, variance, norm)
|
||||
|
||||
PROVIDE psi_coef
|
||||
PROVIDE psi_det
|
||||
PROVIDE psi_det_sorted
|
||||
|
||||
call diagonalize_CI
|
||||
call save_wavefunction
|
||||
call ezfio_set_fci_energy(psi_energy_with_nucl_rep(1:N_states))
|
||||
enddo
|
||||
|
||||
if (N_det < N_det_max) then
|
||||
call diagonalize_CI
|
||||
call save_wavefunction
|
||||
call ezfio_set_fci_energy(psi_energy_with_nucl_rep(1:N_states))
|
||||
call ezfio_set_fci_energy_pt2(psi_energy_with_nucl_rep(1:N_states)+pt2)
|
||||
endif
|
||||
|
||||
if (do_pt2) then
|
||||
pt2 = 0.d0
|
||||
variance = 0.d0
|
||||
norm = 0.d0
|
||||
threshold_selectors = 1.d0
|
||||
threshold_generators = 1d0
|
||||
SOFT_TOUCH threshold_selectors threshold_generators
|
||||
call ZMQ_pt2(psi_energy_with_nucl_rep, pt2,relative_error,error,variance,norm) ! Stochastic PT2
|
||||
threshold_selectors = threshold_selectors_save
|
||||
threshold_generators = threshold_generators_save
|
||||
SOFT_TOUCH threshold_selectors threshold_generators
|
||||
call ezfio_set_fci_energy(psi_energy_with_nucl_rep(1:N_states))
|
||||
call ezfio_set_fci_energy_pt2(psi_energy_with_nucl_rep(1:N_states)+pt2)
|
||||
endif
|
||||
print *, 'N_det = ', N_det
|
||||
print *, 'N_sop = ', N_occ_pattern
|
||||
print *, 'N_states = ', N_states
|
||||
print*, 'correlation_ratio = ', correlation_energy_ratio
|
||||
|
||||
|
||||
do k=1,N_states
|
||||
rpt2(:) = pt2(:)/(1.d0 + norm(k))
|
||||
enddo
|
||||
|
||||
call print_summary(psi_energy_with_nucl_rep(1:N_states),pt2,error,variance,norm)
|
||||
call save_iterations(psi_energy_with_nucl_rep(1:N_states),rpt2,N_det)
|
||||
call print_extrapolated_energy(psi_energy_with_nucl_rep(1:N_states),rpt2)
|
||||
|
||||
end
|
||||
|
|
@ -1,7 +0,0 @@
|
|||
Perturbation
|
||||
Selectors_full
|
||||
Generators_full
|
||||
ZMQ
|
||||
MPI
|
||||
DavidsonUndressed
|
||||
Iterations
|
||||
|
|
@ -1,40 +0,0 @@
|
|||
program pt2_stoch
|
||||
implicit none
|
||||
read_wf = .True.
|
||||
SOFT_TOUCH read_wf
|
||||
PROVIDE mo_bielec_integrals_in_map
|
||||
PROVIDE psi_energy
|
||||
call run
|
||||
end
|
||||
|
||||
subroutine run
|
||||
implicit none
|
||||
integer :: i,j,k
|
||||
logical, external :: detEq
|
||||
|
||||
double precision :: pt2(N_states)
|
||||
integer :: degree
|
||||
integer :: n_det_before, to_select
|
||||
double precision :: threshold_davidson_in
|
||||
|
||||
double precision :: E_CI_before(N_states), relative_error, error(N_states), variance(N_states), norm(N_states), rpt2(N_states)
|
||||
|
||||
pt2(:) = 0.d0
|
||||
|
||||
E_CI_before(:) = psi_energy(:) + nuclear_repulsion
|
||||
threshold_selectors = 1.d0
|
||||
threshold_generators = 1.d0
|
||||
relative_error=PT2_relative_error
|
||||
|
||||
call ZMQ_pt2(psi_energy_with_nucl_rep,pt2,relative_error,error, variance, norm) ! Stochastic PT2
|
||||
do k=1,N_states
|
||||
rpt2(:) = pt2(:)/(1.d0 + norm(k))
|
||||
enddo
|
||||
|
||||
call print_summary(psi_energy_with_nucl_rep(1:N_states),pt2,error,variance,norm)
|
||||
|
||||
call ezfio_set_fci_energy(E_CI_before)
|
||||
call ezfio_set_fci_energy_pt2(E_CI_before+pt2)
|
||||
end
|
||||
|
||||
|
||||
|
|
@ -1,112 +0,0 @@
|
|||
===
|
||||
FCI
|
||||
===
|
||||
|
||||
Selected Full Configuration Interaction.
|
||||
|
||||
The :command:`FCI` program starts with a single determinant, or with the wave
|
||||
function in the |EZFIO| database if :option:`determinants read_wf` is |true|.
|
||||
Then, it will iteratively:
|
||||
|
||||
* Select the most important determinants from the external space and add them to the
|
||||
internal space
|
||||
* If :option:`determinants s2_eig` is |true|, add all the necessary
|
||||
determinants to allow the eigenstates of |H| to be eigenstates of |S^2|
|
||||
* Diagonalize |H| in the enlarged internal space
|
||||
* Compute (stochastically) the second-order perturbative contribution to the energy
|
||||
* Extrapolate the variational energy by fitting
|
||||
:math:`E=E_\text{FCI} - \alpha\, E_\text{PT2}`
|
||||
|
||||
|
||||
The number of selected determinants at each iteration will be such that the
|
||||
size of the wave function will double at every iteration. If :option:`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 :math:`10^{-4}`.
|
||||
|
||||
The variational and |PT2| energies of the iterations are stored in the
|
||||
|EZFIO| database, in the :ref:`IterativeSave` module.
|
||||
|
||||
|
||||
|
||||
Computation of the |PT2| energy
|
||||
-------------------------------
|
||||
|
||||
At each iteration, the |PT2| energy is computed considering the Epstein-Nesbet
|
||||
zeroth-order Hamiltonian:
|
||||
|
||||
.. math::
|
||||
|
||||
E_{\text{PT2}} = \sum_{ \alpha }
|
||||
\frac{|\langle \Psi_S | \hat{H} | \alpha \rangle|^2}
|
||||
{E - \langle \alpha | \hat{H} | \alpha \rangle}
|
||||
|
||||
where the |kalpha| determinants are generated by applying all the single and
|
||||
double excitation operators to all the determinants of the wave function
|
||||
:math:`\Psi_G`.
|
||||
|
||||
When the hybrid-deterministic/stochastic algorithm is chosen
|
||||
(default), :math:`Psi_G = \Psi_S = \Psi`, the full wavefunction expanded in the
|
||||
internal space.
|
||||
When the deterministic algorithm is chosen (:option:`perturbation do_pt2`
|
||||
is set to |false|), :math:`Psi_G` is a truncation of |Psi| using
|
||||
:option:`determinants threshold_generators`, and :math:`Psi_S` is a truncation
|
||||
of |Psi| using :option:`determinants threshold_selectors`, and re-weighted
|
||||
by :math:`1/\langle \Psi_s | \Psi_s \rangle`.
|
||||
|
||||
At every iteration, while computing the |PT2|, the variance of the wave
|
||||
function is also computed:
|
||||
|
||||
.. math::
|
||||
|
||||
\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.
|
||||
|
||||
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 :option:`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 :command:`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
|
||||
|
||||
.. math::
|
||||
|
||||
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.
|
||||
|
||||
|
|
@ -1,33 +0,0 @@
|
|||
BEGIN_PROVIDER [ logical, initialize_pt2_E0_denominator ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! If true, initialize pt2_E0_denominator
|
||||
END_DOC
|
||||
initialize_pt2_E0_denominator = .True.
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ double precision, pt2_E0_denominator, (N_states) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! E0 in the denominator of the PT2
|
||||
END_DOC
|
||||
if (initialize_pt2_E0_denominator) then
|
||||
if (h0_type == "EN") then
|
||||
pt2_E0_denominator(1:N_states) = psi_energy(1:N_states)
|
||||
else if (h0_type == "Barycentric") then
|
||||
pt2_E0_denominator(1:N_states) = barycentric_electronic_energy(1:N_states)
|
||||
else if (h0_type == "Variance") then
|
||||
pt2_E0_denominator(1:N_states) = psi_energy(1:N_states) !1.d0-nuclear_repulsion
|
||||
else if (h0_type == "SOP") then
|
||||
pt2_E0_denominator(1:N_states) = psi_energy(1:N_states)
|
||||
else
|
||||
print *, h0_type, ' not implemented'
|
||||
stop
|
||||
endif
|
||||
call write_double(6,pt2_E0_denominator(1)+nuclear_repulsion, 'PT2 Energy denominator')
|
||||
else
|
||||
pt2_E0_denominator = -huge(1.d0)
|
||||
endif
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
|
|
@ -1,669 +0,0 @@
|
|||
BEGIN_PROVIDER [ integer, pt2_stoch_istate ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! State for stochatsic PT2
|
||||
END_DOC
|
||||
pt2_stoch_istate = 1
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ integer, pt2_F, (N_det_generators) ]
|
||||
&BEGIN_PROVIDER [ integer, pt2_n_tasks_max ]
|
||||
implicit none
|
||||
logical, external :: testTeethBuilding
|
||||
integer :: i
|
||||
integer :: e
|
||||
e = elec_num - n_core_orb * 2
|
||||
pt2_n_tasks_max = 1+min((e*(e-1))/2, int(dsqrt(dble(N_det_selectors)))/10)
|
||||
do i=1,N_det_generators
|
||||
pt2_F(i) = 1 + int(dble(pt2_n_tasks_max)*maxval(dsqrt(dabs(psi_coef_sorted_gen(i,1:N_states)))))
|
||||
enddo
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ integer, pt2_N_teeth ]
|
||||
&BEGIN_PROVIDER [ integer, pt2_minDetInFirstTeeth ]
|
||||
implicit none
|
||||
logical, external :: testTeethBuilding
|
||||
|
||||
if(N_det_generators < 1024) then
|
||||
pt2_minDetInFirstTeeth = 1
|
||||
pt2_N_teeth = 1
|
||||
else
|
||||
pt2_minDetInFirstTeeth = min(5, N_det_generators)
|
||||
do pt2_N_teeth=100,2,-1
|
||||
if(testTeethBuilding(pt2_minDetInFirstTeeth, pt2_N_teeth)) exit
|
||||
end do
|
||||
end if
|
||||
call write_int(6,pt2_N_teeth,'Number of comb teeth')
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
logical function testTeethBuilding(minF, N)
|
||||
implicit none
|
||||
integer, intent(in) :: minF, N
|
||||
integer :: n0, i
|
||||
double precision :: u0, Wt, r
|
||||
|
||||
double precision, allocatable :: tilde_w(:), tilde_cW(:)
|
||||
integer, external :: dress_find_sample
|
||||
|
||||
double precision :: rss
|
||||
double precision, external :: memory_of_double, memory_of_int
|
||||
|
||||
rss = memory_of_double(2*N_det_generators+1)
|
||||
call check_mem(rss,irp_here)
|
||||
|
||||
allocate(tilde_w(N_det_generators), tilde_cW(0:N_det_generators))
|
||||
|
||||
do i=1,N_det_generators
|
||||
tilde_w(i) = psi_coef_sorted_gen(i,pt2_stoch_istate)**2 !+ 1.d-20
|
||||
enddo
|
||||
|
||||
double precision :: norm
|
||||
norm = 0.d0
|
||||
do i=N_det_generators,1,-1
|
||||
norm += tilde_w(i)
|
||||
enddo
|
||||
|
||||
tilde_w(:) = tilde_w(:) / norm
|
||||
|
||||
tilde_cW(0) = -1.d0
|
||||
do i=1,N_det_generators
|
||||
tilde_cW(i) = tilde_cW(i-1) + tilde_w(i)
|
||||
enddo
|
||||
tilde_cW(:) = tilde_cW(:) + 1.d0
|
||||
|
||||
n0 = 0
|
||||
testTeethBuilding = .false.
|
||||
do
|
||||
u0 = tilde_cW(n0)
|
||||
r = tilde_cW(n0 + minF)
|
||||
Wt = (1d0 - u0) / dble(N)
|
||||
if (dabs(Wt) <= 1.d-3) then
|
||||
return
|
||||
endif
|
||||
if(Wt >= r - u0) then
|
||||
testTeethBuilding = .true.
|
||||
return
|
||||
end if
|
||||
n0 += 1
|
||||
if(N_det_generators - n0 < minF * N) then
|
||||
return
|
||||
end if
|
||||
end do
|
||||
stop "exited testTeethBuilding"
|
||||
end function
|
||||
|
||||
|
||||
|
||||
subroutine ZMQ_pt2(E, pt2,relative_error, error, variance, norm)
|
||||
use f77_zmq
|
||||
use selection_types
|
||||
|
||||
implicit none
|
||||
|
||||
integer(ZMQ_PTR) :: zmq_to_qp_run_socket, zmq_socket_pull
|
||||
integer, external :: omp_get_thread_num
|
||||
double precision, intent(in) :: relative_error, E(N_states)
|
||||
double precision, intent(out) :: pt2(N_states),error(N_states)
|
||||
double precision, intent(out) :: variance(N_states),norm(N_states)
|
||||
|
||||
|
||||
integer :: i
|
||||
|
||||
double precision, external :: omp_get_wtime
|
||||
double precision :: state_average_weight_save(N_states), w(N_states,4)
|
||||
integer(ZMQ_PTR), external :: new_zmq_to_qp_run_socket
|
||||
|
||||
if (N_det < max(10,N_states)) then
|
||||
pt2=0.d0
|
||||
variance=0.d0
|
||||
norm=0.d0
|
||||
call ZMQ_selection(0, pt2, variance, norm)
|
||||
error(:) = 0.d0
|
||||
else
|
||||
|
||||
state_average_weight_save(:) = state_average_weight(:)
|
||||
do pt2_stoch_istate=1,N_states
|
||||
state_average_weight(:) = 0.d0
|
||||
state_average_weight(pt2_stoch_istate) = 1.d0
|
||||
TOUCH state_average_weight pt2_stoch_istate
|
||||
|
||||
provide nproc pt2_F mo_bielec_integrals_in_map mo_mono_elec_integral pt2_w psi_selectors
|
||||
call new_parallel_job(zmq_to_qp_run_socket, zmq_socket_pull, 'pt2')
|
||||
|
||||
integer, external :: zmq_put_psi
|
||||
integer, external :: zmq_put_N_det_generators
|
||||
integer, external :: zmq_put_N_det_selectors
|
||||
integer, external :: zmq_put_dvector
|
||||
integer, external :: zmq_put_ivector
|
||||
if (zmq_put_psi(zmq_to_qp_run_socket,1) == -1) then
|
||||
stop 'Unable to put psi on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_N_det_generators(zmq_to_qp_run_socket, 1) == -1) then
|
||||
stop 'Unable to put N_det_generators on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_N_det_selectors(zmq_to_qp_run_socket, 1) == -1) then
|
||||
stop 'Unable to put N_det_selectors on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_dvector(zmq_to_qp_run_socket,1,'energy',pt2_e0_denominator,size(pt2_e0_denominator)) == -1) then
|
||||
stop 'Unable to put energy on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_dvector(zmq_to_qp_run_socket,1,'state_average_weight',state_average_weight,N_states) == -1) then
|
||||
stop 'Unable to put state_average_weight on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_ivector(zmq_to_qp_run_socket,1,'pt2_stoch_istate',pt2_stoch_istate,1) == -1) then
|
||||
stop 'Unable to put pt2_stoch_istate on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_dvector(zmq_to_qp_run_socket,1,'threshold_selectors',threshold_selectors,1) == -1) then
|
||||
stop 'Unable to put threshold_selectors on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_dvector(zmq_to_qp_run_socket,1,'threshold_generators',threshold_generators,1) == -1) then
|
||||
stop 'Unable to put threshold_generators on ZMQ server'
|
||||
endif
|
||||
|
||||
|
||||
|
||||
integer, external :: add_task_to_taskserver
|
||||
character(400000) :: task
|
||||
|
||||
integer :: j,k,ipos,ifirst
|
||||
ifirst=0
|
||||
|
||||
ipos=0
|
||||
do i=1,N_det_generators
|
||||
if (pt2_F(i) > 1) then
|
||||
ipos += 1
|
||||
endif
|
||||
enddo
|
||||
call write_int(6,sum(pt2_F),'Number of tasks')
|
||||
call write_int(6,ipos,'Number of fragmented tasks')
|
||||
|
||||
ipos=1
|
||||
do i= 1, N_det_generators
|
||||
do j=1,pt2_F(pt2_J(i))
|
||||
write(task(ipos:ipos+20),'(I9,1X,I9,''|'')') j, pt2_J(i)
|
||||
ipos += 20
|
||||
if (ipos > 400000-20) then
|
||||
if (add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos))) == -1) then
|
||||
stop 'Unable to add task to task server'
|
||||
endif
|
||||
ipos=1
|
||||
if (ifirst == 0) then
|
||||
ifirst=1
|
||||
if (zmq_set_running(zmq_to_qp_run_socket) == -1) then
|
||||
print *, irp_here, ': Failed in zmq_set_running'
|
||||
endif
|
||||
endif
|
||||
endif
|
||||
end do
|
||||
enddo
|
||||
if (ipos > 1) then
|
||||
if (add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos))) == -1) then
|
||||
stop 'Unable to add task to task server'
|
||||
endif
|
||||
endif
|
||||
|
||||
integer, external :: zmq_set_running
|
||||
if (zmq_set_running(zmq_to_qp_run_socket) == -1) then
|
||||
print *, irp_here, ': Failed in zmq_set_running'
|
||||
endif
|
||||
|
||||
|
||||
integer :: nproc_target
|
||||
nproc_target = nproc
|
||||
double precision :: mem
|
||||
mem = 8.d0 * N_det * (N_int * 2.d0 * 3.d0 + 3.d0 + 5.d0) / (1024.d0**3)
|
||||
call write_double(6,mem,'Estimated memory/thread (Gb)')
|
||||
if (qp_max_mem > 0) then
|
||||
nproc_target = max(1,int(dble(qp_max_mem)/mem))
|
||||
nproc_target = min(nproc_target,nproc)
|
||||
endif
|
||||
|
||||
call omp_set_nested(.false.)
|
||||
|
||||
|
||||
print '(A)', '========== ================= =========== =============== =============== ================='
|
||||
print '(A)', ' Samples Energy Stat. Err Variance Norm Seconds '
|
||||
print '(A)', '========== ================= =========== =============== =============== ================='
|
||||
|
||||
!$OMP PARALLEL DEFAULT(shared) NUM_THREADS(nproc_target+1) &
|
||||
!$OMP PRIVATE(i)
|
||||
i = omp_get_thread_num()
|
||||
if (i==0) then
|
||||
|
||||
call pt2_collector(zmq_socket_pull, E(pt2_stoch_istate),relative_error, w(1,1), w(1,2), w(1,3), w(1,4))
|
||||
pt2(pt2_stoch_istate) = w(pt2_stoch_istate,1)
|
||||
error(pt2_stoch_istate) = w(pt2_stoch_istate,2)
|
||||
variance(pt2_stoch_istate) = w(pt2_stoch_istate,3)
|
||||
norm(pt2_stoch_istate) = w(pt2_stoch_istate,4)
|
||||
|
||||
else
|
||||
call pt2_slave_inproc(i)
|
||||
endif
|
||||
!$OMP END PARALLEL
|
||||
call end_parallel_job(zmq_to_qp_run_socket, zmq_socket_pull, 'pt2')
|
||||
|
||||
print '(A)', '========== ================= =========== =============== =============== ================='
|
||||
|
||||
enddo
|
||||
! call omp_set_nested(.false.)
|
||||
|
||||
FREE pt2_stoch_istate
|
||||
state_average_weight(:) = state_average_weight_save(:)
|
||||
TOUCH state_average_weight
|
||||
endif
|
||||
do k=N_det+1,N_states
|
||||
pt2(k) = 0.d0
|
||||
enddo
|
||||
|
||||
end subroutine
|
||||
|
||||
|
||||
subroutine pt2_slave_inproc(i)
|
||||
implicit none
|
||||
integer, intent(in) :: i
|
||||
|
||||
call run_pt2_slave(1,i,pt2_e0_denominator)
|
||||
end
|
||||
|
||||
|
||||
subroutine pt2_collector(zmq_socket_pull, E, relative_error, pt2, error, variance, norm)
|
||||
use f77_zmq
|
||||
use selection_types
|
||||
use bitmasks
|
||||
implicit none
|
||||
|
||||
|
||||
integer(ZMQ_PTR), intent(in) :: zmq_socket_pull
|
||||
double precision, intent(in) :: relative_error, E
|
||||
double precision, intent(out) :: pt2(N_states), error(N_states)
|
||||
double precision, intent(out) :: variance(N_states), norm(N_states)
|
||||
|
||||
|
||||
double precision, allocatable :: eI(:,:), eI_task(:,:), S(:), S2(:)
|
||||
double precision, allocatable :: vI(:,:), vI_task(:,:), T2(:)
|
||||
double precision, allocatable :: nI(:,:), nI_task(:,:), T3(:)
|
||||
integer(ZMQ_PTR),external :: new_zmq_to_qp_run_socket
|
||||
integer(ZMQ_PTR) :: zmq_to_qp_run_socket
|
||||
integer, external :: zmq_delete_tasks
|
||||
integer, external :: zmq_abort
|
||||
integer, external :: pt2_find_sample_lr
|
||||
|
||||
integer :: more, n, i, p, c, t, n_tasks, U
|
||||
integer, allocatable :: task_id(:)
|
||||
integer, allocatable :: index(:)
|
||||
|
||||
double precision, external :: omp_get_wtime
|
||||
double precision :: v, x, x2, x3, avg, avg2, avg3, eqt, E0, v0, n0
|
||||
double precision :: time, time0
|
||||
|
||||
integer, allocatable :: f(:)
|
||||
logical, allocatable :: d(:)
|
||||
logical :: do_exit
|
||||
|
||||
|
||||
double precision :: rss
|
||||
double precision, external :: memory_of_double, memory_of_int
|
||||
|
||||
rss = memory_of_int(pt2_n_tasks_max*2+N_det_generators*2)
|
||||
rss += memory_of_double(N_states*N_det_generators)*3.d0
|
||||
rss += memory_of_double(N_states*pt2_n_tasks_max)*3.d0
|
||||
rss += memory_of_double(pt2_N_teeth+1)*4.d0
|
||||
call check_mem(rss,irp_here)
|
||||
|
||||
allocate(task_id(pt2_n_tasks_max), index(pt2_n_tasks_max), f(N_det_generators))
|
||||
allocate(d(N_det_generators+1))
|
||||
allocate(eI(N_states, N_det_generators), eI_task(N_states, pt2_n_tasks_max))
|
||||
allocate(vI(N_states, N_det_generators), vI_task(N_states, pt2_n_tasks_max))
|
||||
allocate(nI(N_states, N_det_generators), nI_task(N_states, pt2_n_tasks_max))
|
||||
allocate(S(pt2_N_teeth+1), S2(pt2_N_teeth+1))
|
||||
allocate(T2(pt2_N_teeth+1), T3(pt2_N_teeth+1))
|
||||
|
||||
|
||||
|
||||
zmq_to_qp_run_socket = new_zmq_to_qp_run_socket()
|
||||
|
||||
pt2(:) = -huge(1.)
|
||||
error(:) = huge(1.)
|
||||
variance(:) = huge(1.)
|
||||
norm(:) = 0.d0
|
||||
S(:) = 0d0
|
||||
S2(:) = 0d0
|
||||
T2(:) = 0d0
|
||||
T3(:) = 0d0
|
||||
n = 1
|
||||
t = 0
|
||||
U = 0
|
||||
eI(:,:) = 0d0
|
||||
vI(:,:) = 0d0
|
||||
nI(:,:) = 0d0
|
||||
f(:) = pt2_F(:)
|
||||
d(:) = .false.
|
||||
n_tasks = 0
|
||||
E0 = E
|
||||
v0 = 0.d0
|
||||
n0 = 0.d0
|
||||
more = 1
|
||||
time0 = omp_get_wtime()
|
||||
|
||||
do_exit = .false.
|
||||
do while (n <= N_det_generators)
|
||||
if(f(pt2_J(n)) == 0) then
|
||||
d(pt2_J(n)) = .true.
|
||||
do while(d(U+1))
|
||||
U += 1
|
||||
end do
|
||||
|
||||
! Deterministic part
|
||||
do while(t <= pt2_N_teeth)
|
||||
if(U >= pt2_n_0(t+1)) then
|
||||
t=t+1
|
||||
E0 = 0.d0
|
||||
v0 = 0.d0
|
||||
n0 = 0.d0
|
||||
do i=pt2_n_0(t),1,-1
|
||||
E0 += eI(pt2_stoch_istate, i)
|
||||
v0 += vI(pt2_stoch_istate, i)
|
||||
n0 += nI(pt2_stoch_istate, i)
|
||||
end do
|
||||
else
|
||||
exit
|
||||
end if
|
||||
end do
|
||||
|
||||
! Add Stochastic part
|
||||
c = pt2_R(n)
|
||||
if(c > 0) then
|
||||
x = 0d0
|
||||
x2 = 0d0
|
||||
x3 = 0d0
|
||||
do p=pt2_N_teeth, 1, -1
|
||||
v = pt2_u_0 + pt2_W_T * (pt2_u(c) + dble(p-1))
|
||||
i = pt2_find_sample_lr(v, pt2_cW,pt2_n_0(p),pt2_n_0(p+1))
|
||||
x += eI(pt2_stoch_istate, i) * pt2_W_T / pt2_w(i)
|
||||
x2 += vI(pt2_stoch_istate, i) * pt2_W_T / pt2_w(i)
|
||||
x3 += nI(pt2_stoch_istate, i) * pt2_W_T / pt2_w(i)
|
||||
S(p) += x
|
||||
S2(p) += x*x
|
||||
T2(p) += x2
|
||||
T3(p) += x3
|
||||
end do
|
||||
avg = E0 + S(t) / dble(c)
|
||||
avg2 = v0 + T2(t) / dble(c)
|
||||
avg3 = n0 + T3(t) / dble(c)
|
||||
if ((avg /= 0.d0) .or. (n == N_det_generators) ) then
|
||||
do_exit = .true.
|
||||
endif
|
||||
pt2(pt2_stoch_istate) = avg
|
||||
variance(pt2_stoch_istate) = avg2
|
||||
norm(pt2_stoch_istate) = avg3
|
||||
! 1/(N-1.5) : see Brugger, The American Statistician (23) 4 p. 32 (1969)
|
||||
if(c > 2) then
|
||||
eqt = dabs((S2(t) / c) - (S(t)/c)**2) ! dabs for numerical stability
|
||||
eqt = sqrt(eqt / (dble(c) - 1.5d0))
|
||||
error(pt2_stoch_istate) = eqt
|
||||
if(mod(c,10)==0 .or. n==N_det_generators) then
|
||||
print '(G10.3, 2X, F16.10, 2X, G10.3, 2X, F14.10, 2X, F14.10, 2X, F10.4, A10)', c, avg+E, eqt, avg2, avg3, time-time0, ''
|
||||
if(do_exit .and. (dabs(error(pt2_stoch_istate)) / (1.d-20 + dabs(pt2(pt2_stoch_istate)) ) <= relative_error)) then
|
||||
if (zmq_abort(zmq_to_qp_run_socket) == -1) then
|
||||
call sleep(10)
|
||||
if (zmq_abort(zmq_to_qp_run_socket) == -1) then
|
||||
print *, irp_here, ': Error in sending abort signal (2)'
|
||||
endif
|
||||
endif
|
||||
endif
|
||||
endif
|
||||
endif
|
||||
time = omp_get_wtime()
|
||||
end if
|
||||
n += 1
|
||||
else if(more == 0) then
|
||||
exit
|
||||
else
|
||||
call pull_pt2_results(zmq_socket_pull, index, eI_task, vI_task, nI_task, task_id, n_tasks)
|
||||
if (zmq_delete_tasks(zmq_to_qp_run_socket,zmq_socket_pull,task_id,n_tasks,more) == -1) then
|
||||
stop 'Unable to delete tasks'
|
||||
endif
|
||||
do i=1,n_tasks
|
||||
eI(:, index(i)) += eI_task(:,i)
|
||||
vI(:, index(i)) += vI_task(:,i)
|
||||
nI(:, index(i)) += nI_task(:,i)
|
||||
f(index(i)) -= 1
|
||||
end do
|
||||
end if
|
||||
end do
|
||||
call end_zmq_to_qp_run_socket(zmq_to_qp_run_socket)
|
||||
end subroutine
|
||||
|
||||
|
||||
integer function pt2_find_sample(v, w)
|
||||
implicit none
|
||||
double precision, intent(in) :: v, w(0:N_det_generators)
|
||||
integer, external :: pt2_find_sample_lr
|
||||
|
||||
pt2_find_sample = pt2_find_sample_lr(v, w, 0, N_det_generators)
|
||||
end function
|
||||
|
||||
|
||||
integer function pt2_find_sample_lr(v, w, l_in, r_in)
|
||||
implicit none
|
||||
double precision, intent(in) :: v, w(0:N_det_generators)
|
||||
integer, intent(in) :: l_in,r_in
|
||||
integer :: i,l,r
|
||||
|
||||
l=l_in
|
||||
r=r_in
|
||||
|
||||
do while(r-l > 1)
|
||||
i = shiftr(r+l,1)
|
||||
if(w(i) < v) then
|
||||
l = i
|
||||
else
|
||||
r = i
|
||||
end if
|
||||
end do
|
||||
i = r
|
||||
do r=i+1,N_det_generators
|
||||
if (w(r) /= w(i)) then
|
||||
exit
|
||||
endif
|
||||
enddo
|
||||
pt2_find_sample_lr = r-1
|
||||
end function
|
||||
|
||||
|
||||
BEGIN_PROVIDER [ integer, pt2_n_tasks ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Number of parallel tasks for the Monte Carlo
|
||||
END_DOC
|
||||
pt2_n_tasks = N_det_generators
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER[ double precision, pt2_u, (N_det_generators)]
|
||||
implicit none
|
||||
integer, allocatable :: seed(:)
|
||||
integer :: m,i
|
||||
call random_seed(size=m)
|
||||
allocate(seed(m))
|
||||
do i=1,m
|
||||
seed(i) = i
|
||||
enddo
|
||||
call random_seed(put=seed)
|
||||
deallocate(seed)
|
||||
|
||||
call RANDOM_NUMBER(pt2_u)
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER[ integer, pt2_J, (N_det_generators)]
|
||||
&BEGIN_PROVIDER[ integer, pt2_R, (N_det_generators)]
|
||||
implicit none
|
||||
integer :: N_c, N_j
|
||||
integer :: U, t, i
|
||||
double precision :: v
|
||||
integer, external :: pt2_find_sample_lr
|
||||
|
||||
logical, allocatable :: pt2_d(:)
|
||||
integer :: m,l,r,k
|
||||
integer :: ncache
|
||||
integer, allocatable :: ii(:,:)
|
||||
double precision :: dt
|
||||
|
||||
ncache = min(N_det_generators,10000)
|
||||
|
||||
double precision :: rss
|
||||
double precision, external :: memory_of_double, memory_of_int
|
||||
rss = memory_of_int(ncache)*dble(pt2_N_teeth) + memory_of_int(N_det_generators)
|
||||
call check_mem(rss,irp_here)
|
||||
|
||||
allocate(ii(pt2_N_teeth,ncache),pt2_d(N_det_generators))
|
||||
|
||||
pt2_R(:) = 0
|
||||
pt2_d(:) = .false.
|
||||
N_c = 0
|
||||
N_j = pt2_n_0(1)
|
||||
do i=1,N_j
|
||||
pt2_d(i) = .true.
|
||||
pt2_J(i) = i
|
||||
end do
|
||||
|
||||
U = 0
|
||||
do while(N_j < pt2_n_tasks)
|
||||
|
||||
if (N_c+ncache > N_det_generators) then
|
||||
ncache = N_det_generators - N_c
|
||||
endif
|
||||
|
||||
!$OMP PARALLEL DO DEFAULT(SHARED) PRIVATE(dt,v,t,k)
|
||||
do k=1, ncache
|
||||
dt = pt2_u_0
|
||||
do t=1, pt2_N_teeth
|
||||
v = dt + pt2_W_T *pt2_u(N_c+k)
|
||||
dt = dt + pt2_W_T
|
||||
ii(t,k) = pt2_find_sample_lr(v, pt2_cW,pt2_n_0(t),pt2_n_0(t+1))
|
||||
end do
|
||||
enddo
|
||||
!$OMP END PARALLEL DO
|
||||
|
||||
do k=1,ncache
|
||||
!ADD_COMB
|
||||
N_c = N_c+1
|
||||
do t=1, pt2_N_teeth
|
||||
i = ii(t,k)
|
||||
if(.not. pt2_d(i)) then
|
||||
N_j += 1
|
||||
pt2_J(N_j) = i
|
||||
pt2_d(i) = .true.
|
||||
end if
|
||||
end do
|
||||
|
||||
pt2_R(N_j) = N_c
|
||||
|
||||
!FILL_TOOTH
|
||||
do while(U < N_det_generators)
|
||||
U += 1
|
||||
if(.not. pt2_d(U)) then
|
||||
N_j += 1
|
||||
pt2_J(N_j) = U
|
||||
pt2_d(U) = .true.
|
||||
exit
|
||||
end if
|
||||
end do
|
||||
if (N_j >= pt2_n_tasks) exit
|
||||
end do
|
||||
enddo
|
||||
|
||||
if(N_det_generators > 1) then
|
||||
pt2_R(N_det_generators-1) = 0
|
||||
pt2_R(N_det_generators) = N_c
|
||||
end if
|
||||
|
||||
deallocate(ii,pt2_d)
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
|
||||
BEGIN_PROVIDER [ double precision, pt2_w, (N_det_generators) ]
|
||||
&BEGIN_PROVIDER [ double precision, pt2_cW, (0:N_det_generators) ]
|
||||
&BEGIN_PROVIDER [ double precision, pt2_W_T ]
|
||||
&BEGIN_PROVIDER [ double precision, pt2_u_0 ]
|
||||
&BEGIN_PROVIDER [ integer, pt2_n_0, (pt2_N_teeth+1) ]
|
||||
implicit none
|
||||
integer :: i, t
|
||||
double precision, allocatable :: tilde_w(:), tilde_cW(:)
|
||||
double precision :: r, tooth_width
|
||||
integer, external :: pt2_find_sample
|
||||
|
||||
double precision :: rss
|
||||
double precision, external :: memory_of_double, memory_of_int
|
||||
rss = memory_of_double(2*N_det_generators+1)
|
||||
call check_mem(rss,irp_here)
|
||||
|
||||
allocate(tilde_w(N_det_generators), tilde_cW(0:N_det_generators))
|
||||
|
||||
tilde_cW(0) = 0d0
|
||||
|
||||
do i=1,N_det_generators
|
||||
tilde_w(i) = psi_coef_sorted_gen(i,pt2_stoch_istate)**2 !+ 1.d-20
|
||||
enddo
|
||||
|
||||
double precision :: norm
|
||||
norm = 0.d0
|
||||
do i=N_det_generators,1,-1
|
||||
norm += tilde_w(i)
|
||||
enddo
|
||||
|
||||
tilde_w(:) = tilde_w(:) / norm
|
||||
|
||||
tilde_cW(0) = -1.d0
|
||||
do i=1,N_det_generators
|
||||
tilde_cW(i) = tilde_cW(i-1) + tilde_w(i)
|
||||
enddo
|
||||
tilde_cW(:) = tilde_cW(:) + 1.d0
|
||||
|
||||
pt2_n_0(1) = 0
|
||||
do
|
||||
pt2_u_0 = tilde_cW(pt2_n_0(1))
|
||||
r = tilde_cW(pt2_n_0(1) + pt2_minDetInFirstTeeth)
|
||||
pt2_W_T = (1d0 - pt2_u_0) / dble(pt2_N_teeth)
|
||||
if(pt2_W_T >= r - pt2_u_0) then
|
||||
exit
|
||||
end if
|
||||
pt2_n_0(1) += 1
|
||||
if(N_det_generators - pt2_n_0(1) < pt2_minDetInFirstTeeth * pt2_N_teeth) then
|
||||
stop "teeth building failed"
|
||||
end if
|
||||
end do
|
||||
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
||||
|
||||
do t=2, pt2_N_teeth
|
||||
r = pt2_u_0 + pt2_W_T * dble(t-1)
|
||||
pt2_n_0(t) = pt2_find_sample(r, tilde_cW)
|
||||
end do
|
||||
pt2_n_0(pt2_N_teeth+1) = N_det_generators
|
||||
|
||||
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
|
||||
pt2_w(:pt2_n_0(1)) = tilde_w(:pt2_n_0(1))
|
||||
do t=1, pt2_N_teeth
|
||||
tooth_width = tilde_cW(pt2_n_0(t+1)) - tilde_cW(pt2_n_0(t))
|
||||
if (tooth_width == 0.d0) then
|
||||
tooth_width = sum(tilde_w(pt2_n_0(t):pt2_n_0(t+1)))
|
||||
endif
|
||||
ASSERT(tooth_width > 0.d0)
|
||||
do i=pt2_n_0(t)+1, pt2_n_0(t+1)
|
||||
pt2_w(i) = tilde_w(i) * pt2_W_T / tooth_width
|
||||
end do
|
||||
end do
|
||||
|
||||
pt2_cW(0) = 0d0
|
||||
do i=1,N_det_generators
|
||||
pt2_cW(i) = pt2_cW(i-1) + pt2_w(i)
|
||||
end do
|
||||
pt2_n_0(pt2_N_teeth+1) = N_det_generators
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
|
|
@ -1,247 +0,0 @@
|
|||
|
||||
subroutine run_pt2_slave(thread,iproc,energy)
|
||||
use f77_zmq
|
||||
use selection_types
|
||||
implicit none
|
||||
|
||||
double precision, intent(in) :: energy(N_states_diag)
|
||||
integer, intent(in) :: thread, iproc
|
||||
integer :: rc, i
|
||||
|
||||
integer :: worker_id, ctask, ltask
|
||||
character*(512), allocatable :: task(:)
|
||||
integer, allocatable :: task_id(:)
|
||||
|
||||
integer(ZMQ_PTR),external :: new_zmq_to_qp_run_socket
|
||||
integer(ZMQ_PTR) :: zmq_to_qp_run_socket
|
||||
|
||||
integer(ZMQ_PTR), external :: new_zmq_push_socket
|
||||
integer(ZMQ_PTR) :: zmq_socket_push
|
||||
|
||||
type(selection_buffer) :: buf
|
||||
logical :: done
|
||||
|
||||
double precision,allocatable :: pt2(:,:), variance(:,:), norm(:,:)
|
||||
integer :: n_tasks, k
|
||||
integer, allocatable :: i_generator(:), subset(:)
|
||||
|
||||
double precision :: rss
|
||||
double precision, external :: memory_of_double, memory_of_int
|
||||
rss = memory_of_int(pt2_n_tasks_max)*67.d0
|
||||
rss += memory_of_double(pt2_n_tasks_max)*(N_states*3)
|
||||
call check_mem(rss,irp_here)
|
||||
|
||||
allocate(task_id(pt2_n_tasks_max), task(pt2_n_tasks_max))
|
||||
allocate(pt2(N_states,pt2_n_tasks_max), i_generator(pt2_n_tasks_max), subset(pt2_n_tasks_max))
|
||||
allocate(variance(N_states,pt2_n_tasks_max))
|
||||
allocate(norm(N_states,pt2_n_tasks_max))
|
||||
|
||||
zmq_to_qp_run_socket = new_zmq_to_qp_run_socket()
|
||||
|
||||
integer, external :: connect_to_taskserver
|
||||
if (connect_to_taskserver(zmq_to_qp_run_socket,worker_id,thread) == -1) then
|
||||
call end_zmq_to_qp_run_socket(zmq_to_qp_run_socket)
|
||||
return
|
||||
endif
|
||||
|
||||
zmq_socket_push = new_zmq_push_socket(thread)
|
||||
|
||||
buf%N = 0
|
||||
n_tasks = 1
|
||||
call create_selection_buffer(0, 0, buf)
|
||||
|
||||
done = .False.
|
||||
n_tasks = 1
|
||||
do while (.not.done)
|
||||
|
||||
n_tasks = max(1,n_tasks)
|
||||
n_tasks = min(pt2_n_tasks_max,n_tasks)
|
||||
|
||||
integer, external :: get_tasks_from_taskserver
|
||||
if (get_tasks_from_taskserver(zmq_to_qp_run_socket,worker_id, task_id, task, n_tasks) == -1) then
|
||||
exit
|
||||
endif
|
||||
done = task_id(n_tasks) == 0
|
||||
if (done) n_tasks = n_tasks-1
|
||||
if (n_tasks == 0) exit
|
||||
|
||||
do k=1,n_tasks
|
||||
read (task(k),*) subset(k), i_generator(k)
|
||||
enddo
|
||||
|
||||
double precision :: time0, time1
|
||||
call wall_time(time0)
|
||||
do k=1,n_tasks
|
||||
pt2(:,k) = 0.d0
|
||||
variance(:,k) = 0.d0
|
||||
norm(:,k) = 0.d0
|
||||
buf%cur = 0
|
||||
!double precision :: time2
|
||||
!call wall_time(time2)
|
||||
call select_connected(i_generator(k),energy,pt2(1,k),variance(1,k),norm(1,k),buf,subset(k),pt2_F(i_generator(k)))
|
||||
!call wall_time(time1)
|
||||
!print *, i_generator(1), time1-time2, n_tasks, pt2_F(i_generator(1))
|
||||
enddo
|
||||
call wall_time(time1)
|
||||
!print *, i_generator(1), time1-time0, n_tasks
|
||||
|
||||
integer, external :: tasks_done_to_taskserver
|
||||
if (tasks_done_to_taskserver(zmq_to_qp_run_socket,worker_id,task_id,n_tasks) == -1) then
|
||||
done = .true.
|
||||
endif
|
||||
call push_pt2_results(zmq_socket_push, i_generator, pt2, variance, norm, task_id, n_tasks)
|
||||
|
||||
! Try to adjust n_tasks around nproc/8 seconds per job
|
||||
n_tasks = min(2*n_tasks,int( dble(n_tasks * nproc/8) / (time1 - time0 + 1.d0)))
|
||||
end do
|
||||
|
||||
integer, external :: disconnect_from_taskserver
|
||||
do i=1,300
|
||||
if (disconnect_from_taskserver(zmq_to_qp_run_socket,worker_id) /= -2) exit
|
||||
call sleep(1)
|
||||
print *, 'Retry disconnect...'
|
||||
end do
|
||||
|
||||
call end_zmq_push_socket(zmq_socket_push,thread)
|
||||
call end_zmq_to_qp_run_socket(zmq_to_qp_run_socket)
|
||||
call delete_selection_buffer(buf)
|
||||
end subroutine
|
||||
|
||||
|
||||
subroutine push_pt2_results(zmq_socket_push, index, pt2, variance, norm, task_id, n_tasks)
|
||||
use f77_zmq
|
||||
use selection_types
|
||||
implicit none
|
||||
|
||||
integer(ZMQ_PTR), intent(in) :: zmq_socket_push
|
||||
double precision, intent(in) :: pt2(N_states,n_tasks)
|
||||
double precision, intent(in) :: variance(N_states,n_tasks)
|
||||
double precision, intent(in) :: norm(N_states,n_tasks)
|
||||
integer, intent(in) :: n_tasks, index(n_tasks), task_id(n_tasks)
|
||||
integer :: rc
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, n_tasks, 4, ZMQ_SNDMORE)
|
||||
if (rc == -1) then
|
||||
return
|
||||
endif
|
||||
if(rc /= 4) stop 'push'
|
||||
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, index, 4*n_tasks, ZMQ_SNDMORE)
|
||||
if (rc == -1) then
|
||||
return
|
||||
endif
|
||||
if(rc /= 4*n_tasks) stop 'push'
|
||||
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, pt2, 8*N_states*n_tasks, ZMQ_SNDMORE)
|
||||
if (rc == -1) then
|
||||
return
|
||||
endif
|
||||
if(rc /= 8*N_states*n_tasks) stop 'push'
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, variance, 8*N_states*n_tasks, ZMQ_SNDMORE)
|
||||
if (rc == -1) then
|
||||
return
|
||||
endif
|
||||
if(rc /= 8*N_states*n_tasks) stop 'push'
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, norm, 8*N_states*n_tasks, ZMQ_SNDMORE)
|
||||
if (rc == -1) then
|
||||
return
|
||||
endif
|
||||
if(rc /= 8*N_states*n_tasks) stop 'push'
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, task_id, n_tasks*4, 0)
|
||||
if (rc == -1) then
|
||||
return
|
||||
endif
|
||||
if(rc /= 4*n_tasks) stop 'push'
|
||||
|
||||
! Activate is zmq_socket_push is a REQ
|
||||
IRP_IF ZMQ_PUSH
|
||||
IRP_ELSE
|
||||
character*(2) :: ok
|
||||
rc = f77_zmq_recv( zmq_socket_push, ok, 2, 0)
|
||||
if (rc == -1) then
|
||||
return
|
||||
endif
|
||||
if ((rc /= 2).and.(ok(1:2) /= 'ok')) then
|
||||
print *, irp_here//': error in receiving ok'
|
||||
stop -1
|
||||
endif
|
||||
IRP_ENDIF
|
||||
|
||||
end subroutine
|
||||
|
||||
|
||||
subroutine pull_pt2_results(zmq_socket_pull, index, pt2, variance, norm, task_id, n_tasks)
|
||||
use f77_zmq
|
||||
use selection_types
|
||||
implicit none
|
||||
integer(ZMQ_PTR), intent(in) :: zmq_socket_pull
|
||||
double precision, intent(inout) :: pt2(N_states,*)
|
||||
double precision, intent(inout) :: variance(N_states,*)
|
||||
double precision, intent(inout) :: norm(N_states,*)
|
||||
integer, intent(out) :: index(*)
|
||||
integer, intent(out) :: n_tasks, task_id(*)
|
||||
integer :: rc, rn, i
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, n_tasks, 4, 0)
|
||||
if (rc == -1) then
|
||||
n_tasks = 1
|
||||
task_id(1) = 0
|
||||
endif
|
||||
if(rc /= 4) stop 'pull'
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, index, 4*n_tasks, 0)
|
||||
if (rc == -1) then
|
||||
n_tasks = 1
|
||||
task_id(1) = 0
|
||||
endif
|
||||
if(rc /= 4*n_tasks) stop 'pull'
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, pt2, N_states*8*n_tasks, 0)
|
||||
if (rc == -1) then
|
||||
n_tasks = 1
|
||||
task_id(1) = 0
|
||||
endif
|
||||
if(rc /= 8*N_states*n_tasks) stop 'pull'
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, variance, N_states*8*n_tasks, 0)
|
||||
if (rc == -1) then
|
||||
n_tasks = 1
|
||||
task_id(1) = 0
|
||||
endif
|
||||
if(rc /= 8*N_states*n_tasks) stop 'pull'
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, norm, N_states*8*n_tasks, 0)
|
||||
if (rc == -1) then
|
||||
n_tasks = 1
|
||||
task_id(1) = 0
|
||||
endif
|
||||
if(rc /= 8*N_states*n_tasks) stop 'pull'
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, task_id, n_tasks*4, 0)
|
||||
if (rc == -1) then
|
||||
n_tasks = 1
|
||||
task_id(1) = 0
|
||||
endif
|
||||
if(rc /= 4*n_tasks) stop 'pull'
|
||||
|
||||
! Activate is zmq_socket_pull is a REP
|
||||
IRP_IF ZMQ_PUSH
|
||||
IRP_ELSE
|
||||
rc = f77_zmq_send( zmq_socket_pull, 'ok', 2, 0)
|
||||
if (rc == -1) then
|
||||
n_tasks = 1
|
||||
task_id(1) = 0
|
||||
endif
|
||||
if (rc /= 2) then
|
||||
print *, irp_here//': error in sending ok'
|
||||
stop -1
|
||||
endif
|
||||
IRP_ENDIF
|
||||
|
||||
end subroutine
|
||||
|
||||
|
|
@ -1,254 +0,0 @@
|
|||
subroutine run_selection_slave(thread,iproc,energy)
|
||||
use f77_zmq
|
||||
use selection_types
|
||||
implicit none
|
||||
|
||||
double precision, intent(in) :: energy(N_states)
|
||||
integer, intent(in) :: thread, iproc
|
||||
integer :: rc, i
|
||||
|
||||
integer :: worker_id, task_id(1), ctask, ltask
|
||||
character*(512) :: task
|
||||
|
||||
integer(ZMQ_PTR),external :: new_zmq_to_qp_run_socket
|
||||
integer(ZMQ_PTR) :: zmq_to_qp_run_socket
|
||||
|
||||
integer(ZMQ_PTR), external :: new_zmq_push_socket
|
||||
integer(ZMQ_PTR) :: zmq_socket_push
|
||||
|
||||
type(selection_buffer) :: buf, buf2
|
||||
logical :: done, buffer_ready
|
||||
double precision :: pt2(N_states)
|
||||
double precision :: variance(N_states)
|
||||
double precision :: norm(N_states)
|
||||
|
||||
PROVIDE psi_bilinear_matrix_columns_loc psi_det_alpha_unique psi_det_beta_unique
|
||||
PROVIDE psi_bilinear_matrix_rows psi_det_sorted_order psi_bilinear_matrix_order
|
||||
PROVIDE psi_bilinear_matrix_transp_rows_loc psi_bilinear_matrix_transp_columns
|
||||
PROVIDE psi_bilinear_matrix_transp_order N_int pt2_F
|
||||
|
||||
zmq_to_qp_run_socket = new_zmq_to_qp_run_socket()
|
||||
|
||||
integer, external :: connect_to_taskserver
|
||||
if (connect_to_taskserver(zmq_to_qp_run_socket,worker_id,thread) == -1) then
|
||||
call end_zmq_to_qp_run_socket(zmq_to_qp_run_socket)
|
||||
return
|
||||
endif
|
||||
|
||||
zmq_socket_push = new_zmq_push_socket(thread)
|
||||
|
||||
buf%N = 0
|
||||
buffer_ready = .False.
|
||||
ctask = 1
|
||||
pt2(:) = 0d0
|
||||
variance(:) = 0d0
|
||||
norm(:) = 0.d0
|
||||
|
||||
do
|
||||
integer, external :: get_task_from_taskserver
|
||||
if (get_task_from_taskserver(zmq_to_qp_run_socket,worker_id, task_id(ctask), task) == -1) then
|
||||
exit
|
||||
endif
|
||||
done = task_id(ctask) == 0
|
||||
if (done) then
|
||||
ctask = ctask - 1
|
||||
else
|
||||
integer :: i_generator, N, subset
|
||||
read(task,*) subset, i_generator, N
|
||||
if(buf%N == 0) then
|
||||
! Only first time
|
||||
call create_selection_buffer(N, N*2, buf)
|
||||
call create_selection_buffer(N, N*2, buf2)
|
||||
buffer_ready = .True.
|
||||
else
|
||||
ASSERT (N == buf%N)
|
||||
end if
|
||||
call select_connected(i_generator,energy,pt2,variance,norm,buf,subset,pt2_F(i_generator))
|
||||
endif
|
||||
|
||||
integer, external :: task_done_to_taskserver
|
||||
|
||||
if(done .or. ctask == size(task_id)) then
|
||||
do i=1, ctask
|
||||
if (task_done_to_taskserver(zmq_to_qp_run_socket,worker_id,task_id(i)) == -1) then
|
||||
call sleep(1)
|
||||
if (task_done_to_taskserver(zmq_to_qp_run_socket,worker_id,task_id(i)) == -1) then
|
||||
ctask = 0
|
||||
done = .true.
|
||||
exit
|
||||
endif
|
||||
endif
|
||||
end do
|
||||
if(ctask > 0) then
|
||||
call sort_selection_buffer(buf)
|
||||
call merge_selection_buffers(buf,buf2)
|
||||
call push_selection_results(zmq_socket_push, pt2, variance, norm, buf, task_id(1), ctask)
|
||||
buf%mini = buf2%mini
|
||||
pt2(:) = 0d0
|
||||
variance(:) = 0d0
|
||||
norm(:) = 0d0
|
||||
buf%cur = 0
|
||||
end if
|
||||
ctask = 0
|
||||
end if
|
||||
|
||||
if(done) exit
|
||||
ctask = ctask + 1
|
||||
end do
|
||||
|
||||
|
||||
integer, external :: disconnect_from_taskserver
|
||||
if (disconnect_from_taskserver(zmq_to_qp_run_socket,worker_id) == -1) then
|
||||
continue
|
||||
endif
|
||||
|
||||
call end_zmq_to_qp_run_socket(zmq_to_qp_run_socket)
|
||||
call end_zmq_push_socket(zmq_socket_push,thread)
|
||||
if (buffer_ready) then
|
||||
call delete_selection_buffer(buf)
|
||||
call delete_selection_buffer(buf2)
|
||||
endif
|
||||
end subroutine
|
||||
|
||||
|
||||
subroutine push_selection_results(zmq_socket_push, pt2, variance, norm, b, task_id, ntask)
|
||||
use f77_zmq
|
||||
use selection_types
|
||||
implicit none
|
||||
|
||||
integer(ZMQ_PTR), intent(in) :: zmq_socket_push
|
||||
double precision, intent(in) :: pt2(N_states)
|
||||
double precision, intent(in) :: variance(N_states)
|
||||
double precision, intent(in) :: norm(N_states)
|
||||
type(selection_buffer), intent(inout) :: b
|
||||
integer, intent(in) :: ntask, task_id(*)
|
||||
integer :: rc
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, b%cur, 4, ZMQ_SNDMORE)
|
||||
if(rc /= 4) then
|
||||
print *, 'f77_zmq_send( zmq_socket_push, b%cur, 4, ZMQ_SNDMORE)'
|
||||
endif
|
||||
|
||||
if (b%cur > 0) then
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, pt2, 8*N_states, ZMQ_SNDMORE)
|
||||
if(rc /= 8*N_states) then
|
||||
print *, 'f77_zmq_send( zmq_socket_push, pt2, 8*N_states, ZMQ_SNDMORE)'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, variance, 8*N_states, ZMQ_SNDMORE)
|
||||
if(rc /= 8*N_states) then
|
||||
print *, 'f77_zmq_send( zmq_socket_push, variance, 8*N_states, ZMQ_SNDMORE)'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, norm, 8*N_states, ZMQ_SNDMORE)
|
||||
if(rc /= 8*N_states) then
|
||||
print *, 'f77_zmq_send( zmq_socket_push, norm, 8*N_states, ZMQ_SNDMORE)'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, b%val(1), 8*b%cur, ZMQ_SNDMORE)
|
||||
if(rc /= 8*b%cur) then
|
||||
print *, 'f77_zmq_send( zmq_socket_push, b%val(1), 8*b%cur, ZMQ_SNDMORE)'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, b%det(1,1,1), bit_kind*N_int*2*b%cur, ZMQ_SNDMORE)
|
||||
if(rc /= bit_kind*N_int*2*b%cur) then
|
||||
print *, 'f77_zmq_send( zmq_socket_push, b%det(1,1,1), bit_kind*N_int*2*b%cur, ZMQ_SNDMORE)'
|
||||
endif
|
||||
|
||||
endif
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, ntask, 4, ZMQ_SNDMORE)
|
||||
if(rc /= 4) then
|
||||
print *, 'f77_zmq_send( zmq_socket_push, ntask, 4, ZMQ_SNDMORE)'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, task_id(1), ntask*4, 0)
|
||||
if(rc /= 4*ntask) then
|
||||
print *, 'f77_zmq_send( zmq_socket_push, task_id(1), ntask*4, 0)'
|
||||
endif
|
||||
|
||||
! Activate is zmq_socket_push is a REQ
|
||||
IRP_IF ZMQ_PUSH
|
||||
IRP_ELSE
|
||||
character*(2) :: ok
|
||||
rc = f77_zmq_recv( zmq_socket_push, ok, 2, 0)
|
||||
if ((rc /= 2).and.(ok(1:2) /= 'ok')) then
|
||||
print *, irp_here//': error in receiving ok'
|
||||
stop -1
|
||||
endif
|
||||
IRP_ENDIF
|
||||
|
||||
end subroutine
|
||||
|
||||
|
||||
subroutine pull_selection_results(zmq_socket_pull, pt2, variance, norm, val, det, N, task_id, ntask)
|
||||
use f77_zmq
|
||||
use selection_types
|
||||
implicit none
|
||||
integer(ZMQ_PTR), intent(in) :: zmq_socket_pull
|
||||
double precision, intent(inout) :: pt2(N_states)
|
||||
double precision, intent(inout) :: variance(N_states)
|
||||
double precision, intent(inout) :: norm(N_states)
|
||||
double precision, intent(out) :: val(*)
|
||||
integer(bit_kind), intent(out) :: det(N_int, 2, *)
|
||||
integer, intent(out) :: N, ntask, task_id(*)
|
||||
integer :: rc, rn, i
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, N, 4, 0)
|
||||
if(rc /= 4) then
|
||||
print *, 'f77_zmq_recv( zmq_socket_pull, N, 4, 0)'
|
||||
endif
|
||||
|
||||
if (N>0) then
|
||||
rc = f77_zmq_recv( zmq_socket_pull, pt2, N_states*8, 0)
|
||||
if(rc /= 8*N_states) then
|
||||
print *, 'f77_zmq_recv( zmq_socket_pull, pt2, N_states*8, 0)'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, variance, N_states*8, 0)
|
||||
if(rc /= 8*N_states) then
|
||||
print *, 'f77_zmq_recv( zmq_socket_pull, variance, N_states*8, 0)'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, norm, N_states*8, 0)
|
||||
if(rc /= 8*N_states) then
|
||||
print *, 'f77_zmq_recv( zmq_socket_pull, norm, N_states*8, 0)'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, val(1), 8*N, 0)
|
||||
if(rc /= 8*N) then
|
||||
print *, 'f77_zmq_recv( zmq_socket_pull, val(1), 8*N, 0)'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, det(1,1,1), bit_kind*N_int*2*N, 0)
|
||||
if(rc /= bit_kind*N_int*2*N) then
|
||||
print *, 'f77_zmq_recv( zmq_socket_pull, det(1,1,1), bit_kind*N_int*2*N, 0)'
|
||||
endif
|
||||
else
|
||||
pt2(:) = 0.d0
|
||||
endif
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, ntask, 4, 0)
|
||||
if(rc /= 4) then
|
||||
print *, 'f77_zmq_recv( zmq_socket_pull, ntask, 4, 0)'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, task_id(1), ntask*4, 0)
|
||||
if(rc /= 4*ntask) then
|
||||
print *, 'f77_zmq_recv( zmq_socket_pull, task_id(1), ntask*4, 0)'
|
||||
endif
|
||||
|
||||
! Activate is zmq_socket_pull is a REP
|
||||
IRP_IF ZMQ_PUSH
|
||||
IRP_ELSE
|
||||
rc = f77_zmq_send( zmq_socket_pull, 'ok', 2, 0)
|
||||
if (rc /= 2) then
|
||||
print *, irp_here//': error in sending ok'
|
||||
stop -1
|
||||
endif
|
||||
IRP_ENDIF
|
||||
end subroutine
|
||||
|
||||
|
||||
|
||||
File diff suppressed because it is too large
Load Diff
|
|
@ -1,303 +0,0 @@
|
|||
|
||||
subroutine create_selection_buffer(N, siz_, res)
|
||||
use selection_types
|
||||
implicit none
|
||||
|
||||
integer, intent(in) :: N, siz_
|
||||
type(selection_buffer), intent(out) :: res
|
||||
|
||||
integer :: siz
|
||||
siz = max(siz_,1)
|
||||
|
||||
double precision :: rss
|
||||
double precision, external :: memory_of_double, memory_of_int
|
||||
rss = memory_of_double(siz)*(N_int*2+1)
|
||||
call check_mem(rss,irp_here)
|
||||
|
||||
allocate(res%det(N_int, 2, siz), res%val(siz))
|
||||
|
||||
res%val(:) = 0d0
|
||||
res%det(:,:,:) = 0_8
|
||||
res%N = N
|
||||
res%mini = 0d0
|
||||
res%cur = 0
|
||||
end subroutine
|
||||
|
||||
subroutine delete_selection_buffer(b)
|
||||
use selection_types
|
||||
implicit none
|
||||
type(selection_buffer), intent(inout) :: b
|
||||
if (associated(b%det)) then
|
||||
deallocate(b%det)
|
||||
endif
|
||||
if (associated(b%val)) then
|
||||
deallocate(b%val)
|
||||
endif
|
||||
end
|
||||
|
||||
|
||||
subroutine add_to_selection_buffer(b, det, val)
|
||||
use selection_types
|
||||
implicit none
|
||||
|
||||
type(selection_buffer), intent(inout) :: b
|
||||
integer(bit_kind), intent(in) :: det(N_int, 2)
|
||||
double precision, intent(in) :: val
|
||||
integer :: i
|
||||
|
||||
if(b%N > 0 .and. val <= b%mini) then
|
||||
b%cur += 1
|
||||
b%det(1:N_int,1:2,b%cur) = det(1:N_int,1:2)
|
||||
b%val(b%cur) = val
|
||||
if(b%cur == size(b%val)) then
|
||||
call sort_selection_buffer(b)
|
||||
end if
|
||||
end if
|
||||
end subroutine
|
||||
|
||||
subroutine merge_selection_buffers(b1, b2)
|
||||
use selection_types
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Merges the selection buffers b1 and b2 into b2
|
||||
END_DOC
|
||||
type(selection_buffer), intent(inout) :: b1
|
||||
type(selection_buffer), intent(inout) :: b2
|
||||
integer(bit_kind), pointer :: detmp(:,:,:)
|
||||
double precision, pointer :: val(:)
|
||||
integer :: i, i1, i2, k, nmwen
|
||||
if (b1%cur == 0) return
|
||||
do while (b1%val(b1%cur) > b2%mini)
|
||||
b1%cur = b1%cur-1
|
||||
if (b1%cur == 0) then
|
||||
return
|
||||
endif
|
||||
enddo
|
||||
nmwen = min(b1%N, b1%cur+b2%cur)
|
||||
double precision :: rss
|
||||
double precision, external :: memory_of_double
|
||||
rss = memory_of_double(size(b1%val)) + 2*N_int*memory_of_double(size(b1%det,3))
|
||||
call check_mem(rss,irp_here)
|
||||
allocate( val(size(b1%val)), detmp(N_int, 2, size(b1%det,3)) )
|
||||
i1=1
|
||||
i2=1
|
||||
do i=1,nmwen
|
||||
if ( (i1 > b1%cur).and.(i2 > b2%cur) ) then
|
||||
exit
|
||||
else if (i1 > b1%cur) then
|
||||
val(i) = b2%val(i2)
|
||||
detmp(1:N_int,1,i) = b2%det(1:N_int,1,i2)
|
||||
detmp(1:N_int,2,i) = b2%det(1:N_int,2,i2)
|
||||
i2=i2+1
|
||||
else if (i2 > b2%cur) then
|
||||
val(i) = b1%val(i1)
|
||||
detmp(1:N_int,1,i) = b1%det(1:N_int,1,i1)
|
||||
detmp(1:N_int,2,i) = b1%det(1:N_int,2,i1)
|
||||
i1=i1+1
|
||||
else
|
||||
if (b1%val(i1) <= b2%val(i2)) then
|
||||
val(i) = b1%val(i1)
|
||||
detmp(1:N_int,1,i) = b1%det(1:N_int,1,i1)
|
||||
detmp(1:N_int,2,i) = b1%det(1:N_int,2,i1)
|
||||
i1=i1+1
|
||||
else
|
||||
val(i) = b2%val(i2)
|
||||
detmp(1:N_int,1,i) = b2%det(1:N_int,1,i2)
|
||||
detmp(1:N_int,2,i) = b2%det(1:N_int,2,i2)
|
||||
i2=i2+1
|
||||
endif
|
||||
endif
|
||||
enddo
|
||||
deallocate(b2%det, b2%val)
|
||||
do i=nmwen+1,b2%N
|
||||
val(i) = 0.d0
|
||||
detmp(1:N_int,1:2,i) = 0_bit_kind
|
||||
enddo
|
||||
b2%det => detmp
|
||||
b2%val => val
|
||||
b2%mini = min(b2%mini,b2%val(b2%N))
|
||||
b2%cur = nmwen
|
||||
end
|
||||
|
||||
|
||||
subroutine sort_selection_buffer(b)
|
||||
use selection_types
|
||||
implicit none
|
||||
|
||||
type(selection_buffer), intent(inout) :: b
|
||||
integer, allocatable :: iorder(:)
|
||||
integer(bit_kind), pointer :: detmp(:,:,:)
|
||||
integer :: i, nmwen
|
||||
logical, external :: detEq
|
||||
if (b%N == 0 .or. b%cur == 0) return
|
||||
nmwen = min(b%N, b%cur)
|
||||
|
||||
double precision :: rss
|
||||
double precision, external :: memory_of_double, memory_of_int
|
||||
rss = memory_of_int(b%cur) + 2*N_int*memory_of_double(size(b%det,3))
|
||||
call check_mem(rss,irp_here)
|
||||
allocate(iorder(b%cur), detmp(N_int, 2, size(b%det,3)))
|
||||
do i=1,b%cur
|
||||
iorder(i) = i
|
||||
end do
|
||||
call dsort(b%val, iorder, b%cur)
|
||||
do i=1, nmwen
|
||||
detmp(1:N_int,1,i) = b%det(1:N_int,1,iorder(i))
|
||||
detmp(1:N_int,2,i) = b%det(1:N_int,2,iorder(i))
|
||||
end do
|
||||
deallocate(b%det,iorder)
|
||||
b%det => detmp
|
||||
b%mini = min(b%mini,b%val(b%N))
|
||||
b%cur = nmwen
|
||||
end subroutine
|
||||
|
||||
subroutine make_selection_buffer_s2(b)
|
||||
use selection_types
|
||||
type(selection_buffer), intent(inout) :: b
|
||||
|
||||
integer(bit_kind), allocatable :: o(:,:,:)
|
||||
double precision, allocatable :: val(:)
|
||||
|
||||
integer :: n_d
|
||||
integer :: i,k,sze,n_alpha,j,n
|
||||
logical :: dup
|
||||
|
||||
! Sort
|
||||
integer, allocatable :: iorder(:)
|
||||
integer*8, allocatable :: bit_tmp(:)
|
||||
integer*8, external :: occ_pattern_search_key
|
||||
integer(bit_kind), allocatable :: tmp_array(:,:,:)
|
||||
logical, allocatable :: duplicate(:)
|
||||
|
||||
n_d = b%cur
|
||||
double precision :: rss
|
||||
double precision, external :: memory_of_double, memory_of_int
|
||||
rss = (4*N_int+4)*memory_of_double(n_d)
|
||||
call check_mem(rss,irp_here)
|
||||
allocate(o(N_int,2,n_d), iorder(n_d), duplicate(n_d), bit_tmp(n_d), &
|
||||
tmp_array(N_int,2,n_d), val(n_d) )
|
||||
|
||||
do i=1,n_d
|
||||
do k=1,N_int
|
||||
o(k,1,i) = ieor(b%det(k,1,i), b%det(k,2,i))
|
||||
o(k,2,i) = iand(b%det(k,1,i), b%det(k,2,i))
|
||||
enddo
|
||||
iorder(i) = i
|
||||
bit_tmp(i) = occ_pattern_search_key(o(1,1,i),N_int)
|
||||
enddo
|
||||
|
||||
deallocate(b%det)
|
||||
|
||||
call i8sort(bit_tmp,iorder,n_d)
|
||||
|
||||
do i=1,n_d
|
||||
do k=1,N_int
|
||||
tmp_array(k,1,i) = o(k,1,iorder(i))
|
||||
tmp_array(k,2,i) = o(k,2,iorder(i))
|
||||
enddo
|
||||
val(i) = b%val(iorder(i))
|
||||
duplicate(i) = .False.
|
||||
enddo
|
||||
|
||||
! Find duplicates
|
||||
do i=1,n_d-1
|
||||
if (duplicate(i)) then
|
||||
cycle
|
||||
endif
|
||||
j = i+1
|
||||
do while (bit_tmp(j)==bit_tmp(i))
|
||||
if (duplicate(j)) then
|
||||
j+=1
|
||||
if (j>n_d) then
|
||||
exit
|
||||
endif
|
||||
cycle
|
||||
endif
|
||||
dup = .True.
|
||||
do k=1,N_int
|
||||
if ( (tmp_array(k,1,i) /= tmp_array(k,1,j)) &
|
||||
.or. (tmp_array(k,2,i) /= tmp_array(k,2,j)) ) then
|
||||
dup = .False.
|
||||
exit
|
||||
endif
|
||||
enddo
|
||||
if (dup) then
|
||||
! val(i) = val(i) + val(j)
|
||||
val(i) = max(val(i), val(j))
|
||||
duplicate(j) = .True.
|
||||
endif
|
||||
j+=1
|
||||
if (j>n_d) then
|
||||
exit
|
||||
endif
|
||||
enddo
|
||||
enddo
|
||||
|
||||
deallocate (b%val)
|
||||
! Copy filtered result
|
||||
integer :: n_p
|
||||
n_p=0
|
||||
do i=1,n_d
|
||||
if (duplicate(i)) then
|
||||
cycle
|
||||
endif
|
||||
n_p = n_p + 1
|
||||
do k=1,N_int
|
||||
o(k,1,n_p) = tmp_array(k,1,i)
|
||||
o(k,2,n_p) = tmp_array(k,2,i)
|
||||
enddo
|
||||
val(n_p) = val(i)
|
||||
enddo
|
||||
|
||||
! Sort by importance
|
||||
do i=1,n_p
|
||||
iorder(i) = i
|
||||
end do
|
||||
call dsort(val,iorder,n_p)
|
||||
do i=1,n_p
|
||||
do k=1,N_int
|
||||
tmp_array(k,1,i) = o(k,1,iorder(i))
|
||||
tmp_array(k,2,i) = o(k,2,iorder(i))
|
||||
enddo
|
||||
enddo
|
||||
do i=1,n_p
|
||||
do k=1,N_int
|
||||
o(k,1,i) = tmp_array(k,1,i)
|
||||
o(k,2,i) = tmp_array(k,2,i)
|
||||
enddo
|
||||
enddo
|
||||
|
||||
! Create determinants
|
||||
n_d = 0
|
||||
do i=1,n_p
|
||||
call occ_pattern_to_dets_size(o(1,1,i),sze,elec_alpha_num,N_int)
|
||||
n_d = n_d + sze
|
||||
if (n_d > b%cur) then
|
||||
! if (n_d - b%cur > b%cur - n_d + sze) then
|
||||
! n_d = n_d - sze
|
||||
! endif
|
||||
exit
|
||||
endif
|
||||
enddo
|
||||
|
||||
rss = (4*N_int+2)*memory_of_double(n_d)
|
||||
call check_mem(rss,irp_here)
|
||||
allocate(b%det(N_int,2,2*n_d), b%val(2*n_d))
|
||||
k=1
|
||||
do i=1,n_p
|
||||
n=n_d
|
||||
call occ_pattern_to_dets_size(o(1,1,i),n,elec_alpha_num,N_int)
|
||||
call occ_pattern_to_dets(o(1,1,i),b%det(1,1,k),n,elec_alpha_num,N_int)
|
||||
do j=k,k+n-1
|
||||
b%val(j) = val(i)
|
||||
enddo
|
||||
k = k+n
|
||||
if (k > n_d) exit
|
||||
enddo
|
||||
deallocate(o)
|
||||
b%N = 2*n_d
|
||||
b%cur = n_d
|
||||
end
|
||||
|
||||
|
||||
|
||||
|
|
@ -1,9 +0,0 @@
|
|||
module selection_types
|
||||
type selection_buffer
|
||||
integer :: N, cur
|
||||
integer(8) , pointer :: det(:,:,:)
|
||||
double precision, pointer :: val(:)
|
||||
double precision :: mini
|
||||
endtype
|
||||
end module
|
||||
|
||||
|
|
@ -1,219 +0,0 @@
|
|||
subroutine ZMQ_selection(N_in, pt2, variance, norm)
|
||||
use f77_zmq
|
||||
use selection_types
|
||||
|
||||
implicit none
|
||||
|
||||
integer(ZMQ_PTR) :: zmq_to_qp_run_socket , zmq_socket_pull
|
||||
integer, intent(in) :: N_in
|
||||
type(selection_buffer) :: b
|
||||
integer :: i, N
|
||||
integer, external :: omp_get_thread_num
|
||||
double precision, intent(out) :: pt2(N_states)
|
||||
double precision, intent(out) :: variance(N_states)
|
||||
double precision, intent(out) :: norm(N_states)
|
||||
|
||||
! PROVIDE psi_det psi_coef N_det qp_max_mem N_states pt2_F s2_eig N_det_generators
|
||||
|
||||
N = max(N_in,1)
|
||||
if (.True.) then
|
||||
PROVIDE pt2_e0_denominator nproc
|
||||
PROVIDE psi_bilinear_matrix_columns_loc psi_det_alpha_unique psi_det_beta_unique
|
||||
PROVIDE psi_bilinear_matrix_rows psi_det_sorted_order psi_bilinear_matrix_order
|
||||
PROVIDE psi_bilinear_matrix_transp_rows_loc psi_bilinear_matrix_transp_columns
|
||||
PROVIDE psi_bilinear_matrix_transp_order
|
||||
|
||||
call new_parallel_job(zmq_to_qp_run_socket,zmq_socket_pull,'selection')
|
||||
|
||||
integer, external :: zmq_put_psi
|
||||
integer, external :: zmq_put_N_det_generators
|
||||
integer, external :: zmq_put_N_det_selectors
|
||||
integer, external :: zmq_put_dvector
|
||||
|
||||
if (zmq_put_psi(zmq_to_qp_run_socket,1) == -1) then
|
||||
stop 'Unable to put psi on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_N_det_generators(zmq_to_qp_run_socket, 1) == -1) then
|
||||
stop 'Unable to put N_det_generators on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_N_det_selectors(zmq_to_qp_run_socket, 1) == -1) then
|
||||
stop 'Unable to put N_det_selectors on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_dvector(zmq_to_qp_run_socket,1,'energy',pt2_e0_denominator,size(pt2_e0_denominator)) == -1) then
|
||||
stop 'Unable to put energy on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_dvector(zmq_to_qp_run_socket,1,'threshold_selectors',threshold_selectors,1) == -1) then
|
||||
stop 'Unable to put threshold_selectors on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_dvector(zmq_to_qp_run_socket,1,'state_average_weight',state_average_weight,N_states) == -1) then
|
||||
stop 'Unable to put state_average_weight on ZMQ server'
|
||||
endif
|
||||
if (zmq_put_dvector(zmq_to_qp_run_socket,1,'threshold_generators',threshold_generators,1) == -1) then
|
||||
stop 'Unable to put threshold_generators on ZMQ server'
|
||||
endif
|
||||
call create_selection_buffer(N, N*2, b)
|
||||
endif
|
||||
|
||||
integer, external :: add_task_to_taskserver
|
||||
character(len=100000) :: task
|
||||
integer :: j,k,ipos
|
||||
ipos=1
|
||||
task = ' '
|
||||
|
||||
do i= 1, N_det_generators
|
||||
do j=1,pt2_F(i)
|
||||
write(task(ipos:ipos+30),'(I9,1X,I9,1X,I9,''|'')') j, i, N
|
||||
ipos += 30
|
||||
if (ipos > 100000-30) then
|
||||
if (add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos))) == -1) then
|
||||
stop 'Unable to add task to task server'
|
||||
endif
|
||||
ipos=1
|
||||
endif
|
||||
end do
|
||||
enddo
|
||||
if (ipos > 1) then
|
||||
if (add_task_to_taskserver(zmq_to_qp_run_socket,trim(task(1:ipos))) == -1) then
|
||||
stop 'Unable to add task to task server'
|
||||
endif
|
||||
endif
|
||||
|
||||
|
||||
ASSERT (associated(b%det))
|
||||
ASSERT (associated(b%val))
|
||||
|
||||
integer, external :: zmq_set_running
|
||||
if (zmq_set_running(zmq_to_qp_run_socket) == -1) then
|
||||
print *, irp_here, ': Failed in zmq_set_running'
|
||||
endif
|
||||
|
||||
integer :: nproc_target
|
||||
nproc_target = nproc
|
||||
double precision :: mem
|
||||
mem = 8.d0 * N_det * (N_int * 2.d0 * 3.d0 + 3.d0 + 5.d0) / (1024.d0**3)
|
||||
call write_double(6,mem,'Estimated memory/thread (Gb)')
|
||||
if (qp_max_mem > 0) then
|
||||
nproc_target = max(1,int(dble(qp_max_mem)/mem))
|
||||
nproc_target = min(nproc_target,nproc)
|
||||
endif
|
||||
|
||||
f(:) = 1.d0
|
||||
if (.not.do_pt2) then
|
||||
double precision :: f(N_states), u_dot_u
|
||||
do k=1,min(N_det,N_states)
|
||||
f(k) = 1.d0 / u_dot_u(psi_selectors_coef(1,k), N_det_selectors)
|
||||
enddo
|
||||
endif
|
||||
|
||||
!$OMP PARALLEL DEFAULT(shared) SHARED(b, pt2, variance, norm) PRIVATE(i) NUM_THREADS(nproc_target+1)
|
||||
i = omp_get_thread_num()
|
||||
if (i==0) then
|
||||
call selection_collector(zmq_socket_pull, b, N, pt2, variance, norm)
|
||||
else
|
||||
call selection_slave_inproc(i)
|
||||
endif
|
||||
!$OMP END PARALLEL
|
||||
call end_parallel_job(zmq_to_qp_run_socket, zmq_socket_pull, 'selection')
|
||||
do i=N_det+1,N_states
|
||||
pt2(i) = 0.d0
|
||||
variance(i) = 0.d0
|
||||
norm(i) = 0.d0
|
||||
enddo
|
||||
if (N_in > 0) then
|
||||
if (s2_eig) then
|
||||
call make_selection_buffer_s2(b)
|
||||
endif
|
||||
call fill_H_apply_buffer_no_selection(b%cur,b%det,N_int,0)
|
||||
call copy_H_apply_buffer_to_wf()
|
||||
call save_wavefunction
|
||||
endif
|
||||
call delete_selection_buffer(b)
|
||||
do k=1,N_states
|
||||
pt2(k) = pt2(k) * f(k)
|
||||
variance(k) = variance(k) * f(k)
|
||||
norm(k) = norm(k) * f(k)
|
||||
enddo
|
||||
|
||||
end subroutine
|
||||
|
||||
|
||||
subroutine selection_slave_inproc(i)
|
||||
implicit none
|
||||
integer, intent(in) :: i
|
||||
|
||||
call run_selection_slave(1,i,pt2_e0_denominator)
|
||||
end
|
||||
|
||||
subroutine selection_collector(zmq_socket_pull, b, N, pt2, variance, norm)
|
||||
use f77_zmq
|
||||
use selection_types
|
||||
use bitmasks
|
||||
implicit none
|
||||
|
||||
|
||||
integer(ZMQ_PTR), intent(in) :: zmq_socket_pull
|
||||
type(selection_buffer), intent(inout) :: b
|
||||
integer, intent(in) :: N
|
||||
double precision, intent(inout) :: pt2(N_states)
|
||||
double precision, intent(inout) :: variance(N_states)
|
||||
double precision, intent(inout) :: norm(N_states)
|
||||
double precision :: pt2_mwen(N_states)
|
||||
double precision :: variance_mwen(N_states)
|
||||
double precision :: norm_mwen(N_states)
|
||||
integer(ZMQ_PTR),external :: new_zmq_to_qp_run_socket
|
||||
integer(ZMQ_PTR) :: zmq_to_qp_run_socket
|
||||
|
||||
integer(ZMQ_PTR), external :: new_zmq_pull_socket
|
||||
|
||||
integer :: msg_size, rc, more
|
||||
integer :: acc, i, j, robin, ntask
|
||||
double precision, pointer :: val(:)
|
||||
integer(bit_kind), pointer :: det(:,:,:)
|
||||
integer, allocatable :: task_id(:)
|
||||
type(selection_buffer) :: b2
|
||||
|
||||
|
||||
|
||||
|
||||
zmq_to_qp_run_socket = new_zmq_to_qp_run_socket()
|
||||
call create_selection_buffer(N, N*2, b2)
|
||||
double precision :: rss
|
||||
double precision, external :: memory_of_int
|
||||
rss = memory_of_int(N_det_generators)
|
||||
call check_mem(rss,irp_here)
|
||||
allocate(task_id(N_det_generators))
|
||||
more = 1
|
||||
pt2(:) = 0d0
|
||||
variance(:) = 0.d0
|
||||
norm(:) = 0.d0
|
||||
pt2_mwen(:) = 0.d0
|
||||
variance_mwen(:) = 0.d0
|
||||
norm_mwen(:) = 0.d0
|
||||
do while (more == 1)
|
||||
call pull_selection_results(zmq_socket_pull, pt2_mwen, variance_mwen, norm_mwen, b2%val(1), b2%det(1,1,1), b2%cur, task_id, ntask)
|
||||
|
||||
pt2(:) += pt2_mwen(:)
|
||||
variance(:) += variance_mwen(:)
|
||||
norm(:) += norm_mwen(:)
|
||||
do i=1, b2%cur
|
||||
call add_to_selection_buffer(b, b2%det(1,1,i), b2%val(i))
|
||||
if (b2%val(i) > b%mini) exit
|
||||
end do
|
||||
|
||||
do i=1, ntask
|
||||
if(task_id(i) == 0) then
|
||||
print *, "Error in collector"
|
||||
endif
|
||||
integer, external :: zmq_delete_task
|
||||
if (zmq_delete_task(zmq_to_qp_run_socket,zmq_socket_pull,task_id(i),more) == -1) then
|
||||
stop 'Unable to delete task'
|
||||
endif
|
||||
end do
|
||||
end do
|
||||
|
||||
|
||||
call delete_selection_buffer(b2)
|
||||
call sort_selection_buffer(b)
|
||||
call end_zmq_to_qp_run_socket(zmq_to_qp_run_socket)
|
||||
end subroutine
|
||||
|
||||
|
|
@ -1,12 +0,0 @@
|
|||
program four_idx
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! 4-index transformation from AO to MO integrals
|
||||
END_DOC
|
||||
|
||||
disk_access_mo_integrals = 'Write'
|
||||
SOFT_TOUCH disk_access_mo_integrals
|
||||
if (.true.) then
|
||||
PROVIDE mo_bielec_integrals_in_map
|
||||
endif
|
||||
end
|
||||
|
|
@ -1 +0,0 @@
|
|||
Integrals_Bielec
|
||||
|
|
@ -1,6 +0,0 @@
|
|||
=======
|
||||
FourIdx
|
||||
=======
|
||||
|
||||
Four-index transformation.
|
||||
|
||||
|
|
@ -1 +0,0 @@
|
|||
Determinants
|
||||
|
|
@ -1,12 +0,0 @@
|
|||
==============
|
||||
Generators_CAS
|
||||
==============
|
||||
|
||||
Module defining the generator determinants as those belonging to a |CAS|.
|
||||
The |MOs| belonging to the |CAS| are those which were set as active with
|
||||
the :ref:`qp_set_mo_class` command.
|
||||
|
||||
This module is intended to be included in the :file:`NEED` file to define
|
||||
generators on a |CAS|.
|
||||
|
||||
|
||||
|
|
@ -1,83 +0,0 @@
|
|||
use bitmasks
|
||||
|
||||
BEGIN_PROVIDER [ integer, N_det_generators ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Number of generator detetrminants
|
||||
END_DOC
|
||||
integer :: i,k,l
|
||||
logical :: good
|
||||
integer, external :: number_of_holes,number_of_particles
|
||||
call write_time(6)
|
||||
N_det_generators = 0
|
||||
do i=1,N_det
|
||||
good = ( number_of_holes(psi_det_sorted(1,1,i)) ==0).and.(number_of_particles(psi_det_sorted(1,1,i))==0 )
|
||||
if (good) then
|
||||
N_det_generators += 1
|
||||
endif
|
||||
enddo
|
||||
N_det_generators = max(N_det_generators,1)
|
||||
call write_int(6,N_det_generators,'Number of generators')
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ integer(bit_kind), psi_det_generators, (N_int,2,psi_det_size) ]
|
||||
&BEGIN_PROVIDER [ double precision, psi_coef_generators, (psi_det_size,N_states) ]
|
||||
&BEGIN_PROVIDER [ integer(bit_kind), psi_det_sorted_gen, (N_int,2,psi_det_size) ]
|
||||
&BEGIN_PROVIDER [ double precision, psi_coef_sorted_gen, (psi_det_size,N_states) ]
|
||||
&BEGIN_PROVIDER [ integer, psi_det_sorted_gen_order, (psi_det_size) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! For Single reference wave functions, the generator is the
|
||||
! Hartree-Fock determinant
|
||||
END_DOC
|
||||
integer :: i, k, l, m
|
||||
logical :: good
|
||||
integer, external :: number_of_holes,number_of_particles
|
||||
integer, allocatable :: nongen(:)
|
||||
integer :: inongen
|
||||
|
||||
allocate(nongen(N_det))
|
||||
|
||||
inongen = 0
|
||||
m=0
|
||||
do i=1,N_det
|
||||
good = ( number_of_holes(psi_det_sorted(1,1,i)) ==0).and.(number_of_particles(psi_det_sorted(1,1,i))==0 )
|
||||
if (good) then
|
||||
m = m+1
|
||||
psi_det_sorted_gen_order(i) = m
|
||||
do k=1,N_int
|
||||
psi_det_generators(k,1,m) = psi_det_sorted(k,1,i)
|
||||
psi_det_generators(k,2,m) = psi_det_sorted(k,2,i)
|
||||
enddo
|
||||
psi_coef_generators(m,:) = psi_coef_sorted(i,:)
|
||||
else
|
||||
inongen += 1
|
||||
nongen(inongen) = i
|
||||
endif
|
||||
enddo
|
||||
|
||||
psi_det_sorted_gen(:,:,:N_det_generators) = psi_det_generators(:,:,:N_det_generators)
|
||||
psi_coef_sorted_gen(:N_det_generators, :) = psi_coef_generators(:N_det_generators, :)
|
||||
do i=1,inongen
|
||||
psi_det_sorted_gen_order(nongen(i)) = N_det_generators+i
|
||||
psi_det_sorted_gen(:,:,N_det_generators+i) = psi_det_sorted(:,:,nongen(i))
|
||||
psi_coef_sorted_gen(N_det_generators+i, :) = psi_coef_sorted(nongen(i),:)
|
||||
end do
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ integer, size_select_max]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Size of the select_max array
|
||||
END_DOC
|
||||
size_select_max = 10000
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ double precision, select_max, (size_select_max) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Memo to skip useless selectors
|
||||
END_DOC
|
||||
select_max = huge(1.d0)
|
||||
END_PROVIDER
|
||||
|
||||
|
|
@ -1,2 +0,0 @@
|
|||
Determinants
|
||||
Hartree_Fock
|
||||
|
|
@ -1,9 +0,0 @@
|
|||
===============
|
||||
Generators_full
|
||||
===============
|
||||
|
||||
Module defining the generator determinants as all the determinants of the
|
||||
variational space.
|
||||
|
||||
This module is intended to be included in the :file:`NEED` file to define
|
||||
a full set of generators.
|
||||
|
|
@ -1,82 +0,0 @@
|
|||
use bitmasks
|
||||
|
||||
BEGIN_PROVIDER [ integer, N_det_generators ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! For Single reference wave functions, the number of generators is 1 : the
|
||||
! Hartree-Fock determinant
|
||||
END_DOC
|
||||
integer :: i
|
||||
double precision :: norm
|
||||
call write_time(6)
|
||||
norm = 1.d0
|
||||
N_det_generators = N_det
|
||||
do i=1,N_det
|
||||
norm = norm - psi_average_norm_contrib_sorted(i)
|
||||
if (norm - 1.d-10 < 1.d0 - threshold_generators) then
|
||||
N_det_generators = i
|
||||
exit
|
||||
endif
|
||||
enddo
|
||||
N_det_generators = max(N_det_generators,1)
|
||||
call write_int(6,N_det_generators,'Number of generators')
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ integer(bit_kind), psi_det_generators, (N_int,2,psi_det_size) ]
|
||||
&BEGIN_PROVIDER [ double precision, psi_coef_generators, (psi_det_size,N_states) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! For Single reference wave functions, the generator is the
|
||||
! Hartree-Fock determinant
|
||||
END_DOC
|
||||
psi_det_generators(1:N_int,1:2,1:N_det) = psi_det_sorted(1:N_int,1:2,1:N_det)
|
||||
psi_coef_generators(1:N_det,1:N_states) = psi_coef_sorted(1:N_det,1:N_states)
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ integer(bit_kind), psi_det_sorted_gen, (N_int,2,psi_det_size) ]
|
||||
&BEGIN_PROVIDER [ double precision, psi_coef_sorted_gen, (psi_det_size,N_states) ]
|
||||
&BEGIN_PROVIDER [ integer, psi_det_sorted_gen_order, (psi_det_size) ]
|
||||
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! For Single reference wave functions, the generator is the
|
||||
! Hartree-Fock determinant
|
||||
END_DOC
|
||||
psi_det_sorted_gen = psi_det_sorted
|
||||
psi_coef_sorted_gen = psi_coef_sorted
|
||||
psi_det_sorted_gen_order = psi_det_sorted_order
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
BEGIN_PROVIDER [integer, degree_max_generators]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Max degree of excitation (respect to HF) of the generators
|
||||
END_DOC
|
||||
integer :: i,degree
|
||||
degree_max_generators = 0
|
||||
do i = 1, N_det_generators
|
||||
call get_excitation_degree(HF_bitmask,psi_det_generators(1,1,i),degree,N_int)
|
||||
if(degree .gt. degree_max_generators)then
|
||||
degree_max_generators = degree
|
||||
endif
|
||||
enddo
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ integer, size_select_max]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Size of the select_max array
|
||||
END_DOC
|
||||
size_select_max = 10000
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ double precision, select_max, (size_select_max) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Memo to skip useless selectors
|
||||
END_DOC
|
||||
select_max = huge(1.d0)
|
||||
END_PROVIDER
|
||||
|
||||
|
|
@ -1,139 +0,0 @@
|
|||
BEGIN_PROVIDER [ double precision, threshold_DIIS_nonzero ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! If threshold_DIIS is zero, choose sqrt(thresh_scf)
|
||||
END_DOC
|
||||
if (threshold_DIIS == 0.d0) then
|
||||
threshold_DIIS_nonzero = dsqrt(thresh_scf)
|
||||
else
|
||||
threshold_DIIS_nonzero = threshold_DIIS
|
||||
endif
|
||||
ASSERT (threshold_DIIS_nonzero >= 0.d0)
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [double precision, FPS_SPF_Matrix_AO, (AO_num, AO_num)]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Commutator FPS - SPF
|
||||
END_DOC
|
||||
double precision, allocatable :: scratch(:,:)
|
||||
allocate( &
|
||||
scratch(AO_num, AO_num) &
|
||||
)
|
||||
|
||||
! Compute FP
|
||||
|
||||
call dgemm('N','N',AO_num,AO_num,AO_num, &
|
||||
1.d0, &
|
||||
Fock_Matrix_AO,Size(Fock_Matrix_AO,1), &
|
||||
HF_Density_Matrix_AO,Size(HF_Density_Matrix_AO,1), &
|
||||
0.d0, &
|
||||
scratch,Size(scratch,1))
|
||||
|
||||
! Compute FPS
|
||||
|
||||
call dgemm('N','N',AO_num,AO_num,AO_num, &
|
||||
1.d0, &
|
||||
scratch,Size(scratch,1), &
|
||||
AO_Overlap,Size(AO_Overlap,1), &
|
||||
0.d0, &
|
||||
FPS_SPF_Matrix_AO,Size(FPS_SPF_Matrix_AO,1))
|
||||
|
||||
! Compute SP
|
||||
|
||||
call dgemm('N','N',AO_num,AO_num,AO_num, &
|
||||
1.d0, &
|
||||
AO_Overlap,Size(AO_Overlap,1), &
|
||||
HF_Density_Matrix_AO,Size(HF_Density_Matrix_AO,1), &
|
||||
0.d0, &
|
||||
scratch,Size(scratch,1))
|
||||
|
||||
! Compute FPS - SPF
|
||||
|
||||
call dgemm('N','N',AO_num,AO_num,AO_num, &
|
||||
-1.d0, &
|
||||
scratch,Size(scratch,1), &
|
||||
Fock_Matrix_AO,Size(Fock_Matrix_AO,1), &
|
||||
1.d0, &
|
||||
FPS_SPF_Matrix_AO,Size(FPS_SPF_Matrix_AO,1))
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
bEGIN_PROVIDER [double precision, FPS_SPF_Matrix_MO, (mo_tot_num, mo_tot_num)]
|
||||
implicit none
|
||||
begin_doc
|
||||
! Commutator FPS - SPF in MO basis
|
||||
end_doc
|
||||
call ao_to_mo(FPS_SPF_Matrix_AO, size(FPS_SPF_Matrix_AO,1), &
|
||||
FPS_SPF_Matrix_MO, size(FPS_SPF_Matrix_MO,1))
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
BEGIN_PROVIDER [ double precision, eigenvalues_Fock_matrix_AO, (AO_num) ]
|
||||
&BEGIN_PROVIDER [ double precision, eigenvectors_Fock_matrix_AO, (AO_num,AO_num) ]
|
||||
|
||||
BEGIN_DOC
|
||||
! Eigenvalues and eigenvectors of the Fock matrix over the AO basis
|
||||
END_DOC
|
||||
|
||||
implicit none
|
||||
|
||||
double precision, allocatable :: scratch(:,:),work(:),Xt(:,:)
|
||||
integer :: lwork,info
|
||||
integer :: i,j
|
||||
|
||||
lwork = 3*AO_num - 1
|
||||
allocate( &
|
||||
scratch(AO_num,AO_num), &
|
||||
work(lwork), &
|
||||
Xt(AO_num,AO_num) &
|
||||
)
|
||||
|
||||
! Calculate Xt
|
||||
|
||||
do i=1,AO_num
|
||||
do j=1,AO_num
|
||||
Xt(i,j) = S_half_inv(j,i)
|
||||
enddo
|
||||
enddo
|
||||
|
||||
! Calculate Fock matrix in orthogonal basis: F' = Xt.F.X
|
||||
|
||||
call dgemm('N','N',AO_num,AO_num,AO_num, &
|
||||
1.d0, &
|
||||
Fock_matrix_AO,size(Fock_matrix_AO,1), &
|
||||
S_half_inv,size(S_half_inv,1), &
|
||||
0.d0, &
|
||||
eigenvectors_Fock_matrix_AO,size(eigenvectors_Fock_matrix_AO,1))
|
||||
|
||||
call dgemm('N','N',AO_num,AO_num,AO_num, &
|
||||
1.d0, &
|
||||
Xt,size(Xt,1), &
|
||||
eigenvectors_Fock_matrix_AO,size(eigenvectors_Fock_matrix_AO,1), &
|
||||
0.d0, &
|
||||
scratch,size(scratch,1))
|
||||
|
||||
! Diagonalize F' to obtain eigenvectors in orthogonal basis C' and eigenvalues
|
||||
|
||||
call dsyev('V','U',AO_num, &
|
||||
scratch,size(scratch,1), &
|
||||
eigenvalues_Fock_matrix_AO, &
|
||||
work,lwork,info)
|
||||
|
||||
if(info /= 0) then
|
||||
print *, irp_here//' failed : ', info
|
||||
stop 1
|
||||
endif
|
||||
|
||||
! Back-transform eigenvectors: C =X.C'
|
||||
|
||||
call dgemm('N','N',AO_num,AO_num,AO_num, &
|
||||
1.d0, &
|
||||
S_half_inv,size(S_half_inv,1), &
|
||||
scratch,size(scratch,1), &
|
||||
0.d0, &
|
||||
eigenvectors_Fock_matrix_AO,size(eigenvectors_Fock_matrix_AO,1))
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
|
|
@ -1,53 +0,0 @@
|
|||
[max_dim_diis]
|
||||
type: integer
|
||||
doc: Maximum size of the |DIIS| extrapolation procedure
|
||||
interface: ezfio,provider,ocaml
|
||||
default: 15
|
||||
|
||||
[threshold_diis]
|
||||
type: Threshold
|
||||
doc: Threshold on the convergence of the |DIIS| error vector during a Hartree-Fock calculation. If 0. is chosen, the square root of thresh_scf will be used.
|
||||
interface: ezfio,provider,ocaml
|
||||
default: 0.
|
||||
|
||||
[thresh_scf]
|
||||
type: Threshold
|
||||
doc: Threshold on the convergence of the Hartree Fock energy.
|
||||
interface: ezfio,provider,ocaml
|
||||
default: 1.e-10
|
||||
|
||||
[n_it_scf_max]
|
||||
type: Strictly_positive_int
|
||||
doc: Maximum number of |SCF| iterations
|
||||
interface: ezfio,provider,ocaml
|
||||
default: 500
|
||||
|
||||
[level_shift]
|
||||
type: Positive_float
|
||||
doc: Initial value of the energy shift on the virtual |MOs|
|
||||
interface: ezfio,provider,ocaml
|
||||
default: 0.0
|
||||
|
||||
[scf_algorithm]
|
||||
type: character*(32)
|
||||
doc: Type of |SCF| algorithm used. Possible choices are [ Simple | DIIS]
|
||||
interface: ezfio,provider,ocaml
|
||||
default: DIIS
|
||||
|
||||
[mo_guess_type]
|
||||
type: MO_guess
|
||||
doc: Initial MO guess. Can be [ Huckel | HCore ]
|
||||
interface: ezfio,provider,ocaml
|
||||
default: Huckel
|
||||
|
||||
[energy]
|
||||
type: double precision
|
||||
doc: Calculated HF energy
|
||||
interface: ezfio
|
||||
|
||||
[no_oa_or_av_opt]
|
||||
type: logical
|
||||
doc: If |true|, skip the (inactive+core) --> (active) and the (active) --> (virtual) orbital rotations within the |SCF| procedure
|
||||
interface: ezfio,provider,ocaml
|
||||
default: False
|
||||
|
||||
|
|
@ -1,320 +0,0 @@
|
|||
BEGIN_PROVIDER [ double precision, Fock_matrix_mo, (mo_tot_num,mo_tot_num) ]
|
||||
&BEGIN_PROVIDER [ double precision, Fock_matrix_diag_mo, (mo_tot_num)]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Fock matrix on the MO basis.
|
||||
! For open shells, the ROHF Fock Matrix is
|
||||
!
|
||||
! | F-K | F + K/2 | F |
|
||||
! |---------------------------------|
|
||||
! | F + K/2 | F | F - K/2 |
|
||||
! |---------------------------------|
|
||||
! | F | F - K/2 | F + K |
|
||||
!
|
||||
! F = 1/2 (Fa + Fb)
|
||||
!
|
||||
! K = Fb - Fa
|
||||
!
|
||||
END_DOC
|
||||
integer :: i,j,n
|
||||
if (elec_alpha_num == elec_beta_num) then
|
||||
Fock_matrix_mo = Fock_matrix_mo_alpha
|
||||
else
|
||||
|
||||
do j=1,elec_beta_num
|
||||
! F-K
|
||||
do i=1,elec_beta_num
|
||||
Fock_matrix_mo(i,j) = 0.5d0*(Fock_matrix_mo_alpha(i,j)+Fock_matrix_mo_beta(i,j))&
|
||||
- (Fock_matrix_mo_beta(i,j) - Fock_matrix_mo_alpha(i,j))
|
||||
enddo
|
||||
! F+K/2
|
||||
do i=elec_beta_num+1,elec_alpha_num
|
||||
Fock_matrix_mo(i,j) = 0.5d0*(Fock_matrix_mo_alpha(i,j)+Fock_matrix_mo_beta(i,j))&
|
||||
+ 0.5d0*(Fock_matrix_mo_beta(i,j) - Fock_matrix_mo_alpha(i,j))
|
||||
enddo
|
||||
! F
|
||||
do i=elec_alpha_num+1, mo_tot_num
|
||||
Fock_matrix_mo(i,j) = 0.5d0*(Fock_matrix_mo_alpha(i,j)+Fock_matrix_mo_beta(i,j))
|
||||
enddo
|
||||
enddo
|
||||
|
||||
do j=elec_beta_num+1,elec_alpha_num
|
||||
! F+K/2
|
||||
do i=1,elec_beta_num
|
||||
Fock_matrix_mo(i,j) = 0.5d0*(Fock_matrix_mo_alpha(i,j)+Fock_matrix_mo_beta(i,j))&
|
||||
+ 0.5d0*(Fock_matrix_mo_beta(i,j) - Fock_matrix_mo_alpha(i,j))
|
||||
enddo
|
||||
! F
|
||||
do i=elec_beta_num+1,elec_alpha_num
|
||||
Fock_matrix_mo(i,j) = 0.5d0*(Fock_matrix_mo_alpha(i,j)+Fock_matrix_mo_beta(i,j))
|
||||
enddo
|
||||
! F-K/2
|
||||
do i=elec_alpha_num+1, mo_tot_num
|
||||
Fock_matrix_mo(i,j) = 0.5d0*(Fock_matrix_mo_alpha(i,j)+Fock_matrix_mo_beta(i,j))&
|
||||
- 0.5d0*(Fock_matrix_mo_beta(i,j) - Fock_matrix_mo_alpha(i,j))
|
||||
enddo
|
||||
enddo
|
||||
|
||||
do j=elec_alpha_num+1, mo_tot_num
|
||||
! F
|
||||
do i=1,elec_beta_num
|
||||
Fock_matrix_mo(i,j) = 0.5d0*(Fock_matrix_mo_alpha(i,j)+Fock_matrix_mo_beta(i,j))
|
||||
enddo
|
||||
! F-K/2
|
||||
do i=elec_beta_num+1,elec_alpha_num
|
||||
Fock_matrix_mo(i,j) = 0.5d0*(Fock_matrix_mo_alpha(i,j)+Fock_matrix_mo_beta(i,j))&
|
||||
- 0.5d0*(Fock_matrix_mo_beta(i,j) - Fock_matrix_mo_alpha(i,j))
|
||||
enddo
|
||||
! F+K
|
||||
do i=elec_alpha_num+1,mo_tot_num
|
||||
Fock_matrix_mo(i,j) = 0.5d0*(Fock_matrix_mo_alpha(i,j)+Fock_matrix_mo_beta(i,j)) &
|
||||
+ (Fock_matrix_mo_beta(i,j) - Fock_matrix_mo_alpha(i,j))
|
||||
enddo
|
||||
enddo
|
||||
|
||||
endif
|
||||
|
||||
do i = 1, mo_tot_num
|
||||
Fock_matrix_diag_mo(i) = Fock_matrix_mo(i,i)
|
||||
enddo
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
|
||||
BEGIN_PROVIDER [ double precision, Fock_matrix_ao_alpha, (ao_num, ao_num) ]
|
||||
&BEGIN_PROVIDER [ double precision, Fock_matrix_ao_beta, (ao_num, ao_num) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Alpha Fock matrix in AO basis set
|
||||
END_DOC
|
||||
|
||||
integer :: i,j
|
||||
do j=1,ao_num
|
||||
do i=1,ao_num
|
||||
Fock_matrix_ao_alpha(i,j) = ao_mono_elec_integral(i,j) + ao_bi_elec_integral_alpha(i,j)
|
||||
Fock_matrix_ao_beta (i,j) = ao_mono_elec_integral(i,j) + ao_bi_elec_integral_beta (i,j)
|
||||
enddo
|
||||
enddo
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
BEGIN_PROVIDER [ double precision, ao_bi_elec_integral_alpha, (ao_num, ao_num) ]
|
||||
&BEGIN_PROVIDER [ double precision, ao_bi_elec_integral_beta , (ao_num, ao_num) ]
|
||||
use map_module
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Alpha Fock matrix in AO basis set
|
||||
END_DOC
|
||||
|
||||
integer :: i,j,k,l,k1,r,s
|
||||
integer :: i0,j0,k0,l0
|
||||
integer*8 :: p,q
|
||||
double precision :: integral, c0, c1, c2
|
||||
double precision :: ao_bielec_integral, local_threshold
|
||||
double precision, allocatable :: ao_bi_elec_integral_alpha_tmp(:,:)
|
||||
double precision, allocatable :: ao_bi_elec_integral_beta_tmp(:,:)
|
||||
|
||||
ao_bi_elec_integral_alpha = 0.d0
|
||||
ao_bi_elec_integral_beta = 0.d0
|
||||
if (do_direct_integrals) then
|
||||
|
||||
!$OMP PARALLEL DEFAULT(NONE) &
|
||||
!$OMP PRIVATE(i,j,l,k1,k,integral,ii,jj,kk,ll,i8,keys,values,p,q,r,s,i0,j0,k0,l0, &
|
||||
!$OMP ao_bi_elec_integral_alpha_tmp,ao_bi_elec_integral_beta_tmp, c0, c1, c2, &
|
||||
!$OMP local_threshold)&
|
||||
!$OMP SHARED(ao_num,HF_density_matrix_ao_alpha,HF_density_matrix_ao_beta,&
|
||||
!$OMP ao_integrals_map,ao_integrals_threshold, ao_bielec_integral_schwartz, &
|
||||
!$OMP ao_overlap_abs, ao_bi_elec_integral_alpha, ao_bi_elec_integral_beta)
|
||||
|
||||
allocate(keys(1), values(1))
|
||||
allocate(ao_bi_elec_integral_alpha_tmp(ao_num,ao_num), &
|
||||
ao_bi_elec_integral_beta_tmp(ao_num,ao_num))
|
||||
ao_bi_elec_integral_alpha_tmp = 0.d0
|
||||
ao_bi_elec_integral_beta_tmp = 0.d0
|
||||
|
||||
q = ao_num*ao_num*ao_num*ao_num
|
||||
!$OMP DO SCHEDULE(dynamic)
|
||||
do p=1_8,q
|
||||
call bielec_integrals_index_reverse(kk,ii,ll,jj,p)
|
||||
if ( (kk(1)>ao_num).or. &
|
||||
(ii(1)>ao_num).or. &
|
||||
(jj(1)>ao_num).or. &
|
||||
(ll(1)>ao_num) ) then
|
||||
cycle
|
||||
endif
|
||||
k = kk(1)
|
||||
i = ii(1)
|
||||
l = ll(1)
|
||||
j = jj(1)
|
||||
|
||||
if (ao_overlap_abs(k,l)*ao_overlap_abs(i,j) &
|
||||
< ao_integrals_threshold) then
|
||||
cycle
|
||||
endif
|
||||
local_threshold = ao_bielec_integral_schwartz(k,l)*ao_bielec_integral_schwartz(i,j)
|
||||
if (local_threshold < ao_integrals_threshold) then
|
||||
cycle
|
||||
endif
|
||||
i0 = i
|
||||
j0 = j
|
||||
k0 = k
|
||||
l0 = l
|
||||
values(1) = 0.d0
|
||||
local_threshold = ao_integrals_threshold/local_threshold
|
||||
do k2=1,8
|
||||
if (kk(k2)==0) then
|
||||
cycle
|
||||
endif
|
||||
i = ii(k2)
|
||||
j = jj(k2)
|
||||
k = kk(k2)
|
||||
l = ll(k2)
|
||||
c0 = HF_density_matrix_ao_alpha(k,l)+HF_density_matrix_ao_beta(k,l)
|
||||
c1 = HF_density_matrix_ao_alpha(k,i)
|
||||
c2 = HF_density_matrix_ao_beta(k,i)
|
||||
if ( dabs(c0)+dabs(c1)+dabs(c2) < local_threshold) then
|
||||
cycle
|
||||
endif
|
||||
if (values(1) == 0.d0) then
|
||||
values(1) = ao_bielec_integral(k0,l0,i0,j0)
|
||||
endif
|
||||
integral = c0 * values(1)
|
||||
ao_bi_elec_integral_alpha_tmp(i,j) += integral
|
||||
ao_bi_elec_integral_beta_tmp (i,j) += integral
|
||||
integral = values(1)
|
||||
ao_bi_elec_integral_alpha_tmp(l,j) -= c1 * integral
|
||||
ao_bi_elec_integral_beta_tmp (l,j) -= c2 * integral
|
||||
enddo
|
||||
enddo
|
||||
!$OMP END DO NOWAIT
|
||||
!$OMP CRITICAL
|
||||
ao_bi_elec_integral_alpha += ao_bi_elec_integral_alpha_tmp
|
||||
!$OMP END CRITICAL
|
||||
!$OMP CRITICAL
|
||||
ao_bi_elec_integral_beta += ao_bi_elec_integral_beta_tmp
|
||||
!$OMP END CRITICAL
|
||||
deallocate(keys,values,ao_bi_elec_integral_alpha_tmp,ao_bi_elec_integral_beta_tmp)
|
||||
!$OMP END PARALLEL
|
||||
else
|
||||
PROVIDE ao_bielec_integrals_in_map
|
||||
|
||||
integer(omp_lock_kind) :: lck(ao_num)
|
||||
integer*8 :: i8
|
||||
integer :: ii(8), jj(8), kk(8), ll(8), k2
|
||||
integer(cache_map_size_kind) :: n_elements_max, n_elements
|
||||
integer(key_kind), allocatable :: keys(:)
|
||||
double precision, allocatable :: values(:)
|
||||
|
||||
!$OMP PARALLEL DEFAULT(NONE) &
|
||||
!$OMP PRIVATE(i,j,l,k1,k,integral,ii,jj,kk,ll,i8,keys,values,n_elements_max, &
|
||||
!$OMP n_elements,ao_bi_elec_integral_alpha_tmp,ao_bi_elec_integral_beta_tmp)&
|
||||
!$OMP SHARED(ao_num,HF_density_matrix_ao_alpha,HF_density_matrix_ao_beta,&
|
||||
!$OMP ao_integrals_map, ao_bi_elec_integral_alpha, ao_bi_elec_integral_beta)
|
||||
|
||||
call get_cache_map_n_elements_max(ao_integrals_map,n_elements_max)
|
||||
allocate(keys(n_elements_max), values(n_elements_max))
|
||||
allocate(ao_bi_elec_integral_alpha_tmp(ao_num,ao_num), &
|
||||
ao_bi_elec_integral_beta_tmp(ao_num,ao_num))
|
||||
ao_bi_elec_integral_alpha_tmp = 0.d0
|
||||
ao_bi_elec_integral_beta_tmp = 0.d0
|
||||
|
||||
!$OMP DO SCHEDULE(dynamic,64)
|
||||
!DIR$ NOVECTOR
|
||||
do i8=0_8,ao_integrals_map%map_size
|
||||
n_elements = n_elements_max
|
||||
call get_cache_map(ao_integrals_map,i8,keys,values,n_elements)
|
||||
do k1=1,n_elements
|
||||
call bielec_integrals_index_reverse(kk,ii,ll,jj,keys(k1))
|
||||
|
||||
do k2=1,8
|
||||
if (kk(k2)==0) then
|
||||
cycle
|
||||
endif
|
||||
i = ii(k2)
|
||||
j = jj(k2)
|
||||
k = kk(k2)
|
||||
l = ll(k2)
|
||||
integral = (HF_density_matrix_ao_alpha(k,l)+HF_density_matrix_ao_beta(k,l)) * values(k1)
|
||||
ao_bi_elec_integral_alpha_tmp(i,j) += integral
|
||||
ao_bi_elec_integral_beta_tmp (i,j) += integral
|
||||
integral = values(k1)
|
||||
ao_bi_elec_integral_alpha_tmp(l,j) -= HF_density_matrix_ao_alpha(k,i) * integral
|
||||
ao_bi_elec_integral_beta_tmp (l,j) -= HF_density_matrix_ao_beta (k,i) * integral
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
!$OMP END DO NOWAIT
|
||||
!$OMP CRITICAL
|
||||
ao_bi_elec_integral_alpha += ao_bi_elec_integral_alpha_tmp
|
||||
!$OMP END CRITICAL
|
||||
!$OMP CRITICAL
|
||||
ao_bi_elec_integral_beta += ao_bi_elec_integral_beta_tmp
|
||||
!$OMP END CRITICAL
|
||||
deallocate(keys,values,ao_bi_elec_integral_alpha_tmp,ao_bi_elec_integral_beta_tmp)
|
||||
!$OMP END PARALLEL
|
||||
|
||||
endif
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ double precision, Fock_matrix_mo_alpha, (mo_tot_num,mo_tot_num) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Fock matrix on the MO basis
|
||||
END_DOC
|
||||
call ao_to_mo(Fock_matrix_ao_alpha,size(Fock_matrix_ao_alpha,1), &
|
||||
Fock_matrix_mo_alpha,size(Fock_matrix_mo_alpha,1))
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
BEGIN_PROVIDER [ double precision, Fock_matrix_mo_beta, (mo_tot_num,mo_tot_num) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Fock matrix on the MO basis
|
||||
END_DOC
|
||||
call ao_to_mo(Fock_matrix_ao_beta,size(Fock_matrix_ao_beta,1), &
|
||||
Fock_matrix_mo_beta,size(Fock_matrix_mo_beta,1))
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ double precision, HF_energy ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Hartree-Fock energy
|
||||
END_DOC
|
||||
HF_energy = nuclear_repulsion
|
||||
|
||||
integer :: i,j
|
||||
do j=1,ao_num
|
||||
do i=1,ao_num
|
||||
HF_energy += 0.5d0 * ( &
|
||||
(ao_mono_elec_integral(i,j) + Fock_matrix_ao_alpha(i,j) ) * HF_density_matrix_ao_alpha(i,j) +&
|
||||
(ao_mono_elec_integral(i,j) + Fock_matrix_ao_beta (i,j) ) * HF_density_matrix_ao_beta (i,j) )
|
||||
enddo
|
||||
enddo
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
BEGIN_PROVIDER [ double precision, Fock_matrix_ao, (ao_num, ao_num) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Fock matrix in AO basis set
|
||||
END_DOC
|
||||
|
||||
if ( (elec_alpha_num == elec_beta_num).and. &
|
||||
(level_shift == 0.) ) &
|
||||
then
|
||||
integer :: i,j
|
||||
do j=1,ao_num
|
||||
do i=1,ao_num
|
||||
Fock_matrix_ao(i,j) = Fock_matrix_ao_alpha(i,j)
|
||||
enddo
|
||||
enddo
|
||||
else
|
||||
call mo_to_ao(Fock_matrix_mo,size(Fock_matrix_mo,1), &
|
||||
Fock_matrix_ao,size(Fock_matrix_ao,1))
|
||||
endif
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
|
|
@ -1,41 +0,0 @@
|
|||
BEGIN_PROVIDER [double precision, HF_density_matrix_ao_alpha, (ao_num,ao_num) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! S^{-1}.P_alpha.S^{-1}
|
||||
END_DOC
|
||||
|
||||
call dgemm('N','T',ao_num,ao_num,elec_alpha_num,1.d0, &
|
||||
mo_coef, size(mo_coef,1), &
|
||||
mo_coef, size(mo_coef,1), 0.d0, &
|
||||
HF_density_matrix_ao_alpha, size(HF_density_matrix_ao_alpha,1))
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ double precision, HF_density_matrix_ao_beta, (ao_num,ao_num) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! S^{-1}.P_beta.S^{-1}
|
||||
END_DOC
|
||||
|
||||
call dgemm('N','T',ao_num,ao_num,elec_beta_num,1.d0, &
|
||||
mo_coef, size(mo_coef,1), &
|
||||
mo_coef, size(mo_coef,1), 0.d0, &
|
||||
HF_density_matrix_ao_beta, size(HF_density_matrix_ao_beta,1))
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ double precision, HF_density_matrix_ao, (ao_num,ao_num) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! S^{-1}.P.S^{-1} where P = C.C^t
|
||||
END_DOC
|
||||
ASSERT (size(HF_density_matrix_ao,1) == size(HF_density_matrix_ao_alpha,1))
|
||||
if (elec_alpha_num== elec_beta_num) then
|
||||
HF_density_matrix_ao = HF_density_matrix_ao_alpha + HF_density_matrix_ao_alpha
|
||||
else
|
||||
ASSERT (size(HF_density_matrix_ao,1) == size(HF_density_matrix_ao_beta ,1))
|
||||
HF_density_matrix_ao = HF_density_matrix_ao_alpha + HF_density_matrix_ao_beta
|
||||
endif
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
|
|
@ -1,4 +0,0 @@
|
|||
Integrals_Bielec
|
||||
AO_one_e_integrals
|
||||
MOGuess
|
||||
Bitmask
|
||||
|
|
@ -1,33 +0,0 @@
|
|||
============
|
||||
Hartree-Fock
|
||||
============
|
||||
|
||||
|
||||
The Hartree-Fock module performs *Restricted* Hartree-Fock calculations (the
|
||||
spatial part of the |MOs| is common for alpha and beta spinorbitals).
|
||||
|
||||
The Hartree-Fock program does the following:
|
||||
|
||||
#. Compute/Read all the one- and two-electron integrals, and store them in memory
|
||||
#. Check in the |EZFIO| database if there is a set of |MOs|. If there is, it
|
||||
will read them as initial guess. Otherwise, it will create a guess.
|
||||
#. Perform the |SCF| iterations
|
||||
|
||||
At each iteration, the |MOs| are saved in the |EZFIO| database. Hence, if the calculation
|
||||
crashes for any unexpected reason, the calculation can be restarted by running again
|
||||
the |SCF| with the same |EZFIO| database.
|
||||
|
||||
The `DIIS`_ algorithm is implemented, as well as the `level-shifting`_ method.
|
||||
If the |SCF| does not converge, try again with a higher value of :option:`level_shift`.
|
||||
|
||||
To start a calculation from scratch, the simplest way is to remove the
|
||||
``mo_basis`` directory from the |EZFIO| database, and run the |SCF| again.
|
||||
|
||||
|
||||
|
||||
|
||||
.. _DIIS: https://en.wikipedia.org/w/index.php?title=DIIS
|
||||
.. _level-shifting: https://doi.org/10.1002/qua.560070407
|
||||
|
||||
|
||||
|
||||
|
|
@ -1,289 +0,0 @@
|
|||
subroutine Roothaan_Hall_SCF
|
||||
|
||||
BEGIN_DOC
|
||||
! Roothaan-Hall algorithm for SCF Hartree-Fock calculation
|
||||
END_DOC
|
||||
|
||||
implicit none
|
||||
|
||||
double precision :: energy_SCF,energy_SCF_previous,Delta_energy_SCF
|
||||
double precision :: max_error_DIIS,max_error_DIIS_alpha,max_error_DIIS_beta
|
||||
double precision, allocatable :: Fock_matrix_DIIS(:,:,:),error_matrix_DIIS(:,:,:)
|
||||
|
||||
integer :: iteration_SCF,dim_DIIS,index_dim_DIIS
|
||||
|
||||
integer :: i,j
|
||||
double precision, allocatable :: mo_coef_save(:,:)
|
||||
|
||||
PROVIDE ao_md5 mo_occ level_shift
|
||||
|
||||
allocate(mo_coef_save(ao_num,mo_tot_num), &
|
||||
Fock_matrix_DIIS (ao_num,ao_num,max_dim_DIIS), &
|
||||
error_matrix_DIIS(ao_num,ao_num,max_dim_DIIS) &
|
||||
)
|
||||
|
||||
call write_time(6)
|
||||
|
||||
write(6,'(A4, 1X, A16, 1X, A16, 1X, A16)') &
|
||||
'====','================','================','================'
|
||||
write(6,'(A4, 1X, A16, 1X, A16, 1X, A16)') &
|
||||
' N ', 'Energy ', 'Energy diff ', 'DIIS error '
|
||||
write(6,'(A4, 1X, A16, 1X, A16, 1X, A16)') &
|
||||
'====','================','================','================'
|
||||
|
||||
! Initialize energies and density matrices
|
||||
|
||||
energy_SCF_previous = HF_energy
|
||||
Delta_energy_SCF = 1.d0
|
||||
iteration_SCF = 0
|
||||
dim_DIIS = 0
|
||||
max_error_DIIS = 1.d0
|
||||
|
||||
|
||||
!
|
||||
! Start of main SCF loop
|
||||
!
|
||||
do while(( (max_error_DIIS > threshold_DIIS_nonzero).or.(dabs(Delta_energy_SCF) > thresh_SCF) ) .and. (iteration_SCF < n_it_SCF_max))
|
||||
|
||||
! Increment cycle number
|
||||
|
||||
iteration_SCF += 1
|
||||
|
||||
! Current size of the DIIS space
|
||||
|
||||
dim_DIIS = min(dim_DIIS+1,max_dim_DIIS)
|
||||
|
||||
if (scf_algorithm == 'DIIS') then
|
||||
|
||||
! Store Fock and error matrices at each iteration
|
||||
do j=1,ao_num
|
||||
do i=1,ao_num
|
||||
index_dim_DIIS = mod(dim_DIIS-1,max_dim_DIIS)+1
|
||||
Fock_matrix_DIIS (i,j,index_dim_DIIS) = Fock_matrix_AO(i,j)
|
||||
error_matrix_DIIS(i,j,index_dim_DIIS) = FPS_SPF_matrix_AO(i,j)
|
||||
enddo
|
||||
enddo
|
||||
|
||||
! Compute the extrapolated Fock matrix
|
||||
|
||||
call extrapolate_Fock_matrix( &
|
||||
error_matrix_DIIS,Fock_matrix_DIIS, &
|
||||
Fock_matrix_AO,size(Fock_matrix_AO,1), &
|
||||
iteration_SCF,dim_DIIS &
|
||||
)
|
||||
|
||||
Fock_matrix_AO_alpha = Fock_matrix_AO*0.5d0
|
||||
Fock_matrix_AO_beta = Fock_matrix_AO*0.5d0
|
||||
TOUCH Fock_matrix_AO_alpha Fock_matrix_AO_beta
|
||||
|
||||
endif
|
||||
|
||||
MO_coef = eigenvectors_Fock_matrix_MO
|
||||
|
||||
TOUCH MO_coef
|
||||
|
||||
! Calculate error vectors
|
||||
|
||||
max_error_DIIS = maxval(Abs(FPS_SPF_Matrix_MO))
|
||||
|
||||
! SCF energy
|
||||
|
||||
energy_SCF = HF_energy
|
||||
Delta_Energy_SCF = energy_SCF - energy_SCF_previous
|
||||
if ( (SCF_algorithm == 'DIIS').and.(Delta_Energy_SCF > 0.d0) ) then
|
||||
Fock_matrix_AO(1:ao_num,1:ao_num) = Fock_matrix_DIIS (1:ao_num,1:ao_num,index_dim_DIIS)
|
||||
Fock_matrix_AO_alpha = Fock_matrix_AO*0.5d0
|
||||
Fock_matrix_AO_beta = Fock_matrix_AO*0.5d0
|
||||
TOUCH Fock_matrix_AO_alpha Fock_matrix_AO_beta
|
||||
endif
|
||||
|
||||
double precision :: level_shift_save
|
||||
level_shift_save = level_shift
|
||||
mo_coef_save(1:ao_num,1:mo_tot_num) = mo_coef(1:ao_num,1:mo_tot_num)
|
||||
do while (Delta_Energy_SCF > 0.d0)
|
||||
mo_coef(1:ao_num,1:mo_tot_num) = mo_coef_save
|
||||
TOUCH mo_coef
|
||||
level_shift = level_shift + 1.0d0
|
||||
mo_coef(1:ao_num,1:mo_tot_num) = eigenvectors_Fock_matrix_MO(1:ao_num,1:mo_tot_num)
|
||||
TOUCH mo_coef level_shift
|
||||
Delta_Energy_SCF = HF_energy - energy_SCF_previous
|
||||
energy_SCF = HF_energy
|
||||
if (level_shift-level_shift_save > 50.d0) then
|
||||
level_shift = level_shift_save
|
||||
SOFT_TOUCH level_shift
|
||||
exit
|
||||
endif
|
||||
dim_DIIS=0
|
||||
enddo
|
||||
level_shift = level_shift * 0.75d0
|
||||
SOFT_TOUCH level_shift
|
||||
energy_SCF_previous = energy_SCF
|
||||
|
||||
! Print results at the end of each iteration
|
||||
|
||||
write(6,'(I4, 1X, F16.10, 1X, F16.10, 1X, F16.10, 1X, I3)') &
|
||||
iteration_SCF, energy_SCF, Delta_energy_SCF, max_error_DIIS, dim_DIIS
|
||||
|
||||
if (Delta_energy_SCF < 0.d0) then
|
||||
call save_mos
|
||||
endif
|
||||
|
||||
enddo
|
||||
|
||||
!
|
||||
! End of Main SCF loop
|
||||
!
|
||||
|
||||
write(6,'(A4, 1X, A16, 1X, A16, 1X, A16)') &
|
||||
'====','================','================','================'
|
||||
write(6,*)
|
||||
|
||||
if(.not.no_oa_or_av_opt)then
|
||||
call mo_as_eigvectors_of_mo_matrix(Fock_matrix_mo,size(Fock_matrix_mo,1),size(Fock_matrix_mo,2),mo_label,1,.true.)
|
||||
endif
|
||||
|
||||
call write_double(6, Energy_SCF, 'Hartree-Fock energy')
|
||||
call ezfio_set_hartree_fock_energy(Energy_SCF)
|
||||
|
||||
call write_time(6)
|
||||
|
||||
end
|
||||
|
||||
subroutine extrapolate_Fock_matrix( &
|
||||
error_matrix_DIIS,Fock_matrix_DIIS, &
|
||||
Fock_matrix_AO_,size_Fock_matrix_AO, &
|
||||
iteration_SCF,dim_DIIS &
|
||||
)
|
||||
|
||||
BEGIN_DOC
|
||||
! Compute the extrapolated Fock matrix using the DIIS procedure
|
||||
END_DOC
|
||||
|
||||
implicit none
|
||||
|
||||
double precision,intent(in) :: Fock_matrix_DIIS(ao_num,ao_num,*),error_matrix_DIIS(ao_num,ao_num,*)
|
||||
integer,intent(in) :: iteration_SCF, size_Fock_matrix_AO
|
||||
double precision,intent(inout):: Fock_matrix_AO_(size_Fock_matrix_AO,ao_num)
|
||||
integer,intent(inout) :: dim_DIIS
|
||||
|
||||
double precision,allocatable :: B_matrix_DIIS(:,:),X_vector_DIIS(:)
|
||||
double precision,allocatable :: C_vector_DIIS(:)
|
||||
|
||||
double precision,allocatable :: scratch(:,:)
|
||||
integer :: i,j,k,i_DIIS,j_DIIS
|
||||
|
||||
allocate( &
|
||||
B_matrix_DIIS(dim_DIIS+1,dim_DIIS+1), &
|
||||
X_vector_DIIS(dim_DIIS+1), &
|
||||
C_vector_DIIS(dim_DIIS+1), &
|
||||
scratch(ao_num,ao_num) &
|
||||
)
|
||||
|
||||
! Compute the matrices B and X
|
||||
do j=1,dim_DIIS
|
||||
do i=1,dim_DIIS
|
||||
|
||||
j_DIIS = mod(iteration_SCF-j,max_dim_DIIS)+1
|
||||
i_DIIS = mod(iteration_SCF-i,max_dim_DIIS)+1
|
||||
|
||||
! Compute product of two errors vectors
|
||||
|
||||
call dgemm('N','N',ao_num,ao_num,ao_num, &
|
||||
1.d0, &
|
||||
error_matrix_DIIS(1,1,i_DIIS),size(error_matrix_DIIS,1), &
|
||||
error_matrix_DIIS(1,1,j_DIIS),size(error_matrix_DIIS,1), &
|
||||
0.d0, &
|
||||
scratch,size(scratch,1))
|
||||
|
||||
! Compute Trace
|
||||
|
||||
B_matrix_DIIS(i,j) = 0.d0
|
||||
do k=1,ao_num
|
||||
B_matrix_DIIS(i,j) = B_matrix_DIIS(i,j) + scratch(k,k)
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
|
||||
! Pad B matrix and build the X matrix
|
||||
|
||||
do i=1,dim_DIIS
|
||||
B_matrix_DIIS(i,dim_DIIS+1) = -1.d0
|
||||
B_matrix_DIIS(dim_DIIS+1,i) = -1.d0
|
||||
C_vector_DIIS(i) = 0.d0
|
||||
enddo
|
||||
B_matrix_DIIS(dim_DIIS+1,dim_DIIS+1) = 0.d0
|
||||
C_vector_DIIS(dim_DIIS+1) = -1.d0
|
||||
|
||||
! Solve the linear system C = B.X
|
||||
|
||||
integer :: info
|
||||
integer,allocatable :: ipiv(:)
|
||||
|
||||
allocate( &
|
||||
ipiv(dim_DIIS+1) &
|
||||
)
|
||||
|
||||
double precision, allocatable :: AF(:,:)
|
||||
allocate (AF(dim_DIIS+1,dim_DIIS+1))
|
||||
double precision :: rcond, ferr, berr
|
||||
integer :: iwork(dim_DIIS+1), lwork
|
||||
|
||||
call dsysvx('N','U',dim_DIIS+1,1, &
|
||||
B_matrix_DIIS,size(B_matrix_DIIS,1), &
|
||||
AF, size(AF,1), &
|
||||
ipiv, &
|
||||
C_vector_DIIS,size(C_vector_DIIS,1), &
|
||||
X_vector_DIIS,size(X_vector_DIIS,1), &
|
||||
rcond, &
|
||||
ferr, &
|
||||
berr, &
|
||||
scratch,-1, &
|
||||
iwork, &
|
||||
info &
|
||||
)
|
||||
lwork = int(scratch(1,1))
|
||||
deallocate(scratch)
|
||||
allocate(scratch(lwork,1))
|
||||
|
||||
call dsysvx('N','U',dim_DIIS+1,1, &
|
||||
B_matrix_DIIS,size(B_matrix_DIIS,1), &
|
||||
AF, size(AF,1), &
|
||||
ipiv, &
|
||||
C_vector_DIIS,size(C_vector_DIIS,1), &
|
||||
X_vector_DIIS,size(X_vector_DIIS,1), &
|
||||
rcond, &
|
||||
ferr, &
|
||||
berr, &
|
||||
scratch,size(scratch), &
|
||||
iwork, &
|
||||
info &
|
||||
)
|
||||
|
||||
if(info < 0) then
|
||||
stop 'bug in DIIS'
|
||||
endif
|
||||
|
||||
if (rcond > 1.d-12) then
|
||||
|
||||
! Compute extrapolated Fock matrix
|
||||
|
||||
|
||||
!$OMP PARALLEL DO PRIVATE(i,j,k) DEFAULT(SHARED)
|
||||
do j=1,ao_num
|
||||
do i=1,ao_num
|
||||
Fock_matrix_AO_(i,j) = 0.d0
|
||||
enddo
|
||||
do k=1,dim_DIIS
|
||||
do i=1,ao_num
|
||||
Fock_matrix_AO_(i,j) = Fock_matrix_AO_(i,j) + &
|
||||
X_vector_DIIS(k)*Fock_matrix_DIIS(i,j,dim_DIIS-k+1)
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
!$OMP END PARALLEL DO
|
||||
|
||||
else
|
||||
dim_DIIS = 0
|
||||
endif
|
||||
|
||||
end
|
||||
|
|
@ -1,60 +0,0 @@
|
|||
program scf
|
||||
BEGIN_DOC
|
||||
! Produce `Hartree_Fock` MO orbital
|
||||
! output: mo_basis.mo_tot_num mo_basis.mo_label mo_basis.ao_md5 mo_basis.mo_coef mo_basis.mo_occ
|
||||
! output: hartree_fock.energy
|
||||
! optional: mo_basis.mo_coef
|
||||
END_DOC
|
||||
call create_guess
|
||||
call orthonormalize_mos
|
||||
call run
|
||||
end
|
||||
|
||||
subroutine create_guess
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Create a MO guess if no MOs are present in the EZFIO directory
|
||||
END_DOC
|
||||
logical :: exists
|
||||
PROVIDE ezfio_filename
|
||||
call ezfio_has_mo_basis_mo_coef(exists)
|
||||
if (.not.exists) then
|
||||
if (mo_guess_type == "HCore") then
|
||||
mo_coef = ao_ortho_lowdin_coef
|
||||
TOUCH mo_coef
|
||||
mo_label = 'Guess'
|
||||
call mo_as_eigvectors_of_mo_matrix(mo_mono_elec_integral,size(mo_mono_elec_integral,1),size(mo_mono_elec_integral,2),mo_label,1,.false.)
|
||||
SOFT_TOUCH mo_coef mo_label
|
||||
else if (mo_guess_type == "Huckel") then
|
||||
call huckel_guess
|
||||
else
|
||||
print *, 'Unrecognized MO guess type : '//mo_guess_type
|
||||
stop 1
|
||||
endif
|
||||
endif
|
||||
end
|
||||
|
||||
subroutine run
|
||||
|
||||
BEGIN_DOC
|
||||
! Run SCF calculation
|
||||
END_DOC
|
||||
|
||||
use bitmasks
|
||||
implicit none
|
||||
|
||||
double precision :: SCF_energy_before,SCF_energy_after,diag_H_mat_elem
|
||||
double precision :: EHF
|
||||
integer :: i_it, i, j, k
|
||||
|
||||
EHF = HF_energy
|
||||
|
||||
mo_label = "Canonical"
|
||||
|
||||
! Choose SCF algorithm
|
||||
|
||||
call Roothaan_Hall_SCF
|
||||
|
||||
end
|
||||
|
||||
|
||||
|
|
@ -1,100 +0,0 @@
|
|||
BEGIN_PROVIDER [ double precision, eigenvectors_Fock_matrix_mo, (ao_num,mo_tot_num) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Eigenvector of the Fock matrix in the MO basis obtained with level shift.
|
||||
END_DOC
|
||||
|
||||
integer :: i,j
|
||||
integer :: liwork, lwork, n, info
|
||||
integer, allocatable :: iwork(:)
|
||||
double precision, allocatable :: work(:), F(:,:), S(:,:)
|
||||
double precision, allocatable :: diag(:)
|
||||
|
||||
|
||||
allocate( F(mo_tot_num,mo_tot_num) )
|
||||
allocate (diag(mo_tot_num) )
|
||||
|
||||
do j=1,mo_tot_num
|
||||
do i=1,mo_tot_num
|
||||
F(i,j) = Fock_matrix_mo(i,j)
|
||||
enddo
|
||||
enddo
|
||||
if(no_oa_or_av_opt)then
|
||||
integer :: iorb,jorb
|
||||
do i = 1, n_act_orb
|
||||
iorb = list_act(i)
|
||||
do j = 1, n_inact_orb
|
||||
jorb = list_inact(j)
|
||||
F(iorb,jorb) = 0.d0
|
||||
F(jorb,iorb) = 0.d0
|
||||
enddo
|
||||
do j = 1, n_virt_orb
|
||||
jorb = list_virt(j)
|
||||
F(iorb,jorb) = 0.d0
|
||||
F(jorb,iorb) = 0.d0
|
||||
enddo
|
||||
do j = 1, n_core_orb
|
||||
jorb = list_core(j)
|
||||
F(iorb,jorb) = 0.d0
|
||||
F(jorb,iorb) = 0.d0
|
||||
enddo
|
||||
enddo
|
||||
endif
|
||||
|
||||
|
||||
! Insert level shift here
|
||||
do i = elec_beta_num+1, elec_alpha_num
|
||||
F(i,i) += 0.5d0*level_shift
|
||||
enddo
|
||||
|
||||
do i = elec_alpha_num+1, mo_tot_num
|
||||
F(i,i) += level_shift
|
||||
enddo
|
||||
|
||||
n = mo_tot_num
|
||||
lwork = 1+6*n + 2*n*n
|
||||
liwork = 3 + 5*n
|
||||
|
||||
allocate(work(lwork))
|
||||
allocate(iwork(liwork) )
|
||||
|
||||
lwork = -1
|
||||
liwork = -1
|
||||
|
||||
call dsyevd( 'V', 'U', mo_tot_num, F, &
|
||||
size(F,1), diag, work, lwork, iwork, liwork, info)
|
||||
|
||||
if (info /= 0) then
|
||||
print *, irp_here//' DSYEVD failed : ', info
|
||||
stop 1
|
||||
endif
|
||||
lwork = int(work(1))
|
||||
liwork = iwork(1)
|
||||
deallocate(iwork)
|
||||
deallocate(work)
|
||||
|
||||
allocate(work(lwork))
|
||||
allocate(iwork(liwork) )
|
||||
call dsyevd( 'V', 'U', mo_tot_num, F, &
|
||||
size(F,1), diag, work, lwork, iwork, liwork, info)
|
||||
deallocate(iwork)
|
||||
|
||||
|
||||
if (info /= 0) then
|
||||
call dsyev( 'V', 'L', mo_tot_num, F, &
|
||||
size(F,1), diag, work, lwork, info)
|
||||
|
||||
if (info /= 0) then
|
||||
print *, irp_here//' DSYEV failed : ', info
|
||||
stop 1
|
||||
endif
|
||||
endif
|
||||
|
||||
call dgemm('N','N',ao_num,mo_tot_num,mo_tot_num, 1.d0, &
|
||||
mo_coef, size(mo_coef,1), F, size(F,1), &
|
||||
0.d0, eigenvectors_Fock_matrix_mo, size(eigenvectors_Fock_matrix_mo,1))
|
||||
deallocate(work, F, diag)
|
||||
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
|
|
@ -1,34 +0,0 @@
|
|||
subroutine huckel_guess
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Build the MOs using the extended Huckel model
|
||||
END_DOC
|
||||
integer :: i,j
|
||||
double precision :: accu
|
||||
double precision :: c
|
||||
character*(64) :: label
|
||||
double precision, allocatable :: A(:,:)
|
||||
label = "Guess"
|
||||
c = 0.5d0 * 1.75d0
|
||||
|
||||
allocate (A(ao_num, ao_num))
|
||||
A = 0.d0
|
||||
do j=1,ao_num
|
||||
do i=1,ao_num
|
||||
A(i,j) = c * ao_overlap(i,j) * (ao_mono_elec_integral_diag(i) + ao_mono_elec_integral_diag(j))
|
||||
enddo
|
||||
A(j,j) = ao_mono_elec_integral_diag(j) + ao_bi_elec_integral_alpha(j,j)
|
||||
enddo
|
||||
|
||||
Fock_matrix_ao_alpha(1:ao_num,1:ao_num) = A(1:ao_num,1:ao_num)
|
||||
Fock_matrix_ao_beta (1:ao_num,1:ao_num) = A(1:ao_num,1:ao_num)
|
||||
|
||||
! TOUCH mo_coef
|
||||
|
||||
TOUCH Fock_matrix_ao_alpha Fock_matrix_ao_beta
|
||||
mo_coef = eigenvectors_fock_matrix_mo
|
||||
SOFT_TOUCH mo_coef
|
||||
call save_mos
|
||||
deallocate(A)
|
||||
|
||||
end
|
||||
|
|
@ -1,55 +0,0 @@
|
|||
[disk_access_mo_integrals]
|
||||
type: Disk_access
|
||||
doc: Read/Write |MO| integrals from/to disk [ Write | Read | None ]
|
||||
interface: ezfio,provider,ocaml
|
||||
default: None
|
||||
|
||||
[disk_access_ao_integrals]
|
||||
type: Disk_access
|
||||
doc: Read/Write |AO| integrals from/to disk [ Write | Read | None ]
|
||||
interface: ezfio,provider,ocaml
|
||||
default: None
|
||||
|
||||
[ao_integrals_threshold]
|
||||
type: Threshold
|
||||
doc: If | (pq|rs) | < `ao_integrals_threshold` then (pq|rs) is zero
|
||||
interface: ezfio,provider,ocaml
|
||||
default: 1.e-15
|
||||
ezfio_name: threshold_ao
|
||||
|
||||
[mo_integrals_threshold]
|
||||
type: Threshold
|
||||
doc: If | <ij|kl> | < `mo_integrals_threshold` then <ij|kl> is zero
|
||||
interface: ezfio,provider,ocaml
|
||||
default: 1.e-15
|
||||
ezfio_name: threshold_mo
|
||||
|
||||
|
||||
[no_vvvv_integrals]
|
||||
type: logical
|
||||
doc: If `True`, computes all integrals except for the integrals having 4 virtual indices
|
||||
interface: ezfio,provider,ocaml
|
||||
default: False
|
||||
ezfio_name: no_vvvv_integrals
|
||||
|
||||
[no_ivvv_integrals]
|
||||
type: logical
|
||||
doc: Can be switched on only if `no_vvvv_integrals` is `True`, then does not compute the integrals with 3 virtual indices and 1 belonging to the core inactive active orbitals
|
||||
interface: ezfio,provider,ocaml
|
||||
default: False
|
||||
ezfio_name: no_ivvv_integrals
|
||||
|
||||
[no_vvv_integrals]
|
||||
type: logical
|
||||
doc: Can be switched on only if `no_vvvv_integrals` is `True`, then does not compute the integrals with 3 virtual orbitals
|
||||
interface: ezfio,provider,ocaml
|
||||
default: False
|
||||
ezfio_name: no_vvv_integrals
|
||||
|
||||
[do_direct_integrals]
|
||||
type: logical
|
||||
doc: Compute integrals on the fly (very slow, only for debugging)
|
||||
interface: ezfio,provider,ocaml
|
||||
default: False
|
||||
ezfio_name: direct
|
||||
|
||||
|
|
@ -1,7 +0,0 @@
|
|||
AO_one_e_integrals
|
||||
MO_one_e_integrals
|
||||
Pseudo
|
||||
Bitmask
|
||||
ZMQ
|
||||
AO_Basis
|
||||
MO_Basis
|
||||
|
|
@ -1,21 +0,0 @@
|
|||
================
|
||||
Integrals_Bielec
|
||||
================
|
||||
|
||||
Here, all two-electron integrals (:math:`1/r_{12}`) are computed.
|
||||
As they have 4 indices and many are zero, they are stored in a map, as defined
|
||||
in :file:`Utils/map_module.f90`.
|
||||
|
||||
To fetch an |AO| integral, use the
|
||||
`get_ao_bielec_integral(i,j,k,l,ao_integrals_map)` function, and
|
||||
to fetch an |MO| integral, use
|
||||
`get_mo_bielec_integral(i,j,k,l,mo_integrals_map)` or
|
||||
`mo_bielec_integral(i,j,k,l)`.
|
||||
|
||||
The conventions are:
|
||||
|
||||
* For |AO| integrals : (ik|jl) = (11|22)
|
||||
* For |MO| integrals : <ij|kl> = <12|12>
|
||||
|
||||
|
||||
|
||||
File diff suppressed because it is too large
Load Diff
|
|
@ -1,244 +0,0 @@
|
|||
subroutine ao_bielec_integrals_in_map_slave_tcp(i)
|
||||
implicit none
|
||||
integer, intent(in) :: i
|
||||
BEGIN_DOC
|
||||
! Computes a buffer of integrals. i is the ID of the current thread.
|
||||
END_DOC
|
||||
call ao_bielec_integrals_in_map_slave(0,i)
|
||||
end
|
||||
|
||||
|
||||
subroutine ao_bielec_integrals_in_map_slave_inproc(i)
|
||||
implicit none
|
||||
integer, intent(in) :: i
|
||||
BEGIN_DOC
|
||||
! Computes a buffer of integrals. i is the ID of the current thread.
|
||||
END_DOC
|
||||
call ao_bielec_integrals_in_map_slave(1,i)
|
||||
end
|
||||
|
||||
|
||||
subroutine push_integrals(zmq_socket_push, n_integrals, buffer_i, buffer_value, task_id)
|
||||
use f77_zmq
|
||||
use map_module
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Push integrals in the push socket
|
||||
END_DOC
|
||||
integer(ZMQ_PTR), intent(in) :: zmq_socket_push
|
||||
integer, intent(in) :: n_integrals
|
||||
integer(key_kind), intent(in) :: buffer_i(*)
|
||||
real(integral_kind), intent(in) :: buffer_value(*)
|
||||
integer, intent(in) :: task_id
|
||||
integer :: rc
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, n_integrals, 4, ZMQ_SNDMORE)
|
||||
if (rc /= 4) then
|
||||
print *, irp_here, ': f77_zmq_send( zmq_socket_push, n_integrals, 4, ZMQ_SNDMORE)'
|
||||
stop 'error'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, buffer_i, key_kind*n_integrals, ZMQ_SNDMORE)
|
||||
if (rc /= key_kind*n_integrals) then
|
||||
print *, irp_here, ': f77_zmq_send( zmq_socket_push, buffer_i, key_kind*n_integrals, ZMQ_SNDMORE)'
|
||||
stop 'error'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, buffer_value, integral_kind*n_integrals, ZMQ_SNDMORE)
|
||||
if (rc /= integral_kind*n_integrals) then
|
||||
print *, irp_here, ': f77_zmq_send( zmq_socket_push, buffer_value, integral_kind*n_integrals, 0)'
|
||||
stop 'error'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_send( zmq_socket_push, task_id, 4, 0)
|
||||
if (rc /= 4) then
|
||||
print *, irp_here, ': f77_zmq_send( zmq_socket_push, task_id, 4, 0)'
|
||||
stop 'error'
|
||||
endif
|
||||
|
||||
IRP_IF ZMQ_PUSH
|
||||
IRP_ELSE
|
||||
integer :: idummy
|
||||
rc = f77_zmq_recv( zmq_socket_push, idummy, 4, 0)
|
||||
if (rc /= 4) then
|
||||
print *, irp_here, ': f77_zmq_send( zmq_socket_push, idummy, 4, 0)'
|
||||
stop 'error'
|
||||
endif
|
||||
IRP_ENDIF
|
||||
end
|
||||
|
||||
|
||||
|
||||
|
||||
|
||||
subroutine ao_bielec_integrals_in_map_slave(thread,iproc)
|
||||
use map_module
|
||||
use f77_zmq
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Computes a buffer of integrals
|
||||
END_DOC
|
||||
|
||||
integer, intent(in) :: thread, iproc
|
||||
|
||||
integer :: j,l,n_integrals
|
||||
integer :: rc
|
||||
real(integral_kind), allocatable :: buffer_value(:)
|
||||
integer(key_kind), allocatable :: buffer_i(:)
|
||||
|
||||
integer :: worker_id, task_id
|
||||
character*(512) :: task
|
||||
|
||||
integer(ZMQ_PTR),external :: new_zmq_to_qp_run_socket
|
||||
integer(ZMQ_PTR) :: zmq_to_qp_run_socket
|
||||
|
||||
integer(ZMQ_PTR), external :: new_zmq_push_socket
|
||||
integer(ZMQ_PTR) :: zmq_socket_push
|
||||
|
||||
character*(64) :: state
|
||||
|
||||
zmq_to_qp_run_socket = new_zmq_to_qp_run_socket()
|
||||
|
||||
integer, external :: connect_to_taskserver
|
||||
if (connect_to_taskserver(zmq_to_qp_run_socket,worker_id,thread) == -1) then
|
||||
call end_zmq_to_qp_run_socket(zmq_to_qp_run_socket)
|
||||
return
|
||||
endif
|
||||
|
||||
zmq_socket_push = new_zmq_push_socket(thread)
|
||||
|
||||
allocate ( buffer_i(ao_num*ao_num), buffer_value(ao_num*ao_num) )
|
||||
|
||||
|
||||
do
|
||||
integer, external :: get_task_from_taskserver
|
||||
if (get_task_from_taskserver(zmq_to_qp_run_socket,worker_id, task_id, task) == -1) then
|
||||
exit
|
||||
endif
|
||||
if (task_id == 0) exit
|
||||
read(task,*) j, l
|
||||
integer, external :: task_done_to_taskserver
|
||||
call compute_ao_integrals_jl(j,l,n_integrals,buffer_i,buffer_value)
|
||||
if (task_done_to_taskserver(zmq_to_qp_run_socket,worker_id,task_id) == -1) then
|
||||
stop 'Unable to send task_done'
|
||||
endif
|
||||
call push_integrals(zmq_socket_push, n_integrals, buffer_i, buffer_value, task_id)
|
||||
enddo
|
||||
|
||||
integer, external :: disconnect_from_taskserver
|
||||
if (disconnect_from_taskserver(zmq_to_qp_run_socket,worker_id) == -1) then
|
||||
continue
|
||||
endif
|
||||
deallocate( buffer_i, buffer_value )
|
||||
call end_zmq_to_qp_run_socket(zmq_to_qp_run_socket)
|
||||
call end_zmq_push_socket(zmq_socket_push,thread)
|
||||
|
||||
end
|
||||
|
||||
|
||||
subroutine ao_bielec_integrals_in_map_collector(zmq_socket_pull)
|
||||
use map_module
|
||||
use f77_zmq
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Collects results from the AO integral calculation
|
||||
END_DOC
|
||||
|
||||
integer(ZMQ_PTR), intent(in) :: zmq_socket_pull
|
||||
integer :: j,l,n_integrals
|
||||
integer :: rc
|
||||
|
||||
real(integral_kind), allocatable :: buffer_value(:)
|
||||
integer(key_kind), allocatable :: buffer_i(:)
|
||||
|
||||
integer(ZMQ_PTR),external :: new_zmq_to_qp_run_socket
|
||||
integer(ZMQ_PTR) :: zmq_to_qp_run_socket
|
||||
|
||||
integer(ZMQ_PTR), external :: new_zmq_pull_socket
|
||||
|
||||
integer*8 :: control, accu, sze
|
||||
integer :: task_id, more
|
||||
|
||||
zmq_to_qp_run_socket = new_zmq_to_qp_run_socket()
|
||||
|
||||
sze = ao_num*ao_num
|
||||
allocate ( buffer_i(sze), buffer_value(sze) )
|
||||
|
||||
accu = 0_8
|
||||
more = 1
|
||||
do while (more == 1)
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, n_integrals, 4, 0)
|
||||
if (rc == -1) then
|
||||
n_integrals = 0
|
||||
return
|
||||
endif
|
||||
if (rc /= 4) then
|
||||
print *, irp_here, ': f77_zmq_recv( zmq_socket_pull, n_integrals, 4, 0)'
|
||||
stop 'error'
|
||||
endif
|
||||
|
||||
if (n_integrals >= 0) then
|
||||
|
||||
if (n_integrals > sze) then
|
||||
deallocate (buffer_value, buffer_i)
|
||||
sze = n_integrals
|
||||
allocate (buffer_value(sze), buffer_i(sze))
|
||||
endif
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, buffer_i, key_kind*n_integrals, 0)
|
||||
if (rc /= key_kind*n_integrals) then
|
||||
print *, rc, key_kind, n_integrals
|
||||
print *, irp_here, ': f77_zmq_recv( zmq_socket_pull, buffer_i, key_kind*n_integrals, 0)'
|
||||
stop 'error'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, buffer_value, integral_kind*n_integrals, 0)
|
||||
if (rc /= integral_kind*n_integrals) then
|
||||
print *, irp_here, ': f77_zmq_recv( zmq_socket_pull, buffer_value, integral_kind*n_integrals, 0)'
|
||||
stop 'error'
|
||||
endif
|
||||
|
||||
rc = f77_zmq_recv( zmq_socket_pull, task_id, 4, 0)
|
||||
|
||||
IRP_IF ZMQ_PUSH
|
||||
IRP_ELSE
|
||||
rc = f77_zmq_send( zmq_socket_pull, 0, 4, 0)
|
||||
if (rc /= 4) then
|
||||
print *, irp_here, ' : f77_zmq_send (zmq_socket_pull,...'
|
||||
stop 'error'
|
||||
endif
|
||||
IRP_ENDIF
|
||||
|
||||
|
||||
call insert_into_ao_integrals_map(n_integrals,buffer_i,buffer_value)
|
||||
accu += n_integrals
|
||||
if (task_id /= 0) then
|
||||
integer, external :: zmq_delete_task
|
||||
if (zmq_delete_task(zmq_to_qp_run_socket,zmq_socket_pull,task_id,more) == -1) then
|
||||
stop 'Unable to delete task'
|
||||
endif
|
||||
endif
|
||||
endif
|
||||
|
||||
enddo
|
||||
|
||||
deallocate( buffer_i, buffer_value )
|
||||
|
||||
integer (map_size_kind) :: get_ao_map_size
|
||||
control = get_ao_map_size(ao_integrals_map)
|
||||
|
||||
if (control /= accu) then
|
||||
print *, ''
|
||||
print *, irp_here
|
||||
print *, 'Control : ', control
|
||||
print *, 'Accu : ', accu
|
||||
print *, 'Some integrals were lost during the parallel computation.'
|
||||
print *, 'Try to reduce the number of threads.'
|
||||
stop
|
||||
endif
|
||||
|
||||
call end_zmq_to_qp_run_socket(zmq_to_qp_run_socket)
|
||||
|
||||
end
|
||||
|
||||
|
|
@ -1,35 +0,0 @@
|
|||
BEGIN_PROVIDER [double precision, core_energy]
|
||||
implicit none
|
||||
integer :: i,j,k,l
|
||||
core_energy = 0.d0
|
||||
do i = 1, n_core_orb
|
||||
j = list_core(i)
|
||||
core_energy += 2.d0 * mo_mono_elec_integral(j,j) + mo_bielec_integral_jj(j,j)
|
||||
do k = i+1, n_core_orb
|
||||
l = list_core(k)
|
||||
core_energy += 2.d0 * (2.d0 * mo_bielec_integral_jj(j,l) - mo_bielec_integral_jj_exchange(j,l))
|
||||
enddo
|
||||
enddo
|
||||
core_energy += nuclear_repulsion
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [double precision, core_fock_operator, (mo_tot_num,mo_tot_num)]
|
||||
implicit none
|
||||
integer :: i,j,k,l,m,n
|
||||
double precision :: get_mo_bielec_integral
|
||||
BEGIN_DOC
|
||||
! this is the contribution to the Fock operator from the core electrons
|
||||
END_DOC
|
||||
core_fock_operator = 0.d0
|
||||
do i = 1, n_act_orb
|
||||
j = list_act(i)
|
||||
do k = 1, n_act_orb
|
||||
l = list_act(k)
|
||||
do m = 1, n_core_orb
|
||||
n = list_core(m)
|
||||
core_fock_operator(j,l) += 2.d0 * get_mo_bielec_integral(j,n,l,n,mo_integrals_map) - get_mo_bielec_integral(j,n,n,l,mo_integrals_map)
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
END_PROVIDER
|
||||
|
|
@ -1,57 +0,0 @@
|
|||
BEGIN_PROVIDER [ double precision, gauleg_t2, (n_pt_max_integrals,n_pt_max_integrals/2) ]
|
||||
&BEGIN_PROVIDER [ double precision, gauleg_w, (n_pt_max_integrals,n_pt_max_integrals/2) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! t_w(i,1,k) = w(i)
|
||||
! t_w(i,2,k) = t(i)
|
||||
END_DOC
|
||||
integer :: i,j,l
|
||||
l=0
|
||||
do i = 2,n_pt_max_integrals,2
|
||||
l = l+1
|
||||
call gauleg(0.d0,1.d0,gauleg_t2(1,l),gauleg_w(1,l),i)
|
||||
do j=1,i
|
||||
gauleg_t2(j,l) *= gauleg_t2(j,l)
|
||||
enddo
|
||||
enddo
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
subroutine gauleg(x1,x2,x,w,n)
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Gauss-Legendre
|
||||
END_DOC
|
||||
integer, intent(in) :: n
|
||||
double precision, intent(in) :: x1, x2
|
||||
double precision, intent (out) :: x(n),w(n)
|
||||
double precision, parameter :: eps=3.d-14
|
||||
|
||||
integer :: m,i,j
|
||||
double precision :: xm, xl, z, z1, p1, p2, p3, pp, dn
|
||||
m=(n+1)/2
|
||||
xm=0.5d0*(x2+x1)
|
||||
xl=0.5d0*(x2-x1)
|
||||
dn = dble(n)
|
||||
do i=1,m
|
||||
z=dcos(3.141592654d0*(dble(i)-.25d0)/(dble(n)+.5d0))
|
||||
z1 = z+1.d0
|
||||
do while (dabs(z-z1) > eps)
|
||||
p1=1.d0
|
||||
p2=0.d0
|
||||
do j=1,n
|
||||
p3=p2
|
||||
p2=p1
|
||||
p1=(dble(j+j-1)*z*p2-dble(j-1)*p3)/j
|
||||
enddo
|
||||
pp=dn*(z*p1-p2)/(z*z-1.d0)
|
||||
z1=z
|
||||
z=z1-p1/pp
|
||||
end do
|
||||
x(i)=xm-xl*z
|
||||
x(n+1-i)=xm+xl*z
|
||||
w(i)=(xl+xl)/((1.d0-z*z)*pp*pp)
|
||||
w(n+1-i)=w(i)
|
||||
enddo
|
||||
end
|
||||
|
||||
|
|
@ -1,22 +0,0 @@
|
|||
BEGIN_PROVIDER [double precision, big_array_coulomb_integrals, (mo_tot_num,mo_tot_num, mo_tot_num)]
|
||||
&BEGIN_PROVIDER [double precision, big_array_exchange_integrals,(mo_tot_num,mo_tot_num, mo_tot_num)]
|
||||
implicit none
|
||||
integer :: i,j,k,l
|
||||
double precision :: get_mo_bielec_integral
|
||||
double precision :: integral
|
||||
|
||||
do k = 1, mo_tot_num
|
||||
do i = 1, mo_tot_num
|
||||
do j = 1, mo_tot_num
|
||||
l = j
|
||||
integral = get_mo_bielec_integral(i,j,k,l,mo_integrals_map)
|
||||
big_array_coulomb_integrals(j,i,k) = integral
|
||||
l = j
|
||||
integral = get_mo_bielec_integral(i,j,l,k,mo_integrals_map)
|
||||
big_array_exchange_integrals(j,i,k) = integral
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
|
||||
|
||||
END_PROVIDER
|
||||
|
|
@ -1,694 +0,0 @@
|
|||
use map_module
|
||||
|
||||
!! AO Map
|
||||
!! ======
|
||||
|
||||
BEGIN_PROVIDER [ type(map_type), ao_integrals_map ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! AO integrals
|
||||
END_DOC
|
||||
integer(key_kind) :: key_max
|
||||
integer(map_size_kind) :: sze
|
||||
call bielec_integrals_index(ao_num,ao_num,ao_num,ao_num,key_max)
|
||||
sze = key_max
|
||||
call map_init(ao_integrals_map,sze)
|
||||
print*, 'AO map initialized : ', sze
|
||||
END_PROVIDER
|
||||
|
||||
subroutine bielec_integrals_index(i,j,k,l,i1)
|
||||
use map_module
|
||||
implicit none
|
||||
integer, intent(in) :: i,j,k,l
|
||||
integer(key_kind), intent(out) :: i1
|
||||
integer(key_kind) :: p,q,r,s,i2
|
||||
p = min(i,k)
|
||||
r = max(i,k)
|
||||
p = p+shiftr(r*r-r,1)
|
||||
q = min(j,l)
|
||||
s = max(j,l)
|
||||
q = q+shiftr(s*s-s,1)
|
||||
i1 = min(p,q)
|
||||
i2 = max(p,q)
|
||||
i1 = i1+shiftr(i2*i2-i2,1)
|
||||
end
|
||||
|
||||
subroutine bielec_integrals_index_reverse(i,j,k,l,i1)
|
||||
use map_module
|
||||
implicit none
|
||||
integer, intent(out) :: i(8),j(8),k(8),l(8)
|
||||
integer(key_kind), intent(in) :: i1
|
||||
integer(key_kind) :: i2,i3
|
||||
i = 0
|
||||
i2 = ceiling(0.5d0*(dsqrt(8.d0*dble(i1)+1.d0)-1.d0))
|
||||
l(1) = ceiling(0.5d0*(dsqrt(8.d0*dble(i2)+1.d0)-1.d0))
|
||||
i3 = i1 - shiftr(i2*i2-i2,1)
|
||||
k(1) = ceiling(0.5d0*(dsqrt(8.d0*dble(i3)+1.d0)-1.d0))
|
||||
j(1) = int(i2 - shiftr(l(1)*l(1)-l(1),1),4)
|
||||
i(1) = int(i3 - shiftr(k(1)*k(1)-k(1),1),4)
|
||||
|
||||
!ijkl
|
||||
i(2) = i(1) !ilkj
|
||||
j(2) = l(1)
|
||||
k(2) = k(1)
|
||||
l(2) = j(1)
|
||||
|
||||
i(3) = k(1) !kjil
|
||||
j(3) = j(1)
|
||||
k(3) = i(1)
|
||||
l(3) = l(1)
|
||||
|
||||
i(4) = k(1) !klij
|
||||
j(4) = l(1)
|
||||
k(4) = i(1)
|
||||
l(4) = j(1)
|
||||
|
||||
i(5) = j(1) !jilk
|
||||
j(5) = i(1)
|
||||
k(5) = l(1)
|
||||
l(5) = k(1)
|
||||
|
||||
i(6) = j(1) !jkli
|
||||
j(6) = k(1)
|
||||
k(6) = l(1)
|
||||
l(6) = i(1)
|
||||
|
||||
i(7) = l(1) !lijk
|
||||
j(7) = i(1)
|
||||
k(7) = j(1)
|
||||
l(7) = k(1)
|
||||
|
||||
i(8) = l(1) !lkji
|
||||
j(8) = k(1)
|
||||
k(8) = j(1)
|
||||
l(8) = i(1)
|
||||
|
||||
integer :: ii, jj
|
||||
do ii=2,8
|
||||
do jj=1,ii-1
|
||||
if ( (i(ii) == i(jj)).and. &
|
||||
(j(ii) == j(jj)).and. &
|
||||
(k(ii) == k(jj)).and. &
|
||||
(l(ii) == l(jj)) ) then
|
||||
i(ii) = 0
|
||||
exit
|
||||
endif
|
||||
enddo
|
||||
enddo
|
||||
do ii=1,8
|
||||
if (i(ii) /= 0) then
|
||||
call bielec_integrals_index(i(ii),j(ii),k(ii),l(ii),i2)
|
||||
if (i1 /= i2) then
|
||||
print *, i1, i2
|
||||
print *, i(ii), j(ii), k(ii), l(ii)
|
||||
stop 'bielec_integrals_index_reverse failed'
|
||||
endif
|
||||
endif
|
||||
enddo
|
||||
|
||||
|
||||
end
|
||||
|
||||
BEGIN_PROVIDER [ integer, ao_integrals_cache_min ]
|
||||
&BEGIN_PROVIDER [ integer, ao_integrals_cache_max ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Min and max values of the AOs for which the integrals are in the cache
|
||||
END_DOC
|
||||
ao_integrals_cache_min = max(1,ao_num - 63)
|
||||
ao_integrals_cache_max = ao_num
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ double precision, ao_integrals_cache, (0:64*64*64*64) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Cache of AO integrals for fast access
|
||||
END_DOC
|
||||
PROVIDE ao_bielec_integrals_in_map
|
||||
integer :: i,j,k,l,ii
|
||||
integer(key_kind) :: idx
|
||||
real(integral_kind) :: integral
|
||||
!$OMP PARALLEL DO PRIVATE (i,j,k,l,idx,ii,integral)
|
||||
do l=ao_integrals_cache_min,ao_integrals_cache_max
|
||||
do k=ao_integrals_cache_min,ao_integrals_cache_max
|
||||
do j=ao_integrals_cache_min,ao_integrals_cache_max
|
||||
do i=ao_integrals_cache_min,ao_integrals_cache_max
|
||||
!DIR$ FORCEINLINE
|
||||
call bielec_integrals_index(i,j,k,l,idx)
|
||||
!DIR$ FORCEINLINE
|
||||
call map_get(ao_integrals_map,idx,integral)
|
||||
ii = l-ao_integrals_cache_min
|
||||
ii = ior( shiftl(ii,6), k-ao_integrals_cache_min)
|
||||
ii = ior( shiftl(ii,6), j-ao_integrals_cache_min)
|
||||
ii = ior( shiftl(ii,6), i-ao_integrals_cache_min)
|
||||
ao_integrals_cache(ii) = integral
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
!$OMP END PARALLEL DO
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
double precision function get_ao_bielec_integral(i,j,k,l,map) result(result)
|
||||
use map_module
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Gets one AO bi-electronic integral from the AO map
|
||||
END_DOC
|
||||
integer, intent(in) :: i,j,k,l
|
||||
integer(key_kind) :: idx
|
||||
type(map_type), intent(inout) :: map
|
||||
integer :: ii
|
||||
real(integral_kind) :: tmp
|
||||
PROVIDE ao_bielec_integrals_in_map ao_integrals_cache ao_integrals_cache_min
|
||||
!DIR$ FORCEINLINE
|
||||
if (ao_overlap_abs(i,k)*ao_overlap_abs(j,l) < ao_integrals_threshold ) then
|
||||
tmp = 0.d0
|
||||
else if (ao_bielec_integral_schwartz(i,k)*ao_bielec_integral_schwartz(j,l) < ao_integrals_threshold) then
|
||||
tmp = 0.d0
|
||||
else
|
||||
ii = l-ao_integrals_cache_min
|
||||
ii = ior(ii, k-ao_integrals_cache_min)
|
||||
ii = ior(ii, j-ao_integrals_cache_min)
|
||||
ii = ior(ii, i-ao_integrals_cache_min)
|
||||
if (iand(ii, -64) /= 0) then
|
||||
!DIR$ FORCEINLINE
|
||||
call bielec_integrals_index(i,j,k,l,idx)
|
||||
!DIR$ FORCEINLINE
|
||||
call map_get(map,idx,tmp)
|
||||
else
|
||||
ii = l-ao_integrals_cache_min
|
||||
ii = ior( shiftl(ii,6), k-ao_integrals_cache_min)
|
||||
ii = ior( shiftl(ii,6), j-ao_integrals_cache_min)
|
||||
ii = ior( shiftl(ii,6), i-ao_integrals_cache_min)
|
||||
tmp = ao_integrals_cache(ii)
|
||||
endif
|
||||
endif
|
||||
result = tmp
|
||||
end
|
||||
|
||||
|
||||
subroutine get_ao_bielec_integrals(j,k,l,sze,out_val)
|
||||
use map_module
|
||||
BEGIN_DOC
|
||||
! Gets multiple AO bi-electronic integral from the AO map .
|
||||
! All i are retrieved for j,k,l fixed.
|
||||
END_DOC
|
||||
implicit none
|
||||
integer, intent(in) :: j,k,l, sze
|
||||
real(integral_kind), intent(out) :: out_val(sze)
|
||||
|
||||
integer :: i
|
||||
integer(key_kind) :: hash
|
||||
double precision :: thresh
|
||||
PROVIDE ao_bielec_integrals_in_map ao_integrals_map
|
||||
thresh = ao_integrals_threshold
|
||||
|
||||
if (ao_overlap_abs(j,l) < thresh) then
|
||||
out_val = 0.d0
|
||||
return
|
||||
endif
|
||||
|
||||
double precision :: get_ao_bielec_integral
|
||||
do i=1,sze
|
||||
out_val(i) = get_ao_bielec_integral(i,j,k,l,ao_integrals_map)
|
||||
enddo
|
||||
|
||||
end
|
||||
|
||||
subroutine get_ao_bielec_integrals_non_zero(j,k,l,sze,out_val,out_val_index,non_zero_int)
|
||||
use map_module
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Gets multiple AO bi-electronic integral from the AO map .
|
||||
! All non-zero i are retrieved for j,k,l fixed.
|
||||
END_DOC
|
||||
integer, intent(in) :: j,k,l, sze
|
||||
real(integral_kind), intent(out) :: out_val(sze)
|
||||
integer, intent(out) :: out_val_index(sze),non_zero_int
|
||||
|
||||
integer :: i
|
||||
integer(key_kind) :: hash
|
||||
double precision :: thresh,tmp
|
||||
PROVIDE ao_bielec_integrals_in_map
|
||||
thresh = ao_integrals_threshold
|
||||
|
||||
non_zero_int = 0
|
||||
if (ao_overlap_abs(j,l) < thresh) then
|
||||
out_val = 0.d0
|
||||
return
|
||||
endif
|
||||
|
||||
non_zero_int = 0
|
||||
do i=1,sze
|
||||
integer, external :: ao_l4
|
||||
double precision, external :: ao_bielec_integral
|
||||
!DIR$ FORCEINLINE
|
||||
if (ao_bielec_integral_schwartz(i,k)*ao_bielec_integral_schwartz(j,l) < thresh) then
|
||||
cycle
|
||||
endif
|
||||
call bielec_integrals_index(i,j,k,l,hash)
|
||||
call map_get(ao_integrals_map, hash,tmp)
|
||||
if (dabs(tmp) < thresh ) cycle
|
||||
non_zero_int = non_zero_int+1
|
||||
out_val_index(non_zero_int) = i
|
||||
out_val(non_zero_int) = tmp
|
||||
enddo
|
||||
|
||||
end
|
||||
|
||||
|
||||
function get_ao_map_size()
|
||||
implicit none
|
||||
integer (map_size_kind) :: get_ao_map_size
|
||||
BEGIN_DOC
|
||||
! Returns the number of elements in the AO map
|
||||
END_DOC
|
||||
get_ao_map_size = ao_integrals_map % n_elements
|
||||
end
|
||||
|
||||
subroutine clear_ao_map
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Frees the memory of the AO map
|
||||
END_DOC
|
||||
call map_deinit(ao_integrals_map)
|
||||
FREE ao_integrals_map
|
||||
end
|
||||
|
||||
|
||||
!! MO Map
|
||||
!! ======
|
||||
|
||||
BEGIN_PROVIDER [ type(map_type), mo_integrals_map ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! MO integrals
|
||||
END_DOC
|
||||
integer(key_kind) :: key_max
|
||||
integer(map_size_kind) :: sze
|
||||
call bielec_integrals_index(mo_tot_num,mo_tot_num,mo_tot_num,mo_tot_num,key_max)
|
||||
sze = key_max
|
||||
call map_init(mo_integrals_map,sze)
|
||||
print*, 'MO map initialized: ', sze
|
||||
END_PROVIDER
|
||||
|
||||
subroutine insert_into_ao_integrals_map(n_integrals,buffer_i, buffer_values)
|
||||
use map_module
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Create new entry into AO map
|
||||
END_DOC
|
||||
|
||||
integer, intent(in) :: n_integrals
|
||||
integer(key_kind), intent(inout) :: buffer_i(n_integrals)
|
||||
real(integral_kind), intent(inout) :: buffer_values(n_integrals)
|
||||
|
||||
call map_append(ao_integrals_map, buffer_i, buffer_values, n_integrals)
|
||||
end
|
||||
|
||||
subroutine insert_into_mo_integrals_map(n_integrals, &
|
||||
buffer_i, buffer_values, thr)
|
||||
use map_module
|
||||
implicit none
|
||||
|
||||
BEGIN_DOC
|
||||
! Create new entry into MO map, or accumulate in an existing entry
|
||||
END_DOC
|
||||
|
||||
integer, intent(in) :: n_integrals
|
||||
integer(key_kind), intent(inout) :: buffer_i(n_integrals)
|
||||
real(integral_kind), intent(inout) :: buffer_values(n_integrals)
|
||||
real(integral_kind), intent(in) :: thr
|
||||
call map_update(mo_integrals_map, buffer_i, buffer_values, n_integrals, thr)
|
||||
end
|
||||
|
||||
BEGIN_PROVIDER [ integer*4, mo_integrals_cache_min ]
|
||||
&BEGIN_PROVIDER [ integer*4, mo_integrals_cache_max ]
|
||||
&BEGIN_PROVIDER [ integer*8, mo_integrals_cache_min_8 ]
|
||||
&BEGIN_PROVIDER [ integer*8, mo_integrals_cache_max_8 ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Min and max values of the MOs for which the integrals are in the cache
|
||||
END_DOC
|
||||
mo_integrals_cache_min_8 = max(1_8,elec_alpha_num - 63_8)
|
||||
mo_integrals_cache_max_8 = min(int(mo_tot_num,8),mo_integrals_cache_min_8+127_8)
|
||||
mo_integrals_cache_min = max(1,elec_alpha_num - 63)
|
||||
mo_integrals_cache_max = min(mo_tot_num,mo_integrals_cache_min+127)
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ double precision, mo_integrals_cache, (0_8:128_8*128_8*128_8*128_8) ]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Cache of MO integrals for fast access
|
||||
END_DOC
|
||||
PROVIDE mo_bielec_integrals_in_map
|
||||
integer*8 :: i,j,k,l
|
||||
integer*4 :: i4,j4,k4,l4
|
||||
integer*8 :: ii
|
||||
integer(key_kind) :: idx
|
||||
real(integral_kind) :: integral
|
||||
FREE ao_integrals_cache
|
||||
!$OMP PARALLEL DO PRIVATE (i,j,k,l,i4,j4,k4,l4,idx,ii,integral)
|
||||
do l=mo_integrals_cache_min_8,mo_integrals_cache_max_8
|
||||
l4 = int(l,4)
|
||||
do k=mo_integrals_cache_min_8,mo_integrals_cache_max_8
|
||||
k4 = int(k,4)
|
||||
do j=mo_integrals_cache_min_8,mo_integrals_cache_max_8
|
||||
j4 = int(j,4)
|
||||
do i=mo_integrals_cache_min_8,mo_integrals_cache_max_8
|
||||
i4 = int(i,4)
|
||||
!DIR$ FORCEINLINE
|
||||
call bielec_integrals_index(i4,j4,k4,l4,idx)
|
||||
!DIR$ FORCEINLINE
|
||||
call map_get(mo_integrals_map,idx,integral)
|
||||
ii = l-mo_integrals_cache_min_8
|
||||
ii = ior( shiftl(ii,7), k-mo_integrals_cache_min_8)
|
||||
ii = ior( shiftl(ii,7), j-mo_integrals_cache_min_8)
|
||||
ii = ior( shiftl(ii,7), i-mo_integrals_cache_min_8)
|
||||
mo_integrals_cache(ii) = integral
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
!$OMP END PARALLEL DO
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
|
||||
double precision function get_mo_bielec_integral(i,j,k,l,map)
|
||||
use map_module
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Returns one integral <ij|kl> in the MO basis
|
||||
END_DOC
|
||||
integer, intent(in) :: i,j,k,l
|
||||
integer(key_kind) :: idx
|
||||
integer :: ii
|
||||
integer*8 :: ii_8
|
||||
type(map_type), intent(inout) :: map
|
||||
real(integral_kind) :: tmp
|
||||
PROVIDE mo_bielec_integrals_in_map mo_integrals_cache
|
||||
ii = l-mo_integrals_cache_min
|
||||
ii = ior(ii, k-mo_integrals_cache_min)
|
||||
ii = ior(ii, j-mo_integrals_cache_min)
|
||||
ii = ior(ii, i-mo_integrals_cache_min)
|
||||
if (iand(ii, -128) /= 0) then
|
||||
!DIR$ FORCEINLINE
|
||||
call bielec_integrals_index(i,j,k,l,idx)
|
||||
!DIR$ FORCEINLINE
|
||||
call map_get(map,idx,tmp)
|
||||
get_mo_bielec_integral = dble(tmp)
|
||||
else
|
||||
ii_8 = int(l,8)-mo_integrals_cache_min_8
|
||||
ii_8 = ior( shiftl(ii_8,7), int(k,8)-mo_integrals_cache_min_8)
|
||||
ii_8 = ior( shiftl(ii_8,7), int(j,8)-mo_integrals_cache_min_8)
|
||||
ii_8 = ior( shiftl(ii_8,7), int(i,8)-mo_integrals_cache_min_8)
|
||||
get_mo_bielec_integral = mo_integrals_cache(ii_8)
|
||||
endif
|
||||
end
|
||||
|
||||
|
||||
double precision function mo_bielec_integral(i,j,k,l)
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Returns one integral <ij|kl> in the MO basis
|
||||
END_DOC
|
||||
integer, intent(in) :: i,j,k,l
|
||||
double precision :: get_mo_bielec_integral
|
||||
PROVIDE mo_bielec_integrals_in_map mo_integrals_cache
|
||||
PROVIDE mo_bielec_integrals_in_map
|
||||
!DIR$ FORCEINLINE
|
||||
mo_bielec_integral = get_mo_bielec_integral(i,j,k,l,mo_integrals_map)
|
||||
return
|
||||
end
|
||||
|
||||
subroutine get_mo_bielec_integrals(j,k,l,sze,out_val,map)
|
||||
use map_module
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Returns multiple integrals <ij|kl> in the MO basis, all
|
||||
! i for j,k,l fixed.
|
||||
END_DOC
|
||||
integer, intent(in) :: j,k,l, sze
|
||||
double precision, intent(out) :: out_val(sze)
|
||||
type(map_type), intent(inout) :: map
|
||||
integer :: i
|
||||
double precision, external :: get_mo_bielec_integral
|
||||
PROVIDE mo_bielec_integrals_in_map mo_integrals_cache
|
||||
|
||||
integer :: ii, ii0
|
||||
integer*8 :: ii_8, ii0_8
|
||||
real(integral_kind) :: tmp
|
||||
integer(key_kind) :: i1, idx
|
||||
integer(key_kind) :: p,q,r,s,i2
|
||||
PROVIDE mo_bielec_integrals_in_map mo_integrals_cache
|
||||
|
||||
ii0 = l-mo_integrals_cache_min
|
||||
ii0 = ior(ii0, k-mo_integrals_cache_min)
|
||||
ii0 = ior(ii0, j-mo_integrals_cache_min)
|
||||
|
||||
ii0_8 = int(l,8)-mo_integrals_cache_min_8
|
||||
ii0_8 = ior( shiftl(ii0_8,7), int(k,8)-mo_integrals_cache_min_8)
|
||||
ii0_8 = ior( shiftl(ii0_8,7), int(j,8)-mo_integrals_cache_min_8)
|
||||
|
||||
q = min(j,l)
|
||||
s = max(j,l)
|
||||
q = q+shiftr(s*s-s,1)
|
||||
|
||||
do i=1,sze
|
||||
ii = ior(ii0, i-mo_integrals_cache_min)
|
||||
if (iand(ii, -128) == 0) then
|
||||
ii_8 = ior( shiftl(ii0_8,7), int(i,8)-mo_integrals_cache_min_8)
|
||||
out_val(i) = mo_integrals_cache(ii_8)
|
||||
else
|
||||
p = min(i,k)
|
||||
r = max(i,k)
|
||||
p = p+shiftr(r*r-r,1)
|
||||
i1 = min(p,q)
|
||||
i2 = max(p,q)
|
||||
idx = i1+shiftr(i2*i2-i2,1)
|
||||
!DIR$ FORCEINLINE
|
||||
call map_get(map,idx,tmp)
|
||||
out_val(i) = dble(tmp)
|
||||
endif
|
||||
enddo
|
||||
end
|
||||
|
||||
subroutine get_mo_bielec_integrals_ij(k,l,sze,out_array,map)
|
||||
use map_module
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Returns multiple integrals <ij|kl> in the MO basis, all
|
||||
! i(1)j(2) 1/r12 k(1)l(2)
|
||||
! i, j for k,l fixed.
|
||||
END_DOC
|
||||
integer, intent(in) :: k,l, sze
|
||||
double precision, intent(out) :: out_array(sze,sze)
|
||||
type(map_type), intent(inout) :: map
|
||||
integer :: i,j,kk,ll,m
|
||||
integer(key_kind),allocatable :: hash(:)
|
||||
integer ,allocatable :: pairs(:,:), iorder(:)
|
||||
real(integral_kind), allocatable :: tmp_val(:)
|
||||
|
||||
PROVIDE mo_bielec_integrals_in_map
|
||||
allocate (hash(sze*sze), pairs(2,sze*sze),iorder(sze*sze), &
|
||||
tmp_val(sze*sze))
|
||||
|
||||
kk=0
|
||||
out_array = 0.d0
|
||||
do j=1,sze
|
||||
do i=1,sze
|
||||
kk += 1
|
||||
!DIR$ FORCEINLINE
|
||||
call bielec_integrals_index(i,j,k,l,hash(kk))
|
||||
pairs(1,kk) = i
|
||||
pairs(2,kk) = j
|
||||
iorder(kk) = kk
|
||||
enddo
|
||||
enddo
|
||||
|
||||
logical :: integral_is_in_map
|
||||
if (key_kind == 8) then
|
||||
call i8radix_sort(hash,iorder,kk,-1)
|
||||
else if (key_kind == 4) then
|
||||
call iradix_sort(hash,iorder,kk,-1)
|
||||
else if (key_kind == 2) then
|
||||
call i2radix_sort(hash,iorder,kk,-1)
|
||||
endif
|
||||
|
||||
call map_get_many(mo_integrals_map, hash, tmp_val, kk)
|
||||
|
||||
do ll=1,kk
|
||||
m = iorder(ll)
|
||||
i=pairs(1,m)
|
||||
j=pairs(2,m)
|
||||
out_array(i,j) = tmp_val(ll)
|
||||
enddo
|
||||
|
||||
deallocate(pairs,hash,iorder,tmp_val)
|
||||
end
|
||||
|
||||
subroutine get_mo_bielec_integrals_coulomb_ii(k,l,sze,out_val,map)
|
||||
use map_module
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Returns multiple integrals <ki|li>
|
||||
! k(1)i(2) 1/r12 l(1)i(2) :: out_val(i1)
|
||||
! for k,l fixed.
|
||||
END_DOC
|
||||
integer, intent(in) :: k,l, sze
|
||||
double precision, intent(out) :: out_val(sze)
|
||||
type(map_type), intent(inout) :: map
|
||||
integer :: i
|
||||
integer(key_kind) :: hash(sze)
|
||||
real(integral_kind) :: tmp_val(sze)
|
||||
PROVIDE mo_bielec_integrals_in_map
|
||||
|
||||
integer :: kk
|
||||
do i=1,sze
|
||||
!DIR$ FORCEINLINE
|
||||
call bielec_integrals_index(k,i,l,i,hash(i))
|
||||
enddo
|
||||
|
||||
if (integral_kind == 8) then
|
||||
call map_get_many(map, hash, out_val, sze)
|
||||
else
|
||||
call map_get_many(map, hash, tmp_val, sze)
|
||||
! Conversion to double precision
|
||||
do i=1,sze
|
||||
out_val(i) = dble(tmp_val(i))
|
||||
enddo
|
||||
endif
|
||||
end
|
||||
|
||||
subroutine get_mo_bielec_integrals_exch_ii(k,l,sze,out_val,map)
|
||||
use map_module
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Returns multiple integrals <ki|il>
|
||||
! k(1)i(2) 1/r12 i(1)l(2) :: out_val(i1)
|
||||
! for k,l fixed.
|
||||
END_DOC
|
||||
integer, intent(in) :: k,l, sze
|
||||
double precision, intent(out) :: out_val(sze)
|
||||
type(map_type), intent(inout) :: map
|
||||
integer :: i
|
||||
integer(key_kind) :: hash(sze)
|
||||
real(integral_kind) :: tmp_val(sze)
|
||||
PROVIDE mo_bielec_integrals_in_map
|
||||
|
||||
integer :: kk
|
||||
do i=1,sze
|
||||
!DIR$ FORCEINLINE
|
||||
call bielec_integrals_index(k,i,i,l,hash(i))
|
||||
enddo
|
||||
|
||||
if (integral_kind == 8) then
|
||||
call map_get_many(map, hash, out_val, sze)
|
||||
else
|
||||
call map_get_many(map, hash, tmp_val, sze)
|
||||
! Conversion to double precision
|
||||
do i=1,sze
|
||||
out_val(i) = dble(tmp_val(i))
|
||||
enddo
|
||||
endif
|
||||
end
|
||||
|
||||
|
||||
integer*8 function get_mo_map_size()
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Return the number of elements in the MO map
|
||||
END_DOC
|
||||
get_mo_map_size = mo_integrals_map % n_elements
|
||||
end
|
||||
|
||||
BEGIN_TEMPLATE
|
||||
|
||||
subroutine dump_$ao_integrals(filename)
|
||||
use map_module
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Save to disk the $ao integrals
|
||||
END_DOC
|
||||
character*(*), intent(in) :: filename
|
||||
integer(cache_key_kind), pointer :: key(:)
|
||||
real(integral_kind), pointer :: val(:)
|
||||
integer*8 :: i,j, n
|
||||
if (.not.mpi_master) then
|
||||
return
|
||||
endif
|
||||
call ezfio_set_work_empty(.False.)
|
||||
open(unit=66,file=filename,FORM='unformatted')
|
||||
write(66) integral_kind, key_kind
|
||||
write(66) $ao_integrals_map%sorted, $ao_integrals_map%map_size, &
|
||||
$ao_integrals_map%n_elements
|
||||
do i=0_8,$ao_integrals_map%map_size
|
||||
write(66) $ao_integrals_map%map(i)%sorted, $ao_integrals_map%map(i)%map_size,&
|
||||
$ao_integrals_map%map(i)%n_elements
|
||||
enddo
|
||||
do i=0_8,$ao_integrals_map%map_size
|
||||
key => $ao_integrals_map%map(i)%key
|
||||
val => $ao_integrals_map%map(i)%value
|
||||
n = $ao_integrals_map%map(i)%n_elements
|
||||
write(66) (key(j), j=1,n), (val(j), j=1,n)
|
||||
enddo
|
||||
close(66)
|
||||
|
||||
end
|
||||
|
||||
|
||||
integer function load_$ao_integrals(filename)
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! Read from disk the $ao integrals
|
||||
END_DOC
|
||||
character*(*), intent(in) :: filename
|
||||
integer*8 :: i
|
||||
integer(cache_key_kind), pointer :: key(:)
|
||||
real(integral_kind), pointer :: val(:)
|
||||
integer :: iknd, kknd
|
||||
integer*8 :: n, j
|
||||
load_$ao_integrals = 1
|
||||
open(unit=66,file=filename,FORM='unformatted',STATUS='UNKNOWN')
|
||||
read(66,err=98,end=98) iknd, kknd
|
||||
if (iknd /= integral_kind) then
|
||||
print *, 'Wrong integrals kind in file :', iknd
|
||||
stop 1
|
||||
endif
|
||||
if (kknd /= key_kind) then
|
||||
print *, 'Wrong key kind in file :', kknd
|
||||
stop 1
|
||||
endif
|
||||
read(66,err=98,end=98) $ao_integrals_map%sorted, $ao_integrals_map%map_size,&
|
||||
$ao_integrals_map%n_elements
|
||||
do i=0_8, $ao_integrals_map%map_size
|
||||
read(66,err=99,end=99) $ao_integrals_map%map(i)%sorted, &
|
||||
$ao_integrals_map%map(i)%map_size, $ao_integrals_map%map(i)%n_elements
|
||||
call cache_map_reallocate($ao_integrals_map%map(i),$ao_integrals_map%map(i)%map_size)
|
||||
enddo
|
||||
do i=0_8, $ao_integrals_map%map_size
|
||||
key => $ao_integrals_map%map(i)%key
|
||||
val => $ao_integrals_map%map(i)%value
|
||||
n = $ao_integrals_map%map(i)%n_elements
|
||||
read(66,err=99,end=99) (key(j), j=1,n), (val(j), j=1,n)
|
||||
enddo
|
||||
call map_sort($ao_integrals_map)
|
||||
load_$ao_integrals = 0
|
||||
return
|
||||
99 continue
|
||||
call map_deinit($ao_integrals_map)
|
||||
98 continue
|
||||
stop 'Problem reading $ao_integrals_map file in work/'
|
||||
|
||||
end
|
||||
|
||||
SUBST [ ao_integrals_map, ao_integrals, ao_num ]
|
||||
ao_integrals_map ; ao_integrals ; ao_num ;;
|
||||
mo_integrals_map ; mo_integrals ; mo_tot_num ;;
|
||||
END_TEMPLATE
|
||||
File diff suppressed because it is too large
Load Diff
|
|
@ -1,47 +0,0 @@
|
|||
BEGIN_PROVIDER [ logical, read_ao_integrals ]
|
||||
&BEGIN_PROVIDER [ logical, read_mo_integrals ]
|
||||
&BEGIN_PROVIDER [ logical, write_ao_integrals ]
|
||||
&BEGIN_PROVIDER [ logical, write_mo_integrals ]
|
||||
|
||||
BEGIN_DOC
|
||||
! One level of abstraction for disk_access_ao_integrals and disk_access_mo_integrals
|
||||
END_DOC
|
||||
implicit none
|
||||
|
||||
if (disk_access_ao_integrals.EQ.'Read') then
|
||||
read_ao_integrals = .True.
|
||||
write_ao_integrals = .False.
|
||||
|
||||
else if (disk_access_ao_integrals.EQ.'Write') then
|
||||
read_ao_integrals = .False.
|
||||
write_ao_integrals = .True.
|
||||
|
||||
else if (disk_access_ao_integrals.EQ.'None') then
|
||||
read_ao_integrals = .False.
|
||||
write_ao_integrals = .False.
|
||||
|
||||
else
|
||||
print *, 'bielec_integrals/disk_access_ao_integrals has a wrong type'
|
||||
stop 1
|
||||
|
||||
endif
|
||||
|
||||
if (disk_access_mo_integrals.EQ.'Read') then
|
||||
read_mo_integrals = .True.
|
||||
write_mo_integrals = .False.
|
||||
|
||||
else if (disk_access_mo_integrals.EQ.'Write') then
|
||||
read_mo_integrals = .False.
|
||||
write_mo_integrals = .True.
|
||||
|
||||
else if (disk_access_mo_integrals.EQ.'None') then
|
||||
read_mo_integrals = .False.
|
||||
write_mo_integrals = .False.
|
||||
|
||||
else
|
||||
print *, 'bielec_integrals/disk_access_mo_integrals has a wrong type'
|
||||
stop 1
|
||||
|
||||
endif
|
||||
|
||||
END_PROVIDER
|
||||
|
|
@ -1,2 +0,0 @@
|
|||
MO_Basis
|
||||
MO_one_e_integrals
|
||||
|
|
@ -1,3 +0,0 @@
|
|||
AO_Basis
|
||||
AO_one_e_integrals
|
||||
Electrons
|
||||
|
|
@ -1,3 +0,0 @@
|
|||
AO_one_e_integrals
|
||||
MO_Basis
|
||||
Pseudo
|
||||
|
|
@ -1,3 +0,0 @@
|
|||
Determinants
|
||||
Davidson
|
||||
Psiref_CAS
|
||||
|
|
@ -1,2 +0,0 @@
|
|||
Ezfio_files
|
||||
Utils
|
||||
|
|
@ -1,4 +0,0 @@
|
|||
Determinants
|
||||
Hartree_Fock
|
||||
Davidson
|
||||
MRPT_Utils
|
||||
|
|
@ -1 +0,0 @@
|
|||
Nuclei
|
||||
|
|
@ -1 +0,0 @@
|
|||
Psiref_Utils
|
||||
|
|
@ -1,2 +0,0 @@
|
|||
Bitmask
|
||||
Determinants
|
||||
1
src/Selectors_CASSD/.gitignore
vendored
1
src/Selectors_CASSD/.gitignore
vendored
|
|
@ -1 +0,0 @@
|
|||
../../data/module_gitignore
|
||||
|
|
@ -1 +0,0 @@
|
|||
Selectors_Utils
|
||||
1
src/Selectors_Utils/.gitignore
vendored
1
src/Selectors_Utils/.gitignore
vendored
|
|
@ -1 +0,0 @@
|
|||
../../data/module_gitignore
|
||||
|
|
@ -1 +0,0 @@
|
|||
Determinants
|
||||
1
src/Selectors_full/.gitignore
vendored
1
src/Selectors_full/.gitignore
vendored
|
|
@ -1 +0,0 @@
|
|||
../../data/module_gitignore
|
||||
|
|
@ -1,3 +0,0 @@
|
|||
Determinants
|
||||
Hartree_Fock
|
||||
Selectors_Utils
|
||||
1
src/SingleRefMethod/.gitignore
vendored
1
src/SingleRefMethod/.gitignore
vendored
|
|
@ -1 +0,0 @@
|
|||
../../data/module_gitignore
|
||||
|
|
@ -1 +0,0 @@
|
|||
Bitmask
|
||||
1
src/Slave/.gitignore
vendored
1
src/Slave/.gitignore
vendored
|
|
@ -1 +0,0 @@
|
|||
../../data/module_gitignore
|
||||
|
|
@ -1,2 +0,0 @@
|
|||
FCI
|
||||
MPI
|
||||
1
src/Tools/.gitignore
vendored
1
src/Tools/.gitignore
vendored
|
|
@ -1 +0,0 @@
|
|||
../../data/module_gitignore
|
||||
|
|
@ -1 +0,0 @@
|
|||
FCI
|
||||
1
src/Utils/.gitignore
vendored
1
src/Utils/.gitignore
vendored
|
|
@ -1 +0,0 @@
|
|||
../../data/module_gitignore
|
||||
1
src/ZMQ/.gitignore
vendored
1
src/ZMQ/.gitignore
vendored
|
|
@ -1 +0,0 @@
|
|||
../../data/module_gitignore
|
||||
|
|
@ -1 +0,0 @@
|
|||
Utils
|
||||
1
src/ao_basis/NEED
Normal file
1
src/ao_basis/NEED
Normal file
|
|
@ -0,0 +1 @@
|
|||
nuclei
|
||||
2
src/ao_one_e_integrals/NEED
Normal file
2
src/ao_one_e_integrals/NEED
Normal file
|
|
@ -0,0 +1,2 @@
|
|||
ao_basis
|
||||
pseudo
|
||||
Some files were not shown because too many files have changed in this diff Show More
Loading…
Reference in New Issue
Block a user