diff --git a/Makefile b/Makefile index 24b8b55..cbc251d 100644 --- a/Makefile +++ b/Makefile @@ -1,8 +1,11 @@ IRPF90 = irpf90 #-a -d -FC = gfortran +FC = gfortran -g FCFLAGS= -O2 -ffree-line-length-none -I . -NINJA = ninja +FC = ifort -C -traceback -g +#FCFLAGS= -O2 -ffree-line-length-none -I . +#NINJA = ninja AR = ar +ARCHIVE= ar crs RANLIB = ranlib SRC= diff --git a/deriv_num.irp.f b/deriv_num.irp.f new file mode 100644 index 0000000..6bacd32 --- /dev/null +++ b/deriv_num.irp.f @@ -0,0 +1,47 @@ +program jastrow + implicit none + print *, 'Nabla J1' + integer::k + double precision :: j1, j2, j0, deriv, dt, lapl + dt = 1.d-4 + +BEGIN_TEMPLATE +lapl = 0.d0 +j0 = $X $Y +do k=1,3 + + elec_coord(1,k) -= dt + TOUCH elec_coord + j1 = $X $Y + + elec_coord(1,k) += 2.d0*dt + TOUCH elec_coord + j2 = $X $Y + + deriv = (j2 - j1)/(2.d0*dt) + lapl += (j2 - 2.d0*j0 + j1)/(dt*dt) + print *, 'deriv $X ' + print *, deriv + print *, $X_deriv_e(k,$Z) + print *, '' + + elec_coord(1,k) -= dt + TOUCH elec_coord + +enddo + print *, 'lapl $X ' + print *, lapl + print *, $X_deriv_e(4 ,$Z) + print *, '' + +SUBST [X,Y,Z] +factor_een ; ; 1 ;; +END_TEMPLATE +!rescale_een_e ; (1,3,1) ; 1,3,1 ;; +!rescale_een_n ; (1,1,2) ; 1,1,2 ;; +!rescale_een_e ; (1,2,2) ; 1,2,2 ;; +!elnuc_dist ; (1,1); 1,1 ;; +!elec_dist ; (1,2); 1,2 ;; + + +end program diff --git a/el_nuc_el.irp.f b/el_nuc_el.irp.f index 6938a27..95332bd 100644 --- a/el_nuc_el.irp.f +++ b/el_nuc_el.irp.f @@ -29,7 +29,8 @@ BEGIN_PROVIDER [ double precision, factor_een ] riam = rescale_een_n(i, a, m) rijk = rescale_een_e(i, j, k) factor_een = factor_een + & - rijk * (rial + rjal) * riam * rjam_cn +! rijk * (rial + rjal) * riam * rjam_cn + rijk * rjam_cn enddo enddo enddo @@ -68,7 +69,6 @@ BEGIN_PROVIDER [ double precision, factor_een_deriv_e, (4, nelec) ] cn = cord_vect_lkp(l, k, p, typenuc_arr(a)) do j = 1, nelec - factor_een_deriv_e(:, j) = 0.d0 rjal = rescale_een_n(j, a, l) rjam_cn = rescale_een_n(j, a, m) * cn @@ -83,7 +83,7 @@ BEGIN_PROVIDER [ double precision, factor_een_deriv_e, (4, nelec) ] rijk = rescale_een_e(i, j, k) do ii = 1, 4 - drijk(ii) = rescale_een_e_deriv_e(ii, i, j, k) + drijk(ii) = rescale_een_e_deriv_e(ii, j, i, k) enddo lap = 0.0d0 @@ -92,15 +92,23 @@ BEGIN_PROVIDER [ double precision, factor_een_deriv_e, (4, nelec) ] v1 = rijk * (rial + rjal) v2 = rjam_cn * riam - do ii = 1, 4 + do ii = 1, 3 d1 = drijk(ii) * (rial + rjal) + rijk * (rial + drjal(ii)) d2 = drjam_cn(ii) * riam factor_een_deriv_e(ii, j) = factor_een_deriv_e(ii, j) + & - v1 * d2 + d1 * v2 + x(ii) * lap - ! v(x) u''(x) + 2 * u'(x) v'(x) + u(x) v''(x) + v1 * d2 + d1 * v2 +! factor_een_deriv_e(ii, j) = factor_een_deriv_e(ii, j) + & +! drijk(ii) lap = lap + d1 * d2 enddo + ii = 4 + d1 = drijk(ii) * (rial + rjal) + rijk * (rial + drjal(ii)) + d2 = drjam_cn(ii) * riam + factor_een_deriv_e(ii, j) = factor_een_deriv_e(ii, j) + & + v1 * d2 + d1 * v2 + x(ii) * lap + ! v(x) u''(x) + 2 * u'(x) v'(x) + u(x) v''(x) + enddo enddo enddo @@ -108,6 +116,4 @@ BEGIN_PROVIDER [ double precision, factor_een_deriv_e, (4, nelec) ] enddo enddo - factor_een_deriv_e = 0.5d0 * factor_een_deriv_e - END_PROVIDER diff --git a/electrons.irp.f b/electrons.irp.f index 0011a93..4dfde8a 100644 --- a/electrons.irp.f +++ b/electrons.irp.f @@ -48,6 +48,7 @@ BEGIN_PROVIDER [ double precision, elec_dist, (nelec, nelec) ] z = elec_coord(i, 3) - elec_coord(j, 3) elec_dist(i, j) = dsqrt( x*x + y*y + z*z ) enddo +! elec_dist(j, j) = 1.d-10 enddo END_PROVIDER diff --git a/jastrow b/jastrow deleted file mode 100755 index b442acd..0000000 Binary files a/jastrow and /dev/null differ diff --git a/rescale.irp.f b/rescale.irp.f index 32fbf2c..f948fe3 100644 --- a/rescale.irp.f +++ b/rescale.irp.f @@ -188,7 +188,7 @@ BEGIN_PROVIDER [double precision, elec_dist_deriv_e, (4, nelec, nelec)] BEGIN_DOC ! Derivative of R_{ij} wrt x ! Dimensions 1-3 : dx, dy, dz - ! Dimension 4 : d2x + d2y + d2z + ! Dimension 4 : d2x + d2y + d2z = 2/rij END_DOC implicit none integer :: i, ii, j @@ -196,17 +196,14 @@ BEGIN_PROVIDER [double precision, elec_dist_deriv_e, (4, nelec, nelec)] do j = 1, nelec do i = 1, nelec - rij_inv = sign(1.0d0, dble(i - j)) / elec_dist(i, j) - lap = 0.0d0 + rij_inv = 1.0d0 / elec_dist(i, j) do ii = 1, 3 ! \frac{x-x0}{\sqrt{c+(x-x0)^2}} elec_dist_deriv_e(ii, i, j) = (elec_coord(i, ii) - elec_coord(j, ii)) * rij_inv - ! 1 / \sqrt{c+(x-x0)^2} - (x-x0)^2 /\left(c+(x-x0)^2\right)^{3/2} - lap = lap + rij_inv - elec_dist_deriv_e(ii, i, j) * elec_dist_deriv_e(ii, i, j) * rij_inv end do - elec_dist_deriv_e(4, i, j) = lap - elec_dist_deriv_e(:, i, i) = 0.0d0 + elec_dist_deriv_e(4, i, j) = 2.d0 * rij_inv end do + elec_dist_deriv_e(:, j, j) = 0.0d0 end do END_PROVIDER