From 2ef05d01c9033abccc5f8f7166b33418741b9420 Mon Sep 17 00:00:00 2001 From: Abdallah Ammar Date: Mon, 19 Jun 2023 23:39:53 +0200 Subject: [PATCH] j1b_type 4/104 modif --- src/ao_many_one_e_ints/listj1b.irp.f | 14 +++--- src/non_h_ints_mu/j12_nucl_utils.irp.f | 8 ++-- src/non_h_ints_mu/jast_deriv.irp.f | 6 +-- src/tc_keywords/EZFIO.cfg | 6 +++ src/tc_keywords/j1b_pen.irp.f | 64 +++++++++++++++++++++----- 5 files changed, 73 insertions(+), 25 deletions(-) diff --git a/src/ao_many_one_e_ints/listj1b.irp.f b/src/ao_many_one_e_ints/listj1b.irp.f index 02963605..33ca8085 100644 --- a/src/ao_many_one_e_ints/listj1b.irp.f +++ b/src/ao_many_one_e_ints/listj1b.irp.f @@ -62,6 +62,7 @@ END_PROVIDER double precision :: tmp_cent_x, tmp_cent_y, tmp_cent_z provide j1b_pen + provide j1b_pen_coef List_all_comb_b2_coef = 0.d0 List_all_comb_b2_expo = 0.d0 @@ -127,8 +128,8 @@ END_PROVIDER List_all_comb_b2_expo( 1) = 0.d0 List_all_comb_b2_cent(1:3,1) = 0.d0 do i = 1, nucl_num - List_all_comb_b2_coef( i+1) = -1.d0 - List_all_comb_b2_expo( i+1) = j1b_pen( i) + List_all_comb_b2_coef( i+1) = -1.d0 * j1b_pen_coef(i) + List_all_comb_b2_expo( i+1) = j1b_pen(i) List_all_comb_b2_cent(1,i+1) = nucl_coord(i,1) List_all_comb_b2_cent(2,i+1) = nucl_coord(i,2) List_all_comb_b2_cent(3,i+1) = nucl_coord(i,3) @@ -225,6 +226,7 @@ END_PROVIDER double precision :: dx, dy, dz, r2 provide j1b_pen + provide j1b_pen_coef List_all_comb_b3_coef = 0.d0 List_all_comb_b3_expo = 0.d0 @@ -296,8 +298,8 @@ END_PROVIDER do i = 1, nucl_num ii = ii + 1 - List_all_comb_b3_coef( ii) = -2.d0 - List_all_comb_b3_expo( ii) = j1b_pen( i) + List_all_comb_b3_coef( ii) = -2.d0 * j1b_pen_coef(i) + List_all_comb_b3_expo( ii) = j1b_pen(i) List_all_comb_b3_cent(1,ii) = nucl_coord(i,1) List_all_comb_b3_cent(2,ii) = nucl_coord(i,2) List_all_comb_b3_cent(3,ii) = nucl_coord(i,3) @@ -305,7 +307,7 @@ END_PROVIDER do i = 1, nucl_num ii = ii + 1 - List_all_comb_b3_coef( ii) = 1.d0 + List_all_comb_b3_coef( ii) = 1.d0 * j1b_pen_coef(i) * j1b_pen_coef(i) List_all_comb_b3_expo( ii) = 2.d0 * j1b_pen(i) List_all_comb_b3_cent(1,ii) = nucl_coord(i,1) List_all_comb_b3_cent(2,ii) = nucl_coord(i,2) @@ -337,7 +339,7 @@ END_PROVIDER ii = ii + 1 ! x 2 to avoid doing integrals twice - List_all_comb_b3_coef( ii) = 2.d0 * dexp(-tmp1*tmp2*tmp4*r2) + List_all_comb_b3_coef( ii) = 2.d0 * dexp(-tmp1*tmp2*tmp4*r2) * j1b_pen_coef(i) * j1b_pen_coef(j) List_all_comb_b3_expo( ii) = tmp3 List_all_comb_b3_cent(1,ii) = tmp4 * (tmp1 * xi + tmp2 * xj) List_all_comb_b3_cent(2,ii) = tmp4 * (tmp1 * yi + tmp2 * yj) diff --git a/src/non_h_ints_mu/j12_nucl_utils.irp.f b/src/non_h_ints_mu/j12_nucl_utils.irp.f index 9b91a8ed..ac077fe0 100644 --- a/src/non_h_ints_mu/j12_nucl_utils.irp.f +++ b/src/non_h_ints_mu/j12_nucl_utils.irp.f @@ -35,7 +35,7 @@ BEGIN_PROVIDER [ double precision, v_1b, (n_points_final_grid)] elseif(j1b_type .eq. 4) then - ! v(r) = 1 - \sum_{a} \exp(-\alpha_a (r - r_a)^2) + ! v(r) = 1 - \sum_{a} \beta_a \exp(-\alpha_a (r - r_a)^2) do ipoint = 1, n_points_final_grid @@ -51,7 +51,7 @@ BEGIN_PROVIDER [ double precision, v_1b, (n_points_final_grid)] dz = z - nucl_coord(j,3) d = dx*dx + dy*dy + dz*dz - fact_r = fact_r - dexp(-a*d) + fact_r = fact_r - j1b_pen_coef(j) * dexp(-a*d) enddo v_1b(ipoint) = fact_r @@ -125,7 +125,7 @@ BEGIN_PROVIDER [double precision, v_1b_grad, (3, n_points_final_grid)] elseif(j1b_type .eq. 4) then - ! v(r) = 1 - \sum_{a} \exp(-\alpha_a (r - r_a)^2) + ! v(r) = 1 - \sum_{a} \beta_a \exp(-\alpha_a (r - r_a)^2) do ipoint = 1, n_points_final_grid @@ -144,7 +144,7 @@ BEGIN_PROVIDER [double precision, v_1b_grad, (3, n_points_final_grid)] r2 = dx*dx + dy*dy + dz*dz a = j1b_pen(j) - e = a * dexp(-a * r2) + e = a * j1b_pen_coef(j) * dexp(-a * r2) ax_der += e * dx ay_der += e * dy diff --git a/src/non_h_ints_mu/jast_deriv.irp.f b/src/non_h_ints_mu/jast_deriv.irp.f index 5e99600e..bd7ff6b7 100644 --- a/src/non_h_ints_mu/jast_deriv.irp.f +++ b/src/non_h_ints_mu/jast_deriv.irp.f @@ -296,7 +296,7 @@ double precision function j1b_nucl(r) d = ( (r(1) - nucl_coord(i,1)) * (r(1) - nucl_coord(i,1)) & + (r(2) - nucl_coord(i,2)) * (r(2) - nucl_coord(i,2)) & + (r(3) - nucl_coord(i,3)) * (r(3) - nucl_coord(i,3)) ) - j1b_nucl = j1b_nucl - dexp(-a*d) + j1b_nucl = j1b_nucl - j1b_pen_coef(i) * dexp(-a*d) enddo elseif((j1b_type .eq. 5) .or. (j1b_type .eq. 105)) then @@ -363,7 +363,7 @@ double precision function j1b_nucl_square(r) d = ( (r(1) - nucl_coord(i,1)) * (r(1) - nucl_coord(i,1)) & + (r(2) - nucl_coord(i,2)) * (r(2) - nucl_coord(i,2)) & + (r(3) - nucl_coord(i,3)) * (r(3) - nucl_coord(i,3)) ) - j1b_nucl_square = j1b_nucl_square - dexp(-a*d) + j1b_nucl_square = j1b_nucl_square - j1b_pen_coef(i) * dexp(-a*d) enddo j1b_nucl_square = j1b_nucl_square * j1b_nucl_square @@ -475,7 +475,7 @@ subroutine grad1_j1b_nucl(r, grad) y = r(2) - nucl_coord(i,2) z = r(3) - nucl_coord(i,3) d = x*x + y*y + z*z - e = a * dexp(-a*d) + e = a * j1b_pen_coef(i) * dexp(-a*d) fact_x += e * x fact_y += e * y diff --git a/src/tc_keywords/EZFIO.cfg b/src/tc_keywords/EZFIO.cfg index a69f5bac..ea1503c3 100644 --- a/src/tc_keywords/EZFIO.cfg +++ b/src/tc_keywords/EZFIO.cfg @@ -130,6 +130,12 @@ doc: exponents of the 1-body Jastrow interface: ezfio size: (nuclei.nucl_num) +[j1b_pen_coef] +type: double precision +doc: coefficients of the 1-body Jastrow +interface: ezfio +size: (nuclei.nucl_num) + [j1b_coeff] type: double precision doc: coeff of the 1-body Jastrow diff --git a/src/tc_keywords/j1b_pen.irp.f b/src/tc_keywords/j1b_pen.irp.f index 57250b52..3f1eb8ac 100644 --- a/src/tc_keywords/j1b_pen.irp.f +++ b/src/tc_keywords/j1b_pen.irp.f @@ -1,17 +1,22 @@ ! --- -BEGIN_PROVIDER [ double precision, j1b_pen, (nucl_num) ] + BEGIN_PROVIDER [ double precision, j1b_pen , (nucl_num) ] +&BEGIN_PROVIDER [ double precision, j1b_pen_coef, (nucl_num) ] BEGIN_DOC - ! exponents of the 1-body Jastrow + ! parameters of the 1-body Jastrow END_DOC implicit none logical :: exists + integer :: i + integer :: ierr PROVIDE ezfio_filename + ! --- + if (mpi_master) then call ezfio_has_tc_keywords_j1b_pen(exists) endif @@ -23,7 +28,6 @@ BEGIN_PROVIDER [ double precision, j1b_pen, (nucl_num) ] IRP_IF MPI include 'mpif.h' - integer :: ierr call MPI_BCAST(j1b_pen, (nucl_num), MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ierr) if (ierr /= MPI_SUCCESS) then stop 'Unable to read j1b_pen with MPI' @@ -31,7 +35,6 @@ BEGIN_PROVIDER [ double precision, j1b_pen, (nucl_num) ] IRP_ENDIF if (exists) then - if (mpi_master) then write(6,'(A)') '.. >>>>> [ IO READ: j1b_pen ] <<<<< ..' call ezfio_get_tc_keywords_j1b_pen(j1b_pen) @@ -42,19 +45,55 @@ BEGIN_PROVIDER [ double precision, j1b_pen, (nucl_num) ] endif IRP_ENDIF endif - else - - integer :: i do i = 1, nucl_num j1b_pen(i) = 1d5 enddo - endif - print*,'parameters for nuclei jastrow' - do i = 1, nucl_num - print*,'i,Z,j1b_pen(i)',i,nucl_charge(i),j1b_pen(i) - enddo + + ! --- + + if (mpi_master) then + call ezfio_has_tc_keywords_j1b_pen_coef(exists) + endif + + IRP_IF MPI_DEBUG + print *, irp_here, mpi_rank + call MPI_BARRIER(MPI_COMM_WORLD, ierr) + IRP_ENDIF + + IRP_IF MPI + include 'mpif.h' + call MPI_BCAST(j1b_pen_coef, (nucl_num), MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ierr) + if (ierr /= MPI_SUCCESS) then + stop 'Unable to read j1b_pen_coef with MPI' + endif + IRP_ENDIF + + if (exists) then + if (mpi_master) then + write(6,'(A)') '.. >>>>> [ IO READ: j1b_pen_coef ] <<<<< ..' + call ezfio_get_tc_keywords_j1b_pen_coef(j1b_pen_coef) + IRP_IF MPI + call MPI_BCAST(j1b_pen_coef, (nucl_num), MPI_DOUBLE_PRECISION, 0, MPI_COMM_WORLD, ierr) + if (ierr /= MPI_SUCCESS) then + stop 'Unable to read j1b_pen_coef with MPI' + endif + IRP_ENDIF + endif + else + do i = 1, nucl_num + j1b_pen_coef(i) = 1d0 + enddo + endif + + ! --- + + print *, ' parameters for nuclei jastrow' + print *, ' i, Z, j1b_pen, j1b_pen_coef' + do i = 1, nucl_num + print *, i, nucl_charge(i), j1b_pen(i), j1b_pen_coef(i) + enddo END_PROVIDER @@ -114,3 +153,4 @@ BEGIN_PROVIDER [ double precision, j1b_coeff, (nucl_num) ] END_PROVIDER ! --- +