mirror of
https://github.com/LCPQ/quantum_package
synced 2024-09-27 03:51:01 +02:00
added density_for_dft, dft_keywords and data_energy_and_density
This commit is contained in:
Emmanuel Giner
parent
a750ed71cc
commit
bc39df1a01
368
src/ao_one_e_integrals/pot_ao_erf_ints.irp.f
Normal file
368
src/ao_one_e_integrals/pot_ao_erf_ints.irp.f
Normal file
@ -0,0 +1,368 @@
|
||||
subroutine give_all_erf_kl_ao(integrals_ao,mu_in,C_center)
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! subroutine that returs all integrals over r of type erf(mu_in * |r-C_center|)/|r-C_center|
|
||||
END_DOC
|
||||
double precision, intent(in) :: mu_in,C_center(3)
|
||||
double precision, intent(out) :: integrals_ao(ao_num,ao_num)
|
||||
double precision :: NAI_pol_mult_erf_ao
|
||||
integer :: i,j,l,k,m
|
||||
do k = 1, ao_num
|
||||
do m = 1, ao_num
|
||||
integrals_ao(m,k) = NAI_pol_mult_erf_ao(m,k,mu_in,C_center)
|
||||
enddo
|
||||
enddo
|
||||
end
|
||||
|
||||
|
||||
double precision function NAI_pol_mult_erf_ao(i_ao,j_ao,mu_in,C_center)
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! computes the following integral :
|
||||
! int[-infty;+infty] dr AO_i_ao (r) AO_j_ao(r) erf(mu_in * |r-C_center|)/|r-C_center|
|
||||
END_DOC
|
||||
integer, intent(in) :: i_ao,j_ao
|
||||
double precision, intent(in) :: mu_in, C_center(3)
|
||||
integer :: i,j,num_A,num_B, power_A(3), power_B(3), n_pt_in
|
||||
double precision :: A_center(3), B_center(3),integral, alpha,beta, NAI_pol_mult_erf
|
||||
num_A = ao_nucl(i_ao)
|
||||
power_A(1:3)= ao_power(i_ao,1:3)
|
||||
A_center(1:3) = nucl_coord(num_A,1:3)
|
||||
num_B = ao_nucl(j_ao)
|
||||
power_B(1:3)= ao_power(j_ao,1:3)
|
||||
B_center(1:3) = nucl_coord(num_B,1:3)
|
||||
n_pt_in = n_pt_max_integrals
|
||||
NAI_pol_mult_erf_ao = 0.d0
|
||||
do i = 1, ao_prim_num(i_ao)
|
||||
alpha = ao_expo_ordered_transp(i,i_ao)
|
||||
do j = 1, ao_prim_num(j_ao)
|
||||
beta = ao_expo_ordered_transp(j,j_ao)
|
||||
integral = NAI_pol_mult_erf(A_center,B_center,power_A,power_B,alpha,beta,C_center,n_pt_in,mu_in)
|
||||
NAI_pol_mult_erf_ao += integral * ao_coef_normalized_ordered_transp(j,j_ao)*ao_coef_normalized_ordered_transp(i,i_ao)
|
||||
enddo
|
||||
enddo
|
||||
|
||||
|
||||
|
||||
end
|
||||
|
||||
|
||||
|
||||
double precision function NAI_pol_mult_erf(A_center,B_center,power_A,power_B,alpha,beta,C_center,n_pt_in,mu_in)
|
||||
! function that computes the folowing integral :
|
||||
! int{dr} of (x-A_x)^ax (x-B_X)^bx exp(-alpha (x-A_x)^2 - beta (x-B_x)^2 ) erf(mu_in*(r-R_c))/(r-R_c)
|
||||
|
||||
implicit none
|
||||
integer, intent(in) :: n_pt_in
|
||||
double precision,intent(in) :: C_center(3),A_center(3),B_center(3),alpha,beta,mu_in
|
||||
integer, intent(in) :: power_A(3),power_B(3)
|
||||
integer :: i,j,k,l,n_pt
|
||||
double precision :: P_center(3)
|
||||
|
||||
double precision :: d(0:n_pt_in),pouet,coeff,dist,const,pouet_2,factor
|
||||
double precision :: I_n_special_exact,integrate_bourrin,I_n_bibi
|
||||
double precision :: V_e_n,const_factor,dist_integral,tmp
|
||||
double precision :: accu,rint,p_inv,p,rho,p_inv_2
|
||||
integer :: n_pt_out,lmax
|
||||
include 'Utils/constants.include.F'
|
||||
p = alpha + beta
|
||||
p_inv = 1.d0/p
|
||||
p_inv_2 = 0.5d0 * p_inv
|
||||
rho = alpha * beta * p_inv
|
||||
|
||||
dist = 0.d0
|
||||
dist_integral = 0.d0
|
||||
do i = 1, 3
|
||||
P_center(i) = (alpha * A_center(i) + beta * B_center(i)) * p_inv
|
||||
dist += (A_center(i) - B_center(i))*(A_center(i) - B_center(i))
|
||||
dist_integral += (P_center(i) - C_center(i))*(P_center(i) - C_center(i))
|
||||
enddo
|
||||
const_factor = dist*rho
|
||||
if(const_factor > 80.d0)then
|
||||
NAI_pol_mult_erf = 0.d0
|
||||
return
|
||||
endif
|
||||
double precision :: p_new
|
||||
p_new = mu_in/dsqrt(p+ mu_in * mu_in)
|
||||
factor = dexp(-const_factor)
|
||||
coeff = dtwo_pi * factor * p_inv * p_new
|
||||
lmax = 20
|
||||
|
||||
! print*, "b"
|
||||
do i = 0, n_pt_in
|
||||
d(i) = 0.d0
|
||||
enddo
|
||||
n_pt = 2 * ( (power_A(1) + power_B(1)) +(power_A(2) + power_B(2)) +(power_A(3) + power_B(3)) )
|
||||
const = p * dist_integral * p_new * p_new
|
||||
if (n_pt == 0) then
|
||||
pouet = rint(0,const)
|
||||
NAI_pol_mult_erf = coeff * pouet
|
||||
return
|
||||
endif
|
||||
|
||||
! call give_polynom_mult_center_mono_elec_erf(A_center,B_center,alpha,beta,power_A,power_B,C_center,n_pt_in,d,n_pt_out,mu_in)
|
||||
p_new = p_new * p_new
|
||||
call give_polynom_mult_center_mono_elec_erf_opt(A_center,B_center,alpha,beta,power_A,power_B,C_center,n_pt_in,d,n_pt_out,mu_in,p,p_inv,p_inv_2,p_new,P_center)
|
||||
|
||||
|
||||
if(n_pt_out<0)then
|
||||
NAI_pol_mult_erf = 0.d0
|
||||
return
|
||||
endif
|
||||
accu = 0.d0
|
||||
|
||||
! sum of integrals of type : int {t,[0,1]} exp-(rho.(P-Q)^2 * t^2) * t^i
|
||||
do i =0 ,n_pt_out,2
|
||||
accu += d(i) * rint(i/2,const)
|
||||
enddo
|
||||
NAI_pol_mult_erf = accu * coeff
|
||||
|
||||
end
|
||||
|
||||
|
||||
subroutine give_polynom_mult_center_mono_elec_erf_opt(A_center,B_center,alpha,beta,power_A,power_B,C_center,n_pt_in,d,n_pt_out,mu_in,p,p_inv,p_inv_2,p_new,P_center)
|
||||
!!!! subroutine that returns the explicit polynom in term of the "t" variable of the following polynomw ::
|
||||
!!!! I_x1(a_x, d_x,p,q) * I_x1(a_y, d_y,p,q) * I_x1(a_z, d_z,p,q)
|
||||
!!!! it is for the nuclear electron atraction
|
||||
implicit none
|
||||
integer, intent(in) :: n_pt_in
|
||||
integer,intent(out) :: n_pt_out
|
||||
double precision, intent(in) :: A_center(3), B_center(3),C_center(3),p,p_inv,p_inv_2,p_new,P_center(3)
|
||||
double precision, intent(in) :: alpha,beta,mu_in
|
||||
integer, intent(in) :: power_A(3), power_B(3)
|
||||
integer :: a_x,b_x,a_y,b_y,a_z,b_z
|
||||
double precision :: d(0:n_pt_in)
|
||||
double precision :: d1(0:n_pt_in)
|
||||
double precision :: d2(0:n_pt_in)
|
||||
double precision :: d3(0:n_pt_in)
|
||||
double precision :: accu
|
||||
accu = 0.d0
|
||||
!COMPTEUR irp_rdtsc1 = irp_rdtsc()
|
||||
ASSERT (n_pt_in > 1)
|
||||
|
||||
double precision :: R1x(0:2), B01(0:2), R1xp(0:2),R2x(0:2)
|
||||
R1x(0) = (P_center(1) - A_center(1))
|
||||
R1x(1) = 0.d0
|
||||
R1x(2) = -(P_center(1) - C_center(1))* p_new
|
||||
! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2
|
||||
R1xp(0) = (P_center(1) - B_center(1))
|
||||
R1xp(1) = 0.d0
|
||||
R1xp(2) =-(P_center(1) - C_center(1))* p_new
|
||||
!R1xp = (P_x - B_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2
|
||||
R2x(0) = p_inv_2
|
||||
R2x(1) = 0.d0
|
||||
R2x(2) = -p_inv_2* p_new
|
||||
!R2x = 0.5 / p - 0.5/p ( t * mu/sqrt(p+mu^2) )^2
|
||||
do i = 0,n_pt_in
|
||||
d(i) = 0.d0
|
||||
enddo
|
||||
do i = 0,n_pt_in
|
||||
d1(i) = 0.d0
|
||||
enddo
|
||||
do i = 0,n_pt_in
|
||||
d2(i) = 0.d0
|
||||
enddo
|
||||
do i = 0,n_pt_in
|
||||
d3(i) = 0.d0
|
||||
enddo
|
||||
integer :: n_pt1,n_pt2,n_pt3,dim,i
|
||||
n_pt1 = n_pt_in
|
||||
n_pt2 = n_pt_in
|
||||
n_pt3 = n_pt_in
|
||||
a_x = power_A(1)
|
||||
b_x = power_B(1)
|
||||
call I_x1_pol_mult_mono_elec(a_x,b_x,R1x,R1xp,R2x,d1,n_pt1,n_pt_in)
|
||||
! print*,'passed the first I_x1'
|
||||
if(n_pt1<0)then
|
||||
n_pt_out = -1
|
||||
do i = 0,n_pt_in
|
||||
d(i) = 0.d0
|
||||
enddo
|
||||
return
|
||||
endif
|
||||
|
||||
R1x(0) = (P_center(2) - A_center(2))
|
||||
R1x(1) = 0.d0
|
||||
R1x(2) = -(P_center(2) - C_center(2))* p_new
|
||||
! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2
|
||||
R1xp(0) = (P_center(2) - B_center(2))
|
||||
R1xp(1) = 0.d0
|
||||
R1xp(2) =-(P_center(2) - C_center(2))* p_new
|
||||
!R1xp = (P_x - B_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2
|
||||
a_y = power_A(2)
|
||||
b_y = power_B(2)
|
||||
call I_x1_pol_mult_mono_elec(a_y,b_y,R1x,R1xp,R2x,d2,n_pt2,n_pt_in)
|
||||
! print*,'passed the second I_x1'
|
||||
if(n_pt2<0)then
|
||||
n_pt_out = -1
|
||||
do i = 0,n_pt_in
|
||||
d(i) = 0.d0
|
||||
enddo
|
||||
return
|
||||
endif
|
||||
|
||||
|
||||
R1x(0) = (P_center(3) - A_center(3))
|
||||
R1x(1) = 0.d0
|
||||
R1x(2) = -(P_center(3) - C_center(3))* p_new
|
||||
! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2
|
||||
R1xp(0) = (P_center(3) - B_center(3))
|
||||
R1xp(1) = 0.d0
|
||||
R1xp(2) =-(P_center(3) - C_center(3))* p_new
|
||||
!R2x = 0.5 / p - 0.5/p ( t * mu/sqrt(p+mu^2) )^2
|
||||
a_z = power_A(3)
|
||||
b_z = power_B(3)
|
||||
|
||||
! print*,'a_z = ',a_z
|
||||
! print*,'b_z = ',b_z
|
||||
call I_x1_pol_mult_mono_elec(a_z,b_z,R1x,R1xp,R2x,d3,n_pt3,n_pt_in)
|
||||
! print*,'passed the third I_x1'
|
||||
if(n_pt3<0)then
|
||||
n_pt_out = -1
|
||||
do i = 0,n_pt_in
|
||||
d(i) = 0.d0
|
||||
enddo
|
||||
return
|
||||
endif
|
||||
integer :: n_pt_tmp
|
||||
n_pt_tmp = 0
|
||||
call multiply_poly(d1,n_pt1,d2,n_pt2,d,n_pt_tmp)
|
||||
do i = 0,n_pt_tmp
|
||||
d1(i) = 0.d0
|
||||
enddo
|
||||
n_pt_out = 0
|
||||
call multiply_poly(d ,n_pt_tmp ,d3,n_pt3,d1,n_pt_out)
|
||||
do i = 0, n_pt_out
|
||||
d(i) = d1(i)
|
||||
enddo
|
||||
|
||||
end
|
||||
|
||||
|
||||
|
||||
|
||||
subroutine give_polynom_mult_center_mono_elec_erf(A_center,B_center,alpha,beta,power_A,power_B,C_center,n_pt_in,d,n_pt_out,mu_in)
|
||||
!!!! subroutine that returns the explicit polynom in term of the "t" variable of the following polynomw ::
|
||||
!!!! I_x1(a_x, d_x,p,q) * I_x1(a_y, d_y,p,q) * I_x1(a_z, d_z,p,q)
|
||||
!!!! it is for the nuclear electron atraction
|
||||
implicit none
|
||||
integer, intent(in) :: n_pt_in
|
||||
integer,intent(out) :: n_pt_out
|
||||
double precision, intent(in) :: A_center(3), B_center(3),C_center(3)
|
||||
double precision, intent(in) :: alpha,beta,mu_in
|
||||
integer, intent(in) :: power_A(3), power_B(3)
|
||||
integer :: a_x,b_x,a_y,b_y,a_z,b_z
|
||||
double precision :: d(0:n_pt_in)
|
||||
double precision :: d1(0:n_pt_in)
|
||||
double precision :: d2(0:n_pt_in)
|
||||
double precision :: d3(0:n_pt_in)
|
||||
double precision :: accu, pq_inv, p10_1, p10_2, p01_1, p01_2
|
||||
double precision :: p,P_center(3),rho,p_inv,p_inv_2
|
||||
accu = 0.d0
|
||||
!COMPTEUR irp_rdtsc1 = irp_rdtsc()
|
||||
ASSERT (n_pt_in > 1)
|
||||
p = alpha+beta
|
||||
p_inv = 1.d0/p
|
||||
p_inv_2 = 0.5d0/p
|
||||
do i =1, 3
|
||||
P_center(i) = (alpha * A_center(i) + beta * B_center(i)) * p_inv
|
||||
enddo
|
||||
|
||||
double precision :: R1x(0:2), B01(0:2), R1xp(0:2),R2x(0:2)
|
||||
R1x(0) = (P_center(1) - A_center(1))
|
||||
R1x(1) = 0.d0
|
||||
R1x(2) = -(P_center(1) - C_center(1))* mu_in**2 / (p+mu_in*mu_in)
|
||||
! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2
|
||||
R1xp(0) = (P_center(1) - B_center(1))
|
||||
R1xp(1) = 0.d0
|
||||
R1xp(2) =-(P_center(1) - C_center(1))* mu_in**2 / (p+mu_in*mu_in)
|
||||
!R1xp = (P_x - B_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2
|
||||
R2x(0) = p_inv_2
|
||||
R2x(1) = 0.d0
|
||||
R2x(2) = -p_inv_2* mu_in**2 / (p+mu_in*mu_in)
|
||||
!R2x = 0.5 / p - 0.5/p ( t * mu/sqrt(p+mu^2) )^2
|
||||
do i = 0,n_pt_in
|
||||
d(i) = 0.d0
|
||||
enddo
|
||||
do i = 0,n_pt_in
|
||||
d1(i) = 0.d0
|
||||
enddo
|
||||
do i = 0,n_pt_in
|
||||
d2(i) = 0.d0
|
||||
enddo
|
||||
do i = 0,n_pt_in
|
||||
d3(i) = 0.d0
|
||||
enddo
|
||||
integer :: n_pt1,n_pt2,n_pt3,dim,i
|
||||
n_pt1 = n_pt_in
|
||||
n_pt2 = n_pt_in
|
||||
n_pt3 = n_pt_in
|
||||
a_x = power_A(1)
|
||||
b_x = power_B(1)
|
||||
call I_x1_pol_mult_mono_elec(a_x,b_x,R1x,R1xp,R2x,d1,n_pt1,n_pt_in)
|
||||
! print*,'passed the first I_x1'
|
||||
if(n_pt1<0)then
|
||||
n_pt_out = -1
|
||||
do i = 0,n_pt_in
|
||||
d(i) = 0.d0
|
||||
enddo
|
||||
return
|
||||
endif
|
||||
|
||||
R1x(0) = (P_center(2) - A_center(2))
|
||||
R1x(1) = 0.d0
|
||||
R1x(2) = -(P_center(2) - C_center(2))* mu_in**2 / (p+mu_in*mu_in)
|
||||
! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2
|
||||
R1xp(0) = (P_center(2) - B_center(2))
|
||||
R1xp(1) = 0.d0
|
||||
R1xp(2) =-(P_center(2) - C_center(2))* mu_in**2 / (p+mu_in*mu_in)
|
||||
!R1xp = (P_x - B_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2
|
||||
a_y = power_A(2)
|
||||
b_y = power_B(2)
|
||||
call I_x1_pol_mult_mono_elec(a_y,b_y,R1x,R1xp,R2x,d2,n_pt2,n_pt_in)
|
||||
! print*,'passed the second I_x1'
|
||||
if(n_pt2<0)then
|
||||
n_pt_out = -1
|
||||
do i = 0,n_pt_in
|
||||
d(i) = 0.d0
|
||||
enddo
|
||||
return
|
||||
endif
|
||||
|
||||
|
||||
R1x(0) = (P_center(3) - A_center(3))
|
||||
R1x(1) = 0.d0
|
||||
R1x(2) = -(P_center(3) - C_center(3))* mu_in**2 / (p+mu_in*mu_in)
|
||||
! R1x = (P_x - A_x) - (P_x - C_x) ( t * mu/sqrt(p+mu^2) )^2
|
||||
R1xp(0) = (P_center(3) - B_center(3))
|
||||
R1xp(1) = 0.d0
|
||||
R1xp(2) =-(P_center(3) - C_center(3))* mu_in**2 / (p+mu_in*mu_in)
|
||||
!R2x = 0.5 / p - 0.5/p ( t * mu/sqrt(p+mu^2) )^2
|
||||
a_z = power_A(3)
|
||||
b_z = power_B(3)
|
||||
|
||||
! print*,'a_z = ',a_z
|
||||
! print*,'b_z = ',b_z
|
||||
call I_x1_pol_mult_mono_elec(a_z,b_z,R1x,R1xp,R2x,d3,n_pt3,n_pt_in)
|
||||
! print*,'passed the third I_x1'
|
||||
if(n_pt3<0)then
|
||||
n_pt_out = -1
|
||||
do i = 0,n_pt_in
|
||||
d(i) = 0.d0
|
||||
enddo
|
||||
return
|
||||
endif
|
||||
integer :: n_pt_tmp
|
||||
n_pt_tmp = 0
|
||||
call multiply_poly(d1,n_pt1,d2,n_pt2,d,n_pt_tmp)
|
||||
do i = 0,n_pt_tmp
|
||||
d1(i) = 0.d0
|
||||
enddo
|
||||
n_pt_out = 0
|
||||
call multiply_poly(d ,n_pt_tmp ,d3,n_pt3,d1,n_pt_out)
|
||||
do i = 0, n_pt_out
|
||||
d(i) = d1(i)
|
||||
enddo
|
||||
|
||||
end
|
||||
|
||||
@ -2,7 +2,14 @@
|
||||
Becke Numerical Grid
|
||||
====================
|
||||
|
||||
Lebedev-Laikov grids, using the code distributed through CCL (http://www.ccl.net/).
|
||||
This module contains all quantities needed to build the Becke's grid used in general for DFT integration. Note that it can be used for whatever integration in R^3 as long as the functions to be integrated are mostly concentrated near the atomic regions.
|
||||
|
||||
This grid is built as the reunion of a spherical grid around each atom. Each spherical grid contains a certain number of radial and angular points.
|
||||
|
||||
For a simple example of how to use the grid, see example.irp.f.
|
||||
|
||||
The spherical integration uses Lebedev-Laikov grids, which was used from the code distributed through CCL (http://www.ccl.net/).
|
||||
See next section for explanations and citation policies.
|
||||
|
||||
.. code-block:: text
|
||||
|
||||
|
||||
@ -1,3 +1,34 @@
|
||||
subroutine set_bit_to_integer(i_physical,key,Nint)
|
||||
use bitmasks
|
||||
BEGIN_DOC
|
||||
! set to 1 the bit number i_physical in the bitstring key
|
||||
END_DOC
|
||||
implicit none
|
||||
integer, intent(in) :: i_physical,Nint
|
||||
integer(bit_kind), intent(inout) :: key(Nint)
|
||||
integer :: k,j,i
|
||||
k = ishft(i_physical-1,-bit_kind_shift)+1
|
||||
j = i_physical-ishft(k-1,bit_kind_shift)-1
|
||||
key(k) = ibset(key(k),j)
|
||||
end
|
||||
|
||||
|
||||
subroutine clear_bit_to_integer(i_physical,key,Nint)
|
||||
use bitmasks
|
||||
BEGIN_DOC
|
||||
! set to 0 the bit number i_physical in the bitstring key
|
||||
END_DOC
|
||||
implicit none
|
||||
integer, intent(in) :: i_physical,Nint
|
||||
integer(bit_kind), intent(inout) :: key(Nint)
|
||||
integer :: k,j,i
|
||||
k = ishft(i_physical-1,-bit_kind_shift)+1
|
||||
j = i_physical-ishft(k-1,bit_kind_shift)-1
|
||||
key(k) = ibclr(key(k),j)
|
||||
end
|
||||
|
||||
|
||||
|
||||
subroutine bitstring_to_list( string, list, n_elements, Nint)
|
||||
use bitmasks
|
||||
implicit none
|
||||
|
||||
26
src/data_energy_and_density/EZFIO.cfg
Normal file
26
src/data_energy_and_density/EZFIO.cfg
Normal file
@ -0,0 +1,26 @@
|
||||
[data_energy_var]
|
||||
type: double precision
|
||||
doc: Variational energy computed with the wave function
|
||||
interface: ezfio, provider
|
||||
size: (determinants.n_states)
|
||||
|
||||
[data_energy_proj]
|
||||
type: double precision
|
||||
doc: Projected energy computed with the wave function
|
||||
interface: ezfio, provider
|
||||
size: (determinants.n_states)
|
||||
|
||||
[data_one_body_alpha_dm_mo]
|
||||
interface: ezfio, provider
|
||||
doc: Alpha one body density matrix on the MO basis computed with the wave function
|
||||
type: double precision
|
||||
size: (mo_basis.mo_tot_num,mo_basis.mo_tot_num,determinants.n_states)
|
||||
|
||||
|
||||
[data_one_body_beta_dm_mo]
|
||||
interface: ezfio, provider
|
||||
doc: Beta one body density matrix on the MO basis computed with the wave function
|
||||
type: double precision
|
||||
size: (mo_basis.mo_tot_num,mo_basis.mo_tot_num,determinants.n_states)
|
||||
|
||||
|
||||
1
src/data_energy_and_density/NEED
Normal file
1
src/data_energy_and_density/NEED
Normal file
@ -0,0 +1 @@
|
||||
determinants
|
||||
6
src/data_energy_and_density/README.rst
Normal file
6
src/data_energy_and_density/README.rst
Normal file
@ -0,0 +1,6 @@
|
||||
============
|
||||
Data energy and density
|
||||
============
|
||||
|
||||
|
||||
This module contains some global variables (such as densities and energies) which are stored in the EZFIO folder in a different place than determinants. This is used in practice to store density matrices which can be obtained from any methods, as long as they are stored in the same MO basis which is used for the calculations. In |RS-DFT| calculations, this can be done to perform damping on the density in order to speed up convergence.
|
||||
20
src/data_energy_and_density/save_one_body_dm.irp.f
Normal file
20
src/data_energy_and_density/save_one_body_dm.irp.f
Normal file
@ -0,0 +1,20 @@
|
||||
program save_one_body_dm
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! programs that computes the one body density on the mo basis for alpha and beta electrons from the wave function stored in the EZFIO folder, and then save it into the EZFIO folder data_energy_and_density.
|
||||
!
|
||||
! Then, the global variable data_one_body_alpha_dm_mo and data_one_body_beta_dm_mo will automatically read the density in a further calculation.
|
||||
!
|
||||
! This can be used to perform dampin on the density in RS-DFT calculation (see the density_for_dft module).
|
||||
END_DOC
|
||||
read_wf = .True.
|
||||
touch read_wf
|
||||
call routine
|
||||
|
||||
end
|
||||
|
||||
subroutine routine
|
||||
|
||||
call ezfio_set_data_energy_and_density_data_one_body_alpha_dm_mo(one_body_dm_mo_alpha)
|
||||
call ezfio_set_data_energy_and_density_data_one_body_beta_dm_mo(one_body_dm_mo_beta)
|
||||
end
|
||||
3
src/density_for_dft/NEED
Normal file
3
src/density_for_dft/NEED
Normal file
@ -0,0 +1,3 @@
|
||||
determinants
|
||||
dft_keywords
|
||||
data_energy_and_density
|
||||
4
src/density_for_dft/README.rst
Normal file
4
src/density_for_dft/README.rst
Normal file
@ -0,0 +1,4 @@
|
||||
==========
|
||||
DM_for_dft
|
||||
==========
|
||||
|
||||
77
src/density_for_dft/density_for_dft.irp.f
Normal file
77
src/density_for_dft/density_for_dft.irp.f
Normal file
@ -0,0 +1,77 @@
|
||||
BEGIN_PROVIDER [double precision, one_body_dm_mo_alpha_for_dft, (mo_tot_num,mo_tot_num, N_states)]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! density matrix for alpha electrons in the MO basis used for all DFT calculations based on the density
|
||||
END_DOC
|
||||
double precision :: delta_alpha(mo_tot_num,mo_tot_num,N_states)
|
||||
if(density_for_dft .EQ. "damping_rs_dft")then
|
||||
delta_alpha = one_body_dm_mo_alpha - data_one_body_alpha_dm_mo
|
||||
one_body_dm_mo_alpha_for_dft = data_one_body_alpha_dm_mo + damping_for_rs_dft * delta_alpha
|
||||
else if (density_for_dft .EQ. "input_density")then
|
||||
one_body_dm_mo_alpha_for_dft = data_one_body_alpha_dm_mo
|
||||
else if (density_for_dft .EQ. "WFT")then
|
||||
provide mo_coef
|
||||
one_body_dm_mo_alpha_for_dft = one_body_dm_mo_alpha
|
||||
endif
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [double precision, one_body_dm_mo_beta_for_dft, (mo_tot_num,mo_tot_num, N_states)]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! density matrix for beta electrons in the MO basis used for all DFT calculations based on the density
|
||||
END_DOC
|
||||
double precision :: delta_beta(mo_tot_num,mo_tot_num,N_states)
|
||||
if(density_for_dft .EQ. "damping_rs_dft")then
|
||||
delta_beta = one_body_dm_mo_beta - data_one_body_beta_dm_mo
|
||||
one_body_dm_mo_beta_for_dft = data_one_body_beta_dm_mo + damping_for_rs_dft * delta_beta
|
||||
else if (density_for_dft .EQ. "input_density")then
|
||||
one_body_dm_mo_beta_for_dft = data_one_body_beta_dm_mo
|
||||
else if (density_for_dft .EQ. "WFT")then
|
||||
provide mo_coef
|
||||
one_body_dm_mo_beta_for_dft = one_body_dm_mo_beta
|
||||
endif
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [double precision, one_body_dm_mo_for_dft, (mo_tot_num,mo_tot_num, N_states)]
|
||||
implicit none
|
||||
one_body_dm_mo_for_dft = one_body_dm_mo_beta_for_dft + one_body_dm_mo_alpha_for_dft
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [double precision, one_body_dm_average_mo_for_dft, (mo_tot_num,mo_tot_num)]
|
||||
implicit none
|
||||
integer :: i
|
||||
one_body_dm_average_mo_for_dft = 0.d0
|
||||
do i = 1, N_states
|
||||
one_body_dm_average_mo_for_dft(:,:) += one_body_dm_mo_for_dft(:,:,i) * state_average_weight(i)
|
||||
enddo
|
||||
END_PROVIDER
|
||||
|
||||
BEGIN_PROVIDER [ double precision, one_body_dm_alpha_ao_for_dft, (ao_num,ao_num,N_states) ]
|
||||
&BEGIN_PROVIDER [ double precision, one_body_dm_beta_ao_for_dft, (ao_num,ao_num,N_states) ]
|
||||
BEGIN_DOC
|
||||
! one body density matrix on the AO basis based on one_body_dm_mo_alpha_for_dft
|
||||
END_DOC
|
||||
implicit none
|
||||
integer :: i,j,k,l,istate
|
||||
double precision :: mo_alpha,mo_beta
|
||||
|
||||
one_body_dm_alpha_ao_for_dft = 0.d0
|
||||
one_body_dm_beta_ao_for_dft = 0.d0
|
||||
do k = 1, ao_num
|
||||
do l = 1, ao_num
|
||||
do i = 1, mo_tot_num
|
||||
do j = 1, mo_tot_num
|
||||
do istate = 1, N_states
|
||||
mo_alpha = one_body_dm_mo_alpha_for_dft(j,i,istate)
|
||||
mo_beta = one_body_dm_mo_beta_for_dft(j,i,istate)
|
||||
one_body_dm_alpha_ao_for_dft(l,k,istate) += mo_coef(k,i) * mo_coef(l,j) * mo_alpha
|
||||
one_body_dm_beta_ao_for_dft(l,k,istate) += mo_coef(k,i) * mo_coef(l,j) * mo_beta
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
enddo
|
||||
|
||||
END_PROVIDER
|
||||
|
||||
@ -249,9 +249,21 @@ integer function connected_to_ref(key,keys,Nint,N_past_in,Ndet)
|
||||
logical :: t
|
||||
double precision :: hij_elec
|
||||
|
||||
! output : 0 : not connected
|
||||
! i : connected to determinant i of the past
|
||||
! -i : is the ith determinant of the reference wf keys
|
||||
BEGIN_DOC
|
||||
! input : key : a given Slater determinant
|
||||
!
|
||||
! : keys: a list of Slater determinants
|
||||
!
|
||||
! : Ndet: the number of Slater determinants in keys
|
||||
!
|
||||
! : N_past_in the number of Slater determinants for the connectivity research
|
||||
!
|
||||
! output : 0 : key not connected to the N_past_in first Slater determinants in keys
|
||||
!
|
||||
! i : key is connected to determinant i of keys
|
||||
!
|
||||
! -i : key is the ith determinant of the reference wf keys
|
||||
END_DOC
|
||||
|
||||
ASSERT (Nint > 0)
|
||||
ASSERT (Nint == N_int)
|
||||
@ -349,9 +361,21 @@ integer function connected_to_ref_by_mono(key,keys,Nint,N_past_in,Ndet)
|
||||
logical :: t
|
||||
double precision :: hij_elec
|
||||
|
||||
! output : 0 : not connected
|
||||
! i : connected to determinant i of the past
|
||||
! -i : is the ith determinant of the refernce wf keys
|
||||
BEGIN_DOC
|
||||
! input : key : a given Slater determinant
|
||||
!
|
||||
! : keys: a list of Slater determinants
|
||||
!
|
||||
! : Ndet: the number of Slater determinants in keys
|
||||
!
|
||||
! : N_past_in the number of Slater determinants for the connectivity research
|
||||
!
|
||||
! output : 0 : key not connected by a MONO EXCITATION to the N_past_in first Slater determinants in keys
|
||||
!
|
||||
! i : key is connected by a MONO EXCITATION to determinant i of keys
|
||||
!
|
||||
! -i : key is the ith determinant of the reference wf keys
|
||||
END_DOC
|
||||
|
||||
ASSERT (Nint > 0)
|
||||
ASSERT (Nint == N_int)
|
||||
|
||||
@ -65,26 +65,3 @@ logical function is_spin_flip_possible(key_in,i_flip,ispin)
|
||||
endif
|
||||
end
|
||||
|
||||
subroutine set_bit_to_integer(i_physical,key,Nint)
|
||||
use bitmasks
|
||||
implicit none
|
||||
integer, intent(in) :: i_physical,Nint
|
||||
integer(bit_kind), intent(inout) :: key(Nint)
|
||||
integer :: k,j,i
|
||||
k = shiftr(i_physical-1,bit_kind_shift)+1
|
||||
j = i_physical-shiftl(k-1,bit_kind_shift)-1
|
||||
key(k) = ibset(key(k),j)
|
||||
end
|
||||
|
||||
|
||||
subroutine clear_bit_to_integer(i_physical,key,Nint)
|
||||
use bitmasks
|
||||
implicit none
|
||||
integer, intent(in) :: i_physical,Nint
|
||||
integer(bit_kind), intent(inout) :: key(Nint)
|
||||
integer :: k,j,i
|
||||
k = shiftr(i_physical-1,bit_kind_shift)+1
|
||||
j = i_physical-shiftl(k-1,bit_kind_shift)-1
|
||||
key(k) = ibclr(key(k),j)
|
||||
end
|
||||
|
||||
|
||||
@ -1,16 +0,0 @@
|
||||
subroutine apply_mono(i_hole,i_particle,ispin_excit,key_in,Nint)
|
||||
implicit none
|
||||
integer, intent(in) :: i_hole,i_particle,ispin_excit,Nint
|
||||
integer(bit_kind), intent(inout) :: key_in(Nint,2)
|
||||
integer :: k,j
|
||||
use bitmasks
|
||||
! hole
|
||||
k = shiftr(i_hole-1,bit_kind_shift)+1
|
||||
j = i_hole-shiftl(k-1,bit_kind_shift)-1
|
||||
key_in(k,ispin_excit) = ibclr(key_in(k,ispin_excit),j)
|
||||
|
||||
k = shiftr(i_particle-1,bit_kind_shift)+1
|
||||
j = i_particle-shiftl(k-1,bit_kind_shift)-1
|
||||
key_in(k,ispin_excit) = ibset(key_in(k,ispin_excit),j)
|
||||
|
||||
end
|
||||
37
src/dft_keywords/EZFIO.cfg
Normal file
37
src/dft_keywords/EZFIO.cfg
Normal file
@ -0,0 +1,37 @@
|
||||
[exchange_functional]
|
||||
type: character*(32)
|
||||
doc: name of the exchange functional
|
||||
interface: ezfio, provider, ocaml
|
||||
default: short_range_LDA
|
||||
|
||||
[correlation_functional]
|
||||
type: character*(32)
|
||||
doc: name of the correlation functional
|
||||
interface: ezfio, provider, ocaml
|
||||
default: short_range_LDA
|
||||
|
||||
[HF_exchange]
|
||||
type: double precision
|
||||
doc: Percentage of HF exchange in the DFT model
|
||||
interface: ezfio,provider,ocaml
|
||||
default: 0.
|
||||
|
||||
[mu_erf]
|
||||
type: double precision
|
||||
doc: cutting of the interaction in the range separated model
|
||||
interface: ezfio,provider,ocaml
|
||||
default: 0.5
|
||||
ezfio_name: mu_erf
|
||||
|
||||
|
||||
[density_for_dft]
|
||||
type: character*(32)
|
||||
doc: Type of density used for DFT calculation. If WFT it uses the density of the WFT stored in terms of determinants. If input_density it uses the one-body dm stored in data_.../ . If damping_rs_dft it uses the damping density between WFT and input_density
|
||||
interface: ezfio, provider, ocaml
|
||||
default: WFT
|
||||
|
||||
[damping_for_rs_dft]
|
||||
type: double precision
|
||||
doc: damping factor for the density used in RSFT.
|
||||
interface: ezfio,provider,ocaml
|
||||
default: 0.5
|
||||
1
src/dft_keywords/NEED
Normal file
1
src/dft_keywords/NEED
Normal file
@ -0,0 +1 @@
|
||||
ezfio_files
|
||||
7
src/dft_keywords/README.rst
Normal file
7
src/dft_keywords/README.rst
Normal file
@ -0,0 +1,7 @@
|
||||
============
|
||||
DFT Keywords
|
||||
============
|
||||
|
||||
This module contains all keywords which are related to a DFT calculation
|
||||
|
||||
|
||||
19
src/dft_keywords/keywords.irp.f
Normal file
19
src/dft_keywords/keywords.irp.f
Normal file
@ -0,0 +1,19 @@
|
||||
BEGIN_PROVIDER [ character*(32), DFT_TYPE]
|
||||
implicit none
|
||||
BEGIN_DOC
|
||||
! defines the type of DFT applied: LDA, GGA etc ...
|
||||
END_DOC
|
||||
logical :: is_lda
|
||||
if(correlation_functional.eq."None")then
|
||||
is_lda = (index(exchange_functional,"LDA") .ne. 0)
|
||||
else if(exchange_functional.eq."None")then
|
||||
is_lda = (index(correlation_functional,"LDA") .ne. 0)
|
||||
else
|
||||
is_lda = (index(correlation_functional,"LDA") .ne. 0) .and. (index(exchange_functional,"LDA") .ne. 0)
|
||||
endif
|
||||
if(is_lda)then
|
||||
DFT_TYPE = "LDA"
|
||||
else
|
||||
DFT_TYPE = "GGA"
|
||||
endif
|
||||
END_PROVIDER
|
||||
Loading…
Reference in New Issue
Block a user