diff --git a/src/bitmask/core_inact_act_virt.irp.f b/src/bitmask/core_inact_act_virt.irp.f index d2efef89..ae00d774 100644 --- a/src/bitmask/core_inact_act_virt.irp.f +++ b/src/bitmask/core_inact_act_virt.irp.f @@ -448,7 +448,7 @@ BEGIN_PROVIDER [ integer, n_core_orb_kpts, (kpt_num)] do k=1,kpt_num n_core_orb_kpts(k) = 0 - kshift = (1-k)*mo_num_per_kpt + kshift = (k-1)*mo_num_per_kpt do i = 1, mo_num_per_kpt if(mo_class(i+kshift) == 'Core')then n_core_orb_kpts(k) += 1 @@ -469,7 +469,7 @@ BEGIN_PROVIDER [ integer, n_inact_orb_kpts, (kpt_num)] do k=1,kpt_num n_inact_orb_kpts(k) = 0 - kshift = (1-k)*mo_num_per_kpt + kshift = (k-1)*mo_num_per_kpt do i = 1, mo_num_per_kpt if(mo_class(i+kshift) == 'Inactive')then n_inact_orb_kpts(k) += 1 @@ -490,7 +490,7 @@ BEGIN_PROVIDER [ integer, n_act_orb_kpts, (kpt_num)] do k=1,kpt_num n_act_orb_kpts(k) = 0 - kshift = (1-k)*mo_num_per_kpt + kshift = (k-1)*mo_num_per_kpt do i = 1, mo_num_per_kpt if(mo_class(i+kshift) == 'Active')then n_act_orb_kpts(k) += 1 @@ -511,7 +511,7 @@ BEGIN_PROVIDER [ integer, n_virt_orb_kpts, (kpt_num)] do k=1,kpt_num n_virt_orb_kpts(k) = 0 - kshift = (1-k)*mo_num_per_kpt + kshift = (k-1)*mo_num_per_kpt do i = 1, mo_num_per_kpt if(mo_class(i+kshift) == 'Virtual')then n_virt_orb_kpts(k) += 1 @@ -532,7 +532,7 @@ BEGIN_PROVIDER [ integer, n_del_orb_kpts, (kpt_num)] do k=1,kpt_num n_del_orb_kpts(k) = 0 - kshift = (1-k)*mo_num_per_kpt + kshift = (k-1)*mo_num_per_kpt do i = 1, mo_num_per_kpt if(mo_class(i+kshift) == 'Deleted')then n_del_orb_kpts(k) += 1 diff --git a/src/mo_two_e_ints/df_mo_ints.irp.f b/src/mo_two_e_ints/df_mo_ints.irp.f index 3a61911e..eba3b3da 100644 --- a/src/mo_two_e_ints/df_mo_ints.irp.f +++ b/src/mo_two_e_ints/df_mo_ints.irp.f @@ -48,7 +48,8 @@ subroutine mo_map_fill_from_df_dot logical :: use_map1 integer(key_kind) :: idx_tmp double precision :: sign - complex*16, external :: zdotc + !complex*16, external :: zdotc + complex*16, external :: zdotu mo_num_kpt_2 = mo_num_per_kpt * mo_num_per_kpt @@ -145,7 +146,8 @@ subroutine mo_map_fill_from_df_dot if ((j==l) .and. (i>k)) exit call idx2_tri_int(i,k,ik2) if (ik2 > jl2) exit - integral = zdotc(df_num,ints_jl(1,ij,il),1,ints_ik(1,ii,ik),1) + !integral = zdotc(df_num,ints_jl(1,ij,il),1,ints_ik(1,ii,ik),1) + integral = zdotu(df_num,ints_jl(1,ij,il),1,ints_ik(1,ii,ik),1) ! print*,i,k,j,l,real(integral),imag(integral) if (cdabs(integral) < mo_integrals_threshold) then cycle diff --git a/src/tools/fcidump.irp.f b/src/tools/fcidump.irp.f index bf4d07fb..de878dc6 100644 --- a/src/tools/fcidump.irp.f +++ b/src/tools/fcidump.irp.f @@ -18,6 +18,97 @@ program fcidump ! electrons ! END_DOC + if (is_complex) then + call fcidump_complex + else + call fcidump_real + endif +end + +subroutine fcidump_complex + implicit none + 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,ik2,jl2 + integer :: ki,kj,kk,kl + integer :: ii,ij,ik,il + 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 + + complex*16 :: get_two_e_integral_complex, integral + + do kl=1,kpt_num + do kj=1,kl + do kk=1,kl + ki=kconserv(kl,kk,kj) + if (ki>kl) cycle + do l1=1,n_act_orb_kpts(kl) + il=list_act_kpts(l1,kl) + l = (kl-1)*mo_num_per_kpt + il + do j1=1,n_act_orb_kpts(kj) + ij=list_act_kpts(j1,kj) + j = (kj-1)*mo_num_per_kpt + ij + if (j>l) exit + call idx2_tri_int(j,l,jl2) + do k1=1,n_act_orb_kpts(kk) + ik=list_act_kpts(k1,kk) + k = (kk-1)*mo_num_per_kpt + ik + if (k>l) exit + do i1=1,n_act_orb_kpts(ki) + ii=list_act_kpts(i1,ki) + i = (ki-1)*mo_num_per_kpt + ii + if ((j==l) .and. (i>k)) exit + call idx2_tri_int(i,k,ik2) + if (ik2 > jl2) exit + integral = get_two_e_integral_complex(i,j,k,l,mo_integrals_map,mo_integrals_map_2) + if (cdabs(integral) > mo_integrals_threshold) then + write(i_unit_output,'(2(E25.15,X),4(I6,X))') dble(integral), dimag(integral),i,k,j,l + endif + enddo + enddo + enddo + enddo + enddo + enddo + enddo + + do kj=1,kpt_num + do j1=1,n_act_orb_kpts(kj) + ij = list_act_kpts(j1,kj) + j = (kj-1)*mo_num_per_kpt + ij + do i1=j1,n_act_orb_kpts(kj) + ii = list_act_kpts(i1,kj) + i = (kj-1)*mo_num_per_kpt + ii + integral = mo_one_e_integrals_kpts(ii,ij,kj) + core_fock_operator_complex(i,j) + if (cdabs(integral) > mo_integrals_threshold) then + write(i_unit_output,'(2(E25.15,X),4(I6,X))') dble(integral),dimag(integral), i,j,0,0 + endif + enddo + enddo + enddo + write(i_unit_output,*) core_energy, 0, 0, 0, 0 +end +subroutine fcidump_real + implicit none character*(128) :: output integer :: i_unit_output,getUnitAndOpen output=trim(ezfio_filename)//'.FCIDUMP' diff --git a/src/utils/constants.include.F b/src/utils/constants.include.F index 7399b4a6..bad68054 100644 --- a/src/utils/constants.include.F +++ b/src/utils/constants.include.F @@ -7,6 +7,7 @@ double precision, parameter :: sqpi = dsqrt(dacos(-1.d0)) double precision, parameter :: pi_5_2 = 34.9868366552d0 double precision, parameter :: dfour_pi = 4.d0*dacos(-1.d0) double precision, parameter :: dtwo_pi = 2.d0*dacos(-1.d0) +double precision, parameter :: inv_pi = 1.d0/dacos(-1.d0) double precision, parameter :: inv_sq_pi = 1.d0/dsqrt(dacos(-1.d0)) double precision, parameter :: inv_sq_pi_2 = 0.5d0/dsqrt(dacos(-1.d0)) double precision, parameter :: thresh = 1.d-15