9
1
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-11-14 17:33:38 +01:00
qp2/src/tools/fcidump.irp.f

79 lines
2.1 KiB
Fortran
Raw Normal View History

2019-01-25 11:39:31 +01:00
program fcidump
implicit none
BEGIN_DOC
2019-01-29 23:10:00 +01:00
! Produce a regular `FCIDUMP` file from the |MOs| stored in the |EZFIO|
! directory.
2019-01-25 11:39:31 +01:00
!
2019-01-29 23:10:00 +01:00
! To specify an active space, the class of the |MOs| have to set in the
! |EZFIO| directory (see :ref:`qp_set_mo_class`).
2019-01-25 11:39:31 +01:00
!
2019-01-29 23:10:00 +01:00
! The :ref:`fcidump` program supports 3 types of |MO| classes :
2019-01-25 11:39:31 +01:00
!
2019-01-29 23:10:00 +01:00
! * the *core* orbitals which are always doubly occupied in the
! calculation
2019-01-25 11:39:31 +01:00
!
2019-01-29 23:10:00 +01:00
! * the *deleted* orbitals that are never occupied in the calculation
!
! * the *active* orbitals that are occupied with a varying number of
! electrons
2019-01-25 11:39:31 +01:00
!
END_DOC
character*(128) :: output
integer :: i_unit_output,getUnitAndOpen
output=trim(ezfio_filename)//'.FCIDUMP'
i_unit_output = getUnitAndOpen(output,'w')
integer :: i,j,k,l
integer :: i1,j1,k1,l1
integer :: i2,j2,k2,l2
integer*8 :: m
character*(2), allocatable :: A(:)
write(i_unit_output,*) '&FCI NORB=', n_act_orb, ', NELEC=', elec_num-n_core_orb*2, &
', MS2=', (elec_alpha_num-elec_beta_num), ','
allocate (A(n_act_orb))
A = '1,'
write(i_unit_output,*) 'ORBSYM=', (A(i), i=1,n_act_orb)
write(i_unit_output,*) 'ISYM=0,'
write(i_unit_output,*) '/'
deallocate(A)
integer(key_kind), allocatable :: keys(:)
double precision, allocatable :: values(:)
integer(cache_map_size_kind) :: n_elements, n_elements_max
PROVIDE mo_two_e_integrals_in_map
double precision :: get_two_e_integral, integral
do l=1,n_act_orb
l1 = list_act(l)
do k=1,n_act_orb
k1 = list_act(k)
do j=l,n_act_orb
j1 = list_act(j)
do i=k,n_act_orb
i1 = list_act(i)
if (i1>=j1) then
integral = get_two_e_integral(i1,j1,k1,l1,mo_integrals_map)
if (dabs(integral) > mo_integrals_threshold) then
write(i_unit_output,*) integral, i,k,j,l
endif
end if
enddo
enddo
enddo
enddo
do j=1,n_act_orb
j1 = list_act(j)
do i=j,n_act_orb
i1 = list_act(i)
integral = mo_one_e_integrals(i1,j1) + core_fock_operator(i1,j1)
if (dabs(integral) > mo_integrals_threshold) then
write(i_unit_output,*) integral, i,j,0,0
endif
enddo
enddo
write(i_unit_output,*) core_energy, 0, 0, 0, 0
end