From 67d2301fde87edbba1e803a7057e5bd1d9d01cb4 Mon Sep 17 00:00:00 2001 From: eginer Date: Wed, 11 Dec 2024 14:30:42 +0100 Subject: [PATCH] beginning the fit on 1s functions --- plugins/local/ao_extra_basis/EZFIO.cfg | 6 + plugins/local/ao_extra_basis/aos_transp.irp.f | 68 ++++++++ .../local/ao_extra_basis/fit_1s_basis.irp.f | 8 + .../local/ao_extra_basis/prov_fit_1s.irp.f | 159 ++++++++++++++++++ plugins/local/extra_basis_int/NEED | 1 + 5 files changed, 242 insertions(+) create mode 100644 plugins/local/ao_extra_basis/aos_transp.irp.f create mode 100644 plugins/local/ao_extra_basis/fit_1s_basis.irp.f create mode 100644 plugins/local/ao_extra_basis/prov_fit_1s.irp.f diff --git a/plugins/local/ao_extra_basis/EZFIO.cfg b/plugins/local/ao_extra_basis/EZFIO.cfg index 8b3a8667..9a73f078 100644 --- a/plugins/local/ao_extra_basis/EZFIO.cfg +++ b/plugins/local/ao_extra_basis/EZFIO.cfg @@ -1,3 +1,9 @@ +[ao_extra_center] +type: double precision +doc: parameter for off-centering s functions to mimick p functions +interface: ezfio, provider +default:0.01 + [ao_extra_basis] type: character*(256) doc: Name of the |ao_extra| basis set diff --git a/plugins/local/ao_extra_basis/aos_transp.irp.f b/plugins/local/ao_extra_basis/aos_transp.irp.f new file mode 100644 index 00000000..ed34835b --- /dev/null +++ b/plugins/local/ao_extra_basis/aos_transp.irp.f @@ -0,0 +1,68 @@ + +! --- + +BEGIN_PROVIDER [ integer, Nucl_ao_extras_transposed, (N_ao_extras_max,nucl_num)] + + BEGIN_DOC + ! List of ao_extras attached on each atom + END_DOC + + implicit none + integer :: i + integer, allocatable :: nucl_tmp(:) + + allocate(nucl_tmp(nucl_num)) + nucl_tmp = 0 + do i = 1, ao_extra_num + nucl_tmp(ao_extra_nucl(i)) += 1 + Nucl_ao_extras_transposed(nucl_tmp(ao_extra_nucl(i)),ao_extra_nucl(i)) = i + enddo + deallocate(nucl_tmp) + +END_PROVIDER + +! --- + +BEGIN_PROVIDER [double precision, ao_extra_expo_ordered_transp_per_nucl, (ao_extra_prim_num_max,N_ao_extras_max,nucl_num) ] + implicit none + integer :: i,j,k,l + do i = 1, nucl_num + do j = 1,Nucl_N_ao_extras(i) + k = Nucl_ao_extras_transposed(j,i) + do l = 1, ao_extra_prim_num(k) + ao_extra_expo_ordered_transp_per_nucl(l,j,i) = ao_extra_expo_ordered_transp(l,k) + enddo + enddo + enddo + +END_PROVIDER + + +BEGIN_PROVIDER [ integer, ao_extra_power_ordered_transp_per_nucl, (3,N_ao_extras_max,nucl_num) ] + implicit none + integer :: i,j,k,l + do i = 1, nucl_num + do j = 1,Nucl_N_ao_extras(i) + k = Nucl_ao_extras_transposed(j,i) + do l = 1, 3 + ao_extra_power_ordered_transp_per_nucl(l,j,i) = ao_extra_power(k,l) + enddo + enddo + enddo + +END_PROVIDER + +BEGIN_PROVIDER [ double precision, ao_extra_coef_normalized_ordered_transp_per_nucl, (ao_extra_prim_num_max,N_ao_extras_max,nucl_num) ] + implicit none + integer :: i,j,k,l + do i = 1, nucl_num + do j = 1,Nucl_N_ao_extras(i) + k = Nucl_ao_extras_transposed(j,i) + do l = 1, ao_extra_prim_num(k) + ao_extra_coef_normalized_ordered_transp_per_nucl(l,j,i) = ao_extra_coef_normalized_ordered_transp(l,k) + enddo + enddo + enddo + +END_PROVIDER + diff --git a/plugins/local/ao_extra_basis/fit_1s_basis.irp.f b/plugins/local/ao_extra_basis/fit_1s_basis.irp.f new file mode 100644 index 00000000..f544ad83 --- /dev/null +++ b/plugins/local/ao_extra_basis/fit_1s_basis.irp.f @@ -0,0 +1,8 @@ +program fit_1s_basis + implicit none + provide lmax_too_big + integer :: i +! print*,'n_func_tot', n_func_tot + provide new_nucl_coord + +end diff --git a/plugins/local/ao_extra_basis/prov_fit_1s.irp.f b/plugins/local/ao_extra_basis/prov_fit_1s.irp.f new file mode 100644 index 00000000..cba0bb05 --- /dev/null +++ b/plugins/local/ao_extra_basis/prov_fit_1s.irp.f @@ -0,0 +1,159 @@ + BEGIN_PROVIDER [ integer, n_func_tot] + implicit none + BEGIN_DOC + ! n_func_tot :: total number of functions in the fitted basis set + ! + ! returned in an uncontracted way + END_DOC + integer :: i,prefact + n_func_tot = 0 + print*,'n_func_tot ' + do i = 1, ao_num + if(ao_l(i) == 0)then + prefact = 1 ! s functions + else + ! p functions are fitted with 2 functions + ! d functions are fitted with 4 functions etc ... + prefact=2*ao_l(i) + endif + n_func_tot += prefact * ao_prim_num(i) + enddo +END_PROVIDER + +BEGIN_PROVIDER [ integer, n_prim_tot_orig] + implicit none + integer :: i + n_prim_tot_orig = 0 + do i = 1, ao_num + n_prim_tot_orig += ao_prim_num(i) + enddo +END_PROVIDER + + +BEGIN_PROVIDER [ logical, lmax_too_big] + implicit none + if (ao_l_max.gt.1)then + lmax_too_big = .True. + else + lmax_too_big = .False. + endif + if(lmax_too_big)then + print*,'STOPPING !! lmax is larger than 1 !' + print*,'Cannot yet fit with 1s functions ...' + stop + endif +END_PROVIDER + + BEGIN_PROVIDER [ integer, n_2p_func_orig] +&BEGIN_PROVIDER [ integer, n_2p_func_tot] + implicit none + integer :: i + BEGIN_DOC + ! n_2p_func_orig :: number of 2p functions in the original basis + ! + ! n_2p_func_tot :: total number of p functions in the fitted basis + END_DOC + n_2p_func_orig= 0 + n_2p_func_tot = 0 + do i = 1, ao_num + if(ao_l(i)==1)then + n_2p_func_orig+= 1 + n_2p_func_tot += ao_prim_num(i) * 2 + endif + enddo +END_PROVIDER + +BEGIN_PROVIDER [ integer, list_2p_functions, (n_2p_func_orig)] + implicit none + BEGIN_DOC + ! list of 2p functions in the original basis + END_DOC + integer :: i,j + j=0 + do i = 1, ao_num + if(ao_l(i)==1)then + j+=1 + list_2p_functions(j) = i + endif + enddo +END_PROVIDER + +BEGIN_PROVIDER [ integer, new_nucl_num] + implicit none + new_nucl_num = n_func_tot +END_PROVIDER + +BEGIN_PROVIDER [ double precision, new_ao_expo_1s , (n_func_tot) ] + implicit none + integer :: i,j,ii,i_ao,k + k = 0 + do i = 1, nucl_num + do ii = 1, Nucl_N_Aos(i) + i_ao = nucl_aos_transposed(ii,i) + if(ao_l(i_ao)==0)then + do j = 1, ao_prim_num(i_ao) + k+=1 + new_ao_expo_1s(k)= ao_expo(i_ao,j) + enddo + else if(ao_l(i_ao)==1)then + ! for 'p' functions + ! you replace the function by 2 functions 's' + do j = 1, ao_prim_num(i_ao) + k+=1 + new_ao_expo_1s(k)= ao_expo(i_ao,j) + k+=1 + new_ao_expo_1s(k)= ao_expo(i_ao,j) + enddo + else + print*,'WARNING ! Lmax value not implemented yet !' + print*,'stopping ...' + stop + endif + enddo + enddo + if(k.ne.n_func_tot)then + print*,'pb !!! k NE n_func_tot !!' + print*,k,n_func_tot + stop + endif + +END_PROVIDER + +BEGIN_PROVIDER [ double precision, new_nucl_coord, (3,new_nucl_num)] + implicit none + integer :: i,ii,j,i_ao,k + k = 0 + do i = 1, nucl_num + do ii = 1, Nucl_N_Aos(i) + i_ao = nucl_aos_transposed(ii,i) + if(ao_l(i_ao)==0)then + do j = 1, ao_prim_num(i_ao) + k+=1 + new_nucl_coord(1:3,k)=nucl_coord_transp(1:3,i) + enddo + else if(ao_l(i_ao)==1)then + ! split the function into 2 s functions + ! one is centered in R_x + d + do j = 1, ao_prim_num(i_ao) + k+=1 + new_nucl_coord(2:3,k)= nucl_coord_transp(2:3,i) + new_nucl_coord(1,k)= nucl_coord_transp(1,i)+ao_extra_center + k+=1 + ! one is centered in R_x - d + new_nucl_coord(2:3,k)= nucl_coord_transp(2:3,i) + new_nucl_coord(1,k)= nucl_coord_transp(1,i)-ao_extra_center + enddo + else + print*,'WARNING ! Lmax value not implemented yet !' + print*,'stopping ...' + stop + endif + enddo + enddo + if(k.ne.n_func_tot)then + print*,'pb !!! k NE n_func_tot !!' + print*,k,n_func_tot + stop + endif + +END_PROVIDER diff --git a/plugins/local/extra_basis_int/NEED b/plugins/local/extra_basis_int/NEED index 7f39af7c..6214e5cd 100644 --- a/plugins/local/extra_basis_int/NEED +++ b/plugins/local/extra_basis_int/NEED @@ -1,2 +1,3 @@ ao_extra_basis ao_one_e_ints +ao_two_e_ints