diff --git a/src/JASTROW/jastrow_mu.irp.f b/src/JASTROW/jastrow_mu.irp.f new file mode 100644 index 0000000..13507f0 --- /dev/null +++ b/src/JASTROW/jastrow_mu.irp.f @@ -0,0 +1,120 @@ +! Mu Jastrow +! -------------- + +BEGIN_PROVIDER [ double precision , jast_elec_Mu_value, (elec_num_8) ] +implicit none + BEGIN_DOC +! J(i) = \sum_j a.rij/(1+b^2.rij) - \sum_A (a.riA/(1+a.riA))^2 + END_DOC + integer :: i,j + double precision :: a, b, rij, tmp + include '../constants.F' + double precision :: mu + mu = mu_erf + + do i=1,elec_num + jast_elec_Mu_value(i) = 0.d0 + enddo + + do j=1,elec_num + !DIR$ LOOP COUNT (50) + do i=1,elec_num + if(j==i)cycle + rij = elec_dist(i,j) + tmp = 0.5d0 * rij * (1.d0 - derf(mu*rij)) - 0.5d0/(dsqpi*mu) * dexp(-mu*mu*rij*rij) + jast_elec_Mu_value(i) += tmp + enddo + enddo + jast_elec_Mu_value = jast_elec_Mu_value * 0.5d0 ! symmetrization +END_PROVIDER + + BEGIN_PROVIDER [ double precision , jast_elec_Mu_grad_x, (elec_num_8) ] +&BEGIN_PROVIDER [ double precision , jast_elec_Mu_grad_y, (elec_num_8) ] +&BEGIN_PROVIDER [ double precision , jast_elec_Mu_grad_z, (elec_num_8) ] + implicit none + BEGIN_DOC +! Gradient of the Jastrow factor + END_DOC + + integer :: i,j + double precision :: a, b, rij, tmp, x, y, z + include '../constants.F' + double precision :: mu + mu = mu_erf + + do i=1,elec_num + jast_elec_Mu_grad_x(i) = 0.d0 + jast_elec_Mu_grad_y(i) = 0.d0 + jast_elec_Mu_grad_z(i) = 0.d0 + !DIR$ LOOP COUNT (100) + enddo + + ! (grad of J(r12) with respect to xi, yi, zi) + do i = 1, elec_num + do j = 1, elec_num + if(i==j)cycle + rij = elec_dist(j,i) + jast_elec_Mu_grad_x(i) += 0.5d0 * ( 1.d0 - derf(mu * rij) ) * elec_dist_inv(j,i) * (-1.d0) * elec_dist_vec_x(j,i) + jast_elec_Mu_grad_y(i) += 0.5d0 * ( 1.d0 - derf(mu * rij) ) * elec_dist_inv(j,i) * (-1.d0) * elec_dist_vec_y(j,i) + jast_elec_Mu_grad_z(i) += 0.5d0 * ( 1.d0 - derf(mu * rij) ) * elec_dist_inv(j,i) * (-1.d0) * elec_dist_vec_z(j,i) + enddo + enddo +END_PROVIDER + +BEGIN_PROVIDER [ double precision , jast_elec_Mu_lapl, (elec_num_8) ] + implicit none + BEGIN_DOC +! Laplacian of the Jastrow factor + END_DOC + + integer :: i,j + double precision :: a, b, rij, tmp, x, y, z + include '../constants.F' + double precision :: mu, x_ij, y_ij, z_ij, rij_inv + mu = mu_erf + + do i=1,elec_num + jast_elec_Mu_lapl(i) = 0.d0 + enddo + + do i=1, elec_num + do j=1, elec_num + if(j==i)cycle + rij = elec_dist(j,i) + rij_inv = elec_dist_inv(j,i) + x_ij = elec_dist_vec_x(j,i) + y_ij = elec_dist_vec_y(j,i) + z_ij = elec_dist_vec_z(j,i) + + jast_elec_Mu_lapl(i) += (1.d0 - derf(mu*rij))*elec_dist_inv(j,i) - mu/dsqpi * dexp(-mu*mu*rij*rij) + enddo + enddo +END_PROVIDER + +BEGIN_PROVIDER [double precision, mu_erf ] + implicit none + mu_erf = 0.5d0 +END_PROVIDER + + BEGIN_PROVIDER [double precision, grad_j_mu_x,(elec_num, elec_num)] +&BEGIN_PROVIDER [double precision, grad_j_mu_y,(elec_num, elec_num)] +&BEGIN_PROVIDER [double precision, grad_j_mu_z,(elec_num, elec_num)] + implicit none + integer :: i,j + double precision :: rij, mu,scal + mu = mu_erf + grad_j_mu_x = 0.d0 + grad_j_mu_y = 0.d0 + grad_j_mu_z = 0.d0 + do j = 1, elec_num + do i = 1, elec_num + if(i==j)cycle + rij = elec_dist(i,j) + scal = 0.5d0 * ( 1.d0 - derf(mu * rij) ) * elec_dist_inv(i,j) + grad_j_mu_x(i,j) = (elec_coord_transp(1,i) - elec_coord_transp(1,j)) * scal + grad_j_mu_y(i,j) = (elec_coord_transp(2,i) - elec_coord_transp(2,j)) * scal + grad_j_mu_z(i,j) = (elec_coord_transp(3,i) - elec_coord_transp(3,j)) * scal + enddo + enddo + +END_PROVIDER diff --git a/src/PROPERTIES/properties_mu.irp.f b/src/PROPERTIES/properties_mu.irp.f new file mode 100644 index 0000000..c477d0d --- /dev/null +++ b/src/PROPERTIES/properties_mu.irp.f @@ -0,0 +1,143 @@ +BEGIN_PROVIDER [ double precision, Energy_mu ] + implicit none + BEGIN_DOC + ! E mu + END_DOC + + integer :: i + energy_mu = E_nucl + !DIR$ VECTOR ALIGNED + !DIR$ LOOP COUNT(200) + do i=1,elec_num + energy_mu += E_kin_elec(i) + enddo + energy_mu += Eff_pot_mu + eff_pot_deriv_mu + E_nucl_elec - three_body_mu + + energy_mu_min = min(energy_mu_min,energy_mu) + energy_mu_max = max(energy_mu_max,energy_mu) + SOFT_TOUCH energy_mu_min energy_mu_max +END_PROVIDER + +BEGIN_PROVIDER [double precision, E_nucl_elec] + implicit none + integer :: i,j + E_nucl_elec = 0.d0 + do i = 1, elec_num +! E_nucl_elec += E_pot_elec_one(i) + E_pot_elec_two(i) + E_nucl_elec += E_pot_elec_one(i) + enddo + E_nucl_elec_min = min(E_nucl_elec_min,E_nucl_elec) + E_nucl_elec_max = max(E_nucl_elec_max,E_nucl_elec) +END_PROVIDER + + + BEGIN_PROVIDER [double precision, Eff_pot_mu_elec, (elec_num)] +&BEGIN_PROVIDER [double precision, Eff_pot_mu_elec_simple, (elec_num)] + implicit none + include '../constants.F' + integer :: i,j + double precision :: rij, mu + mu = mu_erf + Eff_pot_mu_elec = 0.d0 + do i=1,elec_num + !DIR$ VECTOR ALIGNED + !DIR$ LOOP COUNT(50) + do j=1,elec_num + rij = elec_dist(j,i) + if(i==j)cycle + Eff_pot_mu_elec(i) = Eff_pot_mu_elec(i) + 0.5d0 * derf(mu * rij) * elec_dist_inv(j,i) + Eff_pot_mu_elec(i) = Eff_pot_mu_elec(i) + 0.5d0 * mu/dsqpi * dexp(-mu*mu*rij*rij) + Eff_pot_mu_elec_simple(i) = Eff_pot_mu_elec(i) + Eff_pot_mu_elec(i) = Eff_pot_mu_elec(i) + 0.5d0 * (- 0.25d0 * (1.d0 - derf(mu*rij))**2.d0 ) + + enddo + enddo + +END_PROVIDER + + BEGIN_PROVIDER [double precision, Eff_pot_mu ] + implicit none + include '../constants.F' + integer :: i + Eff_pot_mu = 0.d0 + do i=1,elec_num + Eff_pot_mu += eff_pot_mu_elec(i) + enddo + Eff_pot_mu_min = min(Eff_pot_mu_min,Eff_pot_mu) + Eff_pot_mu_max = max(Eff_pot_mu_max,Eff_pot_mu) + SOFT_TOUCH Eff_pot_mu_min Eff_pot_mu_max + +END_PROVIDER + +BEGIN_PROVIDER [double precision, Eff_pot_mu_simple ] + implicit none + include '../constants.F' + integer :: i + Eff_pot_mu_simple = 0.d0 + do i=1,elec_num + Eff_pot_mu_simple += Eff_pot_mu_elec_simple(i) + enddo + Eff_pot_mu_simple_min = min(Eff_pot_mu_simple_min,Eff_pot_mu_simple) + Eff_pot_mu_simple_max = max(Eff_pot_mu_simple_max,Eff_pot_mu_simple) + SOFT_TOUCH Eff_pot_mu_simple_min Eff_pot_mu_simple_max + +END_PROVIDER + +BEGIN_PROVIDER [double precision, eff_pot_deriv_mu_elec, (elec_num) ] + implicit none + integer :: i,j + double precision :: rij, mu + mu = mu_erf + eff_pot_deriv_mu_elec = 0.d0 + do i = 1, elec_num + do j = 1, elec_num + if(i==j)cycle + rij = elec_dist(i,j) + eff_pot_deriv_mu_elec(i) += 0.5d0 * ( derf(mu * rij) - 1.d0 ) * elec_dist_inv(j,i) & + * ( - elec_dist_vec_x(j,i) * psidet_grad_lapl(1,i) & + - elec_dist_vec_y(j,i) * psidet_grad_lapl(2,i) & + - elec_dist_vec_z(j,i) * psidet_grad_lapl(3,i) ) * psidet_inv + enddo + enddo + +END_PROVIDER + + +BEGIN_PROVIDER [double precision, three_body_mu ] + implicit none + integer :: i,j,k + three_body_mu = 0.d0 + do i = 1, elec_num + do j = i+1, elec_num + do k = j+1, elec_num + three_body_mu += grad_j_mu_x(i,j) * grad_j_mu_x(i,k) + three_body_mu += grad_j_mu_y(i,j) * grad_j_mu_y(i,k) + three_body_mu += grad_j_mu_z(i,j) * grad_j_mu_z(i,k) + + three_body_mu += grad_j_mu_x(j,i) * grad_j_mu_x(j,k) + three_body_mu += grad_j_mu_y(j,i) * grad_j_mu_y(j,k) + three_body_mu += grad_j_mu_z(j,i) * grad_j_mu_z(j,k) + + three_body_mu += grad_j_mu_x(k,i) * grad_j_mu_x(k,j) + three_body_mu += grad_j_mu_y(k,i) * grad_j_mu_y(k,j) + three_body_mu += grad_j_mu_z(k,i) * grad_j_mu_z(k,j) + enddo + enddo + enddo + three_body_mu_min = min(three_body_mu_min,three_body_mu) + three_body_mu_max = max(three_body_mu_max,three_body_mu) + SOFT_TOUCH three_body_mu_min three_body_mu_max +END_PROVIDER + +BEGIN_PROVIDER [double precision, eff_pot_deriv_mu] + implicit none + integer :: i + eff_pot_deriv_mu = 0.d0 + do i = 1, elec_num + eff_pot_deriv_mu += eff_pot_deriv_mu_elec(i) + enddo + eff_pot_deriv_mu_min = min(eff_pot_deriv_mu_min,eff_pot_deriv_mu) + eff_pot_deriv_mu_max = max(eff_pot_deriv_mu_max,eff_pot_deriv_mu) + SOFT_TOUCH eff_pot_deriv_mu_min eff_pot_deriv_mu_max + +END_PROVIDER diff --git a/src/ezfio_interface.irp.f b/src/ezfio_interface.irp.f index 09f032f..1b9dc37 100644 --- a/src/ezfio_interface.irp.f +++ b/src/ezfio_interface.irp.f @@ -18,6 +18,7 @@ data = [ \ ("ao_basis_ao_power" , "integer" , "(ao_num,3)" ), ("ao_basis_ao_expo" , "real" , "(ao_num,ao_prim_num_max)" ), ("ao_basis_ao_coef" , "real" , "(ao_num,ao_prim_num_max)" ), +("jastrow_mu_erf" , "real" , "" ), ("jastrow_jast_a_up_up" , "real" , "" ), ("jastrow_jast_a_up_dn" , "real" , "" ), ("jastrow_jast_b_up_up" , "real" , "" ),