diff --git a/ocaml/qp_create_ezfio_from_xyz.ml b/ocaml/qp_create_ezfio_from_xyz.ml index 6ec6e0d0..538c5f63 100644 --- a/ocaml/qp_create_ezfio_from_xyz.ml +++ b/ocaml/qp_create_ezfio_from_xyz.ml @@ -14,7 +14,7 @@ let spec = +> flag "m" (optional_with_default 1 int) ~doc:"int Spin multiplicity (2S+1) of the molecule. Default is 1." +> flag "p" no_arg - ~doc:"Using pseudo." + ~doc:"Using pseudopotentials" +> anon ("xyz_file" %: string) ;; diff --git a/scripts/cache_compile.py b/scripts/cache_compile.py new file mode 100755 index 00000000..9898eae5 --- /dev/null +++ b/scripts/cache_compile.py @@ -0,0 +1,62 @@ +#!/usr/bin/env python + +import os +import sys +import shelve +import hashlib +import re + +r = re.compile(ur'-c\s+(\S+\.[fF]90)\s+-o\s+(\S+\.o)') +p = re.compile(ur'-I IRPF90_temp/\S*\s+') +mod = re.compile(ur'module\s+(?P\S+).+end\s?module\s+(?P=mod)?', re.MULTILINE | re.IGNORECASE) + +TMPDIR="/tmp/qp_compiler/" + +def main(): + # Create temp directory + if "qp_compiler" not in os.listdir("/tmp"): + os.mkdir("/tmp/qp_compiler/") + + line = sys.argv[1:] + command = " ".join(line) + command_clean = p.sub('',command) + + try: + match = r.search(command_clean) + input = match.group(1) + output = match.group(2) + except: + os.system(command) + return + m = hashlib.md5() + + # Fread : read input + with open(input,'r') as file: + fread = file.read() + m.update( " ".join( [ command, fread ] )) + + # Md5 Key containing command + content of Fread + key = TMPDIR+m.hexdigest() + try: + # Try to return the content of the .o file + with open(key,'r') as file: + result = file.read() + except IOError: + # Compile the file -> .o + os.system(command) + # Read the .o + with open(output,'r') as file: + result = file.read() + # Copy the .o in database + if not mod.search(fread.replace('\n',' ')): + with open(key,'w') as file: + file.write(result) + else: + print input+' -> module' + else: + # Write the .o file + with open(output,'w') as file: + file.write(result) + +if __name__ == '__main__': + main() diff --git a/scripts/install/install_m4.sh b/scripts/install/install_m4.sh index 6f012d04..8bf514fd 100755 --- a/scripts/install/install_m4.sh +++ b/scripts/install/install_m4.sh @@ -16,6 +16,7 @@ fi cd ${QPACKAGE_ROOT} +rm -f l${QPACKAGE_ROOT}/bin/m4 if [[ -z ${M4} ]] then rm -f -- bin/m4 diff --git a/scripts/pseudo/put_pseudo_in_ezfio.py b/scripts/pseudo/put_pseudo_in_ezfio.py index 98641a18..d5dfb584 100755 --- a/scripts/pseudo/put_pseudo_in_ezfio.py +++ b/scripts/pseudo/put_pseudo_in_ezfio.py @@ -186,7 +186,7 @@ def add_zero(array, size, type): def make_it_square(matrix, dim, type=float): """ - matix the matrix to squate + matix the matrix to square dim array [lmax, kmax] type the null value you want [[[28.59107316], [19.37583724]], [[50.25646328]]] diff --git a/setup_environment.sh b/setup_environment.sh index 878c23e1..d62bd198 100755 --- a/setup_environment.sh +++ b/setup_environment.sh @@ -28,7 +28,7 @@ EOF source quantum_package.rc - +mkdir -p install_logs echo "${BLUE}===== Installing IRPF90 ===== ${BLACK}" ${QPACKAGE_ROOT}/scripts/install/install_irpf90.sh | tee ${QPACKAGE_ROOT}/install_logs/install_irpf90.log if [[ ! -d ${QPACKAGE_ROOT}/irpf90 ]] || [[ ! -x ${QPACKAGE_ROOT}/bin/irpf90 ]] || [[ ! -x ${QPACKAGE_ROOT}/bin/irpman ]] diff --git a/src/AOs/tree_dependancy.png b/src/AOs/tree_dependancy.png index 35cbcbc3..3bf621ee 100644 Binary files a/src/AOs/tree_dependancy.png and b/src/AOs/tree_dependancy.png differ diff --git a/src/Bitmask/bitmasks_module.f90 b/src/Bitmask/bitmasks_module.f90 index c542fac6..c5a9093f 100644 --- a/src/Bitmask/bitmasks_module.f90 +++ b/src/Bitmask/bitmasks_module.f90 @@ -8,4 +8,4 @@ module bitmasks integer, parameter :: d_part2 = 4 integer, parameter :: s_hole = 5 integer, parameter :: s_part = 6 -end module +end module bitmasks diff --git a/src/Bitmask/tree_dependancy.png b/src/Bitmask/tree_dependancy.png index dd1c0f54..c931c9fa 100644 Binary files a/src/Bitmask/tree_dependancy.png and b/src/Bitmask/tree_dependancy.png differ diff --git a/src/CAS_SD/tree_dependancy.png b/src/CAS_SD/tree_dependancy.png index aae463ee..0942fadc 100644 Binary files a/src/CAS_SD/tree_dependancy.png and b/src/CAS_SD/tree_dependancy.png differ diff --git a/src/CID/tree_dependancy.png b/src/CID/tree_dependancy.png index aee0843d..5e421085 100644 Binary files a/src/CID/tree_dependancy.png and b/src/CID/tree_dependancy.png differ diff --git a/src/CID_SC2_selected/tree_dependancy.png b/src/CID_SC2_selected/tree_dependancy.png index 1b24cd17..86e753a4 100644 Binary files a/src/CID_SC2_selected/tree_dependancy.png and b/src/CID_SC2_selected/tree_dependancy.png differ diff --git a/src/CID_selected/tree_dependancy.png b/src/CID_selected/tree_dependancy.png index f6c4cc3e..f8ad68e3 100644 Binary files a/src/CID_selected/tree_dependancy.png and b/src/CID_selected/tree_dependancy.png differ diff --git a/src/CIS/tree_dependancy.png b/src/CIS/tree_dependancy.png index 87b2e1cd..bf88fde6 100644 Binary files a/src/CIS/tree_dependancy.png and b/src/CIS/tree_dependancy.png differ diff --git a/src/CISD/tree_dependancy.png b/src/CISD/tree_dependancy.png index 2be99e68..2d040122 100644 Binary files a/src/CISD/tree_dependancy.png and b/src/CISD/tree_dependancy.png differ diff --git a/src/CISD_SC2_selected/tree_dependancy.png b/src/CISD_SC2_selected/tree_dependancy.png index 8804d41f..ac57a98d 100644 Binary files a/src/CISD_SC2_selected/tree_dependancy.png and b/src/CISD_SC2_selected/tree_dependancy.png differ diff --git a/src/CISD_selected/tree_dependancy.png b/src/CISD_selected/tree_dependancy.png index 282999d6..c33e2e66 100644 Binary files a/src/CISD_selected/tree_dependancy.png and b/src/CISD_selected/tree_dependancy.png differ diff --git a/src/DDCI_selected/tree_dependancy.png b/src/DDCI_selected/tree_dependancy.png index 4ca5ff35..9975ff6c 100644 Binary files a/src/DDCI_selected/tree_dependancy.png and b/src/DDCI_selected/tree_dependancy.png differ diff --git a/src/DensityFit/NEEDED_CHILDREN_MODULES b/src/DensityFit/NEEDED_CHILDREN_MODULES index 13e425cc..d6315de9 100644 --- a/src/DensityFit/NEEDED_CHILDREN_MODULES +++ b/src/DensityFit/NEEDED_CHILDREN_MODULES @@ -1 +1 @@ -AOs +AOs Pseudo diff --git a/src/DensityFit/README.rst b/src/DensityFit/README.rst index 23d6c06f..9d3c4d9b 100644 --- a/src/DensityFit/README.rst +++ b/src/DensityFit/README.rst @@ -10,6 +10,51 @@ Documentation .. Do not edit this section. It was auto-generated from the .. NEEDED_MODULES file. +`aux_basis_coef `_ + Exponents and coefficients of the auxiliary basis + +`aux_basis_coef_transp `_ + Exponents of the auxiliary basis + +`aux_basis_expo `_ + Exponents and coefficients of the auxiliary basis + +`aux_basis_expo_transp `_ + Exponents of the auxiliary basis + +`aux_basis_idx `_ + aux_basis_idx(k) -> i,j + +`aux_basis_nucl `_ + Exponents of the auxiliary basis + +`aux_basis_num `_ + Number of auxiliary basis functions + +`aux_basis_num_8 `_ + Number of auxiliary basis functions + +`aux_basis_num_sqrt `_ + Number of auxiliary basis functions + +`aux_basis_overlap_matrix `_ + Auxiliary basis set + +`aux_basis_power `_ + Exponents of the auxiliary basis + +`aux_basis_prim_num `_ + Exponents of the auxiliary basis + +`aux_basis_prim_num_max `_ + = ao_prim_num_max + +`save_aux_basis `_ + Undocumented + +`aux_basis_four_overlap `_ + \int \chi_i(r) \chi_j(r) \chi_k(r) \chi_l(r) dr + Needed Modules @@ -21,4 +66,5 @@ Needed Modules .. image:: tree_dependancy.png * `AOs `_ +* `Pseudo `_ diff --git a/src/DensityFit/aux_basis.ezfio_config b/src/DensityFit/aux_basis.ezfio_config new file mode 100644 index 00000000..cd56d1c7 --- /dev/null +++ b/src/DensityFit/aux_basis.ezfio_config @@ -0,0 +1,12 @@ +aux_basis + aux_basis_num integer + aux_basis_num_sqrt integer + aux_basis_idx integer (2,aux_basis_aux_basis_num) + aux_basis_prim_num integer (aux_basis_aux_basis_num_sqrt) + aux_basis_nucl integer (aux_basis_aux_basis_num_sqrt) + aux_basis_power integer (aux_basis_aux_basis_num_sqrt,3) + aux_basis_prim_num_max integer = maxval(aux_basis_aux_basis_prim_num) + aux_basis_coef double precision (aux_basis_aux_basis_num_sqrt,aux_basis_aux_basis_prim_num_max) + aux_basis_expo double precision (aux_basis_aux_basis_num_sqrt,aux_basis_aux_basis_prim_num_max) + + diff --git a/src/DensityFit/aux_basis.irp.f b/src/DensityFit/aux_basis.irp.f index 1f28d985..7ca51aa7 100644 --- a/src/DensityFit/aux_basis.irp.f +++ b/src/DensityFit/aux_basis.irp.f @@ -1,11 +1,19 @@ - BEGIN_PROVIDER [ integer, aux_basis_num ] + BEGIN_PROVIDER [ integer, aux_basis_num_sqrt ] +&BEGIN_PROVIDER [ integer, aux_basis_num ] &BEGIN_PROVIDER [ integer, aux_basis_num_8 ] implicit none BEGIN_DOC ! Number of auxiliary basis functions END_DOC integer :: align_double - aux_basis_num = ao_num * (ao_num+1)/2 + + if (do_pseudo) then + aux_basis_num_sqrt = ao_num + ao_pseudo_num +! aux_basis_num_sqrt = ao_num + else + endif + + aux_basis_num = aux_basis_num_sqrt * (aux_basis_num_sqrt+1)/2 aux_basis_num_8 = align_double(aux_basis_num) END_PROVIDER @@ -16,7 +24,7 @@ BEGIN_PROVIDER [ integer, aux_basis_idx, (2,aux_basis_num) ] END_DOC integer :: i,j,k k=0 - do j=1,ao_num + do j=1,aux_basis_num_sqrt do i=1,j k = k+1 aux_basis_idx(1,k) = i @@ -25,6 +33,38 @@ BEGIN_PROVIDER [ integer, aux_basis_idx, (2,aux_basis_num) ] enddo END_PROVIDER + BEGIN_PROVIDER [ double precision, aux_basis_expo_transp, (ao_prim_num_max_align,aux_basis_num_sqrt) ] +&BEGIN_PROVIDER [ double precision, aux_basis_coef_transp, (ao_prim_num_max_align,aux_basis_num_sqrt) ] +&BEGIN_PROVIDER [ integer, aux_basis_prim_num, (aux_basis_num_sqrt) ] +&BEGIN_PROVIDER [ integer, aux_basis_power, (aux_basis_num_sqrt,3) ] +&BEGIN_PROVIDER [ integer, aux_basis_nucl, (aux_basis_num_sqrt) ] + implicit none + BEGIN_DOC + ! Exponents of the auxiliary basis + END_DOC + integer :: i,j + do j=1,ao_num + do i=1,ao_prim_num_max + aux_basis_expo_transp(i,j) = ao_expo_ordered_transp(i,j) + aux_basis_coef_transp(i,j) = ao_coef_normalized_ordered_transp(i,j) + enddo + enddo + do i=1,ao_num + aux_basis_prim_num(i) = ao_prim_num(i) + aux_basis_nucl(i) = ao_nucl(i) + aux_basis_power(i,1:3) = ao_power(i,1:3) + enddo + + do j=1,ao_pseudo_num + aux_basis_expo_transp(1,ao_num+j) = 0.5d0*pseudo_ao_expo(j) + aux_basis_coef_transp(1,ao_num+j) = 1.d0 + aux_basis_power(ao_num+j,1:3) = 0 + aux_basis_prim_num(ao_num+j) = 1 + aux_basis_nucl(ao_num+j) = pseudo_ao_nucl(j) + enddo + +END_PROVIDER + BEGIN_PROVIDER [ double precision, aux_basis_overlap_matrix, (aux_basis_num_8,aux_basis_num) ] implicit none @@ -32,9 +72,9 @@ BEGIN_PROVIDER [ double precision, aux_basis_overlap_matrix, (aux_basis_num_8,au ! Auxiliary basis set END_DOC integer :: m,n,i,j,k,l - double precision :: ao_four_overlap + double precision :: aux_basis_four_overlap - aux_basis_overlap_matrix(1,1) = ao_four_overlap(1,1,1,1) + aux_basis_overlap_matrix(1,1) = aux_basis_four_overlap(1,1,1,1) !$OMP PARALLEL DO PRIVATE(i,j,k,l,m,n) SCHEDULE(GUIDED) do m=1,aux_basis_num i = aux_basis_idx(1,m) @@ -42,7 +82,7 @@ BEGIN_PROVIDER [ double precision, aux_basis_overlap_matrix, (aux_basis_num_8,au do n=1,m k = aux_basis_idx(1,n) l = aux_basis_idx(2,n) - aux_basis_overlap_matrix(m,n) = ao_four_overlap(i,j,k,l) + aux_basis_overlap_matrix(m,n) = aux_basis_four_overlap(i,j,k,l) aux_basis_overlap_matrix(n,m) = aux_basis_overlap_matrix(m,n) enddo enddo @@ -50,5 +90,41 @@ BEGIN_PROVIDER [ double precision, aux_basis_overlap_matrix, (aux_basis_num_8,au END_PROVIDER + BEGIN_PROVIDER [ double precision, aux_basis_expo, (aux_basis_num_sqrt,aux_basis_prim_num_max) ] +&BEGIN_PROVIDER [ double precision, aux_basis_coef, (aux_basis_num_sqrt,aux_basis_prim_num_max) ] + implicit none + BEGIN_DOC + ! Exponents and coefficients of the auxiliary basis + END_DOC + integer :: i,j + aux_basis_expo = 0.d0 + aux_basis_coef = 0.d0 + do j=1,aux_basis_num_sqrt + do i=1,aux_basis_prim_num(j) + aux_basis_expo(j,i) = aux_basis_expo_transp(i,j) + aux_basis_coef(j,i) = aux_basis_coef_transp(i,j) + enddo + enddo + +END_PROVIDER + +BEGIN_PROVIDER [ integer, aux_basis_prim_num_max ] + implicit none + BEGIN_DOC + ! = ao_prim_num_max + END_DOC + aux_basis_prim_num_max = ao_prim_num_max +END_PROVIDER +subroutine save_aux_basis + implicit none + call ezfio_set_aux_basis_aux_basis_num(aux_basis_num) + call ezfio_set_aux_basis_aux_basis_num_sqrt(aux_basis_num_sqrt) + call ezfio_set_aux_basis_aux_basis_idx(aux_basis_idx) + call ezfio_set_aux_basis_aux_basis_prim_num(aux_basis_prim_num) + call ezfio_set_aux_basis_aux_basis_nucl(aux_basis_nucl) + call ezfio_set_aux_basis_aux_basis_power(aux_basis_power) + call ezfio_set_aux_basis_aux_basis_coef(aux_basis_coef) + call ezfio_set_aux_basis_aux_basis_expo(aux_basis_expo) +end diff --git a/src/DensityFit/overlap.irp.f b/src/DensityFit/overlap.irp.f index 0394dbf9..836df432 100644 --- a/src/DensityFit/overlap.irp.f +++ b/src/DensityFit/overlap.irp.f @@ -1,4 +1,4 @@ -double precision function ao_four_overlap(i,j,k,l) +double precision function aux_basis_four_overlap(i,j,k,l) implicit none BEGIN_DOC ! \int \chi_i(r) \chi_j(r) \chi_k(r) \chi_l(r) dr @@ -15,44 +15,44 @@ double precision function ao_four_overlap(i,j,k,l) dim1 = n_pt_max_integrals - num_i = ao_nucl(i) - num_j = ao_nucl(j) - num_k = ao_nucl(k) - num_l = ao_nucl(l) - ao_four_overlap = 0.d0 + num_i = aux_basis_nucl(i) + num_j = aux_basis_nucl(j) + num_k = aux_basis_nucl(k) + num_l = aux_basis_nucl(l) + aux_basis_four_overlap = 0.d0 do p = 1, 3 - I_power(p) = ao_power(i,p) - J_power(p) = ao_power(j,p) - K_power(p) = ao_power(k,p) - L_power(p) = ao_power(l,p) + I_power(p) = aux_basis_power(i,p) + J_power(p) = aux_basis_power(j,p) + K_power(p) = aux_basis_power(k,p) + L_power(p) = aux_basis_power(l,p) I_center(p) = nucl_coord(num_i,p) J_center(p) = nucl_coord(num_j,p) K_center(p) = nucl_coord(num_k,p) L_center(p) = nucl_coord(num_l,p) enddo - do p = 1, ao_prim_num(i) + do p = 1, aux_basis_prim_num(i) double precision :: coef1 - coef1 = ao_coef_normalized_ordered_transp(p,i) - do q = 1, ao_prim_num(j) + coef1 = aux_basis_coef_transp(p,i) + do q = 1, aux_basis_prim_num(j) call give_explicit_poly_and_gaussian(P_new,P_center,pp,fact_p,iorder_p,& - ao_expo_ordered_transp(p,i),ao_expo_ordered_transp(q,j), & + aux_basis_expo_transp(p,i),aux_basis_expo_transp(q,j), & I_power,J_power,I_center,J_center,dim1) double precision :: coef2 - coef2 = coef1*ao_coef_normalized_ordered_transp(q,j)*fact_p - do r = 1, ao_prim_num(k) + coef2 = coef1*aux_basis_coef_transp(q,j)*fact_p + do r = 1, aux_basis_prim_num(k) double precision :: coef3 - coef3 = coef2*ao_coef_normalized_ordered_transp(r,k) - do s = 1, ao_prim_num(l) + coef3 = coef2*aux_basis_coef_transp(r,k) + do s = 1, aux_basis_prim_num(l) double precision :: general_primitive_integral - call give_explicit_poly_and_gaussian(Q_new,Q_center,qq,fact_q,iorder_q,& - ao_expo_ordered_transp(r,k),ao_expo_ordered_transp(s,l), & + call give_explicit_poly_and_gaussian(Q_new,Q_center,qq,fact_q,iorder_q, & + aux_basis_expo_transp(r,k),aux_basis_expo_transp(s,l), & K_power,L_power,K_center,L_center,dim1) double precision :: coef4 - coef4 = coef3*ao_coef_normalized_ordered_transp(s,l)*fact_q + coef4 = coef3*aux_basis_coef_transp(s,l)*fact_q call overlap_gaussian_xyz(P_center,Q_center,pp,qq,iorder_p,iorder_q,overlap_x,overlap_y,overlap_z,overlap,dim1) - ao_four_overlap += coef4 * overlap + aux_basis_four_overlap += coef4 * overlap enddo ! s enddo ! r enddo ! q diff --git a/src/DensityFit/tree_dependancy.png b/src/DensityFit/tree_dependancy.png new file mode 100644 index 00000000..62bad0ff Binary files /dev/null and b/src/DensityFit/tree_dependancy.png differ diff --git a/src/Determinants/README.rst b/src/Determinants/README.rst index a1f0b330..a9455783 100644 --- a/src/Determinants/README.rst +++ b/src/Determinants/README.rst @@ -532,12 +532,6 @@ Documentation `save_casino `_ Undocumented -`save_dets_qmcchem `_ - Undocumented - -`save_for_qmc `_ - Undocumented - `save_natorb `_ Undocumented diff --git a/src/Determinants/save_for_qmcchem.irp.f b/src/Determinants/save_for_qmcchem.irp.f deleted file mode 100644 index 93a9778d..00000000 --- a/src/Determinants/save_for_qmcchem.irp.f +++ /dev/null @@ -1,51 +0,0 @@ -subroutine save_dets_qmcchem - use bitmasks - implicit none - character :: c(mo_tot_num) - integer :: i,k - - integer, allocatable :: occ(:,:,:), occ_tmp(:,:) - !DIR$ ATTRIBUTES ALIGN : $IRP_ALIGN :: occ, occ_tmp - - call ezfio_set_determinants_det_num(N_det) - call ezfio_set_determinants_det_coef(psi_coef_sorted(1,1)) - - allocate (occ(elec_alpha_num,N_det,2)) - ! OMP PARALLEL DEFAULT(NONE) & - ! OMP PRIVATE(occ_tmp,i,k)& - ! OMP SHARED(N_det,psi_det_sorted,elec_alpha_num, & - ! OMP occ,elec_beta_num,N_int) - allocate (occ_tmp(N_int*bit_kind_size,2)) - occ_tmp = 0 - ! OMP DO - do i=1,N_det - call bitstring_to_list(psi_det_sorted(1,1,i), occ_tmp(1,1), elec_alpha_num, N_int ) - call bitstring_to_list(psi_det_sorted(1,2,i), occ_tmp(1,2), elec_beta_num, N_int ) - do k=1,elec_alpha_num - occ(k,i,1) = occ_tmp(k,1) - occ(k,i,2) = occ_tmp(k,2) - enddo - enddo - ! OMP END DO - deallocate(occ_tmp) - ! OMP END PARALLEL - call ezfio_set_determinants_det_occ(occ) - call write_int(output_determinants,N_det,'Determinants saved for QMC') - deallocate(occ) - open(unit=31,file=trim(ezfio_filename)//'/mo_basis/mo_classif') - write(31,'(I1)') 1 - write(31,*) mo_tot_num - do i=1,mo_tot_num - write(31,'(A)') 'a' - enddo - close(31) - call system('gzip -f '//trim(ezfio_filename)//'/mo_basis/mo_classif') - -end - -program save_for_qmc - read_wf = .True. - TOUCH read_wf -! call save_dets_qmcchem - call write_spindeterminants -end diff --git a/src/Determinants/tree_dependancy.png b/src/Determinants/tree_dependancy.png index fb18c1cd..1227a472 100644 Binary files a/src/Determinants/tree_dependancy.png and b/src/Determinants/tree_dependancy.png differ diff --git a/src/FCIdump/tree_dependancy.png b/src/FCIdump/tree_dependancy.png index 80c25a57..5563cfce 100644 Binary files a/src/FCIdump/tree_dependancy.png and b/src/FCIdump/tree_dependancy.png differ diff --git a/src/Full_CI/tree_dependancy.png b/src/Full_CI/tree_dependancy.png index 1297595b..766eec3a 100644 Binary files a/src/Full_CI/tree_dependancy.png and b/src/Full_CI/tree_dependancy.png differ diff --git a/src/Generators_CAS/tree_dependancy.png b/src/Generators_CAS/tree_dependancy.png index 2e2571d9..49fcefcd 100644 Binary files a/src/Generators_CAS/tree_dependancy.png and b/src/Generators_CAS/tree_dependancy.png differ diff --git a/src/Generators_full/tree_dependancy.png b/src/Generators_full/tree_dependancy.png index 659db178..a489df79 100644 Binary files a/src/Generators_full/tree_dependancy.png and b/src/Generators_full/tree_dependancy.png differ diff --git a/src/Hartree_Fock/tree_dependancy.png b/src/Hartree_Fock/tree_dependancy.png index aa20a902..38a98b49 100644 Binary files a/src/Hartree_Fock/tree_dependancy.png and b/src/Hartree_Fock/tree_dependancy.png differ diff --git a/src/Integrals_Bielec/tree_dependancy.png b/src/Integrals_Bielec/tree_dependancy.png index e9430bf0..d7298f31 100644 Binary files a/src/Integrals_Bielec/tree_dependancy.png and b/src/Integrals_Bielec/tree_dependancy.png differ diff --git a/src/Integrals_Monoelec/README.rst b/src/Integrals_Monoelec/README.rst index 775f45d2..b7054f47 100644 --- a/src/Integrals_Monoelec/README.rst +++ b/src/Integrals_Monoelec/README.rst @@ -97,18 +97,12 @@ Documentation `ao_pseudo_integral `_ Pseudo-potential -`ao_pseudo_integral_local `_ +`ao_pseudo_integral_local `_ Local pseudo-potential -`ao_pseudo_integral_non_local `_ +`ao_pseudo_integral_non_local `_ Local pseudo-potential -`pseudo_grid `_ - Grid points for f(|r-r_A|) = \int Y_{lm}^{C} (|r-r_C|, \Omega_C) \chi_i^{A} (r-r_A) d\Omega_C - .br - - `mo_nucl_elec_integral `_ interaction nuclear electron on the MO basis diff --git a/src/Integrals_Monoelec/pot_ao_pseudo_ints.irp.f b/src/Integrals_Monoelec/pot_ao_pseudo_ints.irp.f index 63828c8c..95023177 100644 --- a/src/Integrals_Monoelec/pot_ao_pseudo_ints.irp.f +++ b/src/Integrals_Monoelec/pot_ao_pseudo_ints.irp.f @@ -4,7 +4,9 @@ BEGIN_PROVIDER [ double precision, ao_pseudo_integral, (ao_num_align,ao_num)] ! Pseudo-potential END_DOC if (do_pseudo) then - ao_pseudo_integral = ao_pseudo_integral_local !+ ao_pseudo_integral_non_local + ao_pseudo_integral = ao_pseudo_integral_local + ao_pseudo_integral_non_local +! ao_pseudo_integral = ao_pseudo_integral_local +! ao_pseudo_integral = ao_pseudo_integral_non_local else ao_pseudo_integral = 0.d0 endif @@ -220,67 +222,5 @@ END_PROVIDER END_PROVIDER -BEGIN_PROVIDER [ double precision, pseudo_grid, (pseudo_grid_size,ao_num,-pseudo_lmax:pseudo_lmax,0:pseudo_lmax,nucl_num) ] - implicit none - BEGIN_DOC -! Grid points for f(|r-r_A|) = \int Y_{lm}^{C} (|r-r_C|, \Omega_C) \chi_i^{A} (r-r_A) d\Omega_C -! -! - END_DOC - ! l,m : Y(l,m) parameters - ! c(3) : pseudopotential center - ! a(3) : Atomic Orbital center - ! n_a(3) : Powers of x,y,z in the Atomic Orbital - ! g_a : Atomic Orbital exponent - ! r : Distance between the Atomic Orbital center and the considered point - double precision, external :: ylm_orb - integer :: n_a(3) - double precision :: a(3), c(3), g_a - integer :: i,j,k,l,m,n,p - double precision :: r(pseudo_grid_size), dr, Ulc - double precision, parameter :: rmax= 10.d0 - - dr = rmax/dble(pseudo_grid_size) - r(1) = 0.d0 - do j=2,pseudo_grid_size - r(j) = r(j-1) + dr - enddo - - pseudo_grid = 0.d0 - do k=1,nucl_num - c(1:3) = nucl_coord(k,1:3) - do l=0,pseudo_lmax - do i=1,ao_num - a(1:3) = nucl_coord(ao_nucl(i),1:3) - n_a(1:3) = ao_power(i,1:3) - do j=1,pseudo_grid_size - do p=1,ao_prim_num(i) - g_a = ao_expo_ordered_transp(p,i) - do m=-l,l - double precision :: y - ! y = ylm_orb(l,m,c,a,n_a,g_a,r(j)) - ! if (y /= 0.d0) then - ! print *, 'y = ', y - ! print *, 'l = ', l - ! print *, 'm = ', m - ! print *, 'c = ', c - ! print *, 'a = ', a - ! print *, 'n = ', n_a - ! print *, 'g = ', g_a - ! print *, 'r = ', r(j) - ! print *, '' - ! endif - y = ylm_orb(l,m,c,a,n_a,g_a,r(j)) - pseudo_grid(j,i,m,l,k) = pseudo_grid(j,i,m,l,k) + & - ao_coef_normalized_ordered_transp(p,i)*y - enddo - enddo - enddo - enddo - enddo - enddo - -END_PROVIDER diff --git a/src/Integrals_Monoelec/pot_mo_pseudo_ints.irp.f b/src/Integrals_Monoelec/pot_mo_pseudo_ints.irp.f index d34a05ff..6c412e4b 100644 --- a/src/Integrals_Monoelec/pot_mo_pseudo_ints.irp.f +++ b/src/Integrals_Monoelec/pot_mo_pseudo_ints.irp.f @@ -28,7 +28,6 @@ BEGIN_PROVIDER [double precision, mo_pseudo_integral, (mo_tot_num_align,mo_tot_n enddo enddo !$OMP END PARALLEL DO - call ezfio_set_pseudo_pseudo_matrix(mo_pseudo_integral(1:mo_tot_num,1:mo_tot_num)) END_PROVIDER diff --git a/src/Integrals_Monoelec/pseudopot.f90 b/src/Integrals_Monoelec/pseudopot.f90 index d8ff7893..b605a45e 100644 --- a/src/Integrals_Monoelec/pseudopot.f90 +++ b/src/Integrals_Monoelec/pseudopot.f90 @@ -1485,7 +1485,6 @@ end a=bessel_mod_exp(n,x) return endif - print *, n,x if(n.eq.0)a=dsinh(x)/x if(n.eq.1)a=(x*dcosh(x)-dsinh(x))/x**2 if(n.ge.2)a=bessel_mod_recur(n-2,x)-(2*n-1)/x*bessel_mod_recur(n-1,x) diff --git a/src/Integrals_Monoelec/tree_dependancy.png b/src/Integrals_Monoelec/tree_dependancy.png index 490d33b5..b12bf22c 100644 Binary files a/src/Integrals_Monoelec/tree_dependancy.png and b/src/Integrals_Monoelec/tree_dependancy.png differ diff --git a/src/MOGuess/tree_dependancy.png b/src/MOGuess/tree_dependancy.png index 14da5742..fa7ab7f1 100644 Binary files a/src/MOGuess/tree_dependancy.png and b/src/MOGuess/tree_dependancy.png differ diff --git a/src/MOs/tree_dependancy.png b/src/MOs/tree_dependancy.png index 7250aa3f..07d88313 100644 Binary files a/src/MOs/tree_dependancy.png and b/src/MOs/tree_dependancy.png differ diff --git a/src/MP2/tree_dependancy.png b/src/MP2/tree_dependancy.png index f0573cdf..f4bd74fa 100644 Binary files a/src/MP2/tree_dependancy.png and b/src/MP2/tree_dependancy.png differ diff --git a/src/MRCC/README.rst b/src/MRCC/README.rst index 20ff8014..ece272a1 100644 --- a/src/MRCC/README.rst +++ b/src/MRCC/README.rst @@ -53,7 +53,7 @@ Documentation `ci_electronic_energy_dressed `_ Eigenvectors/values of the CI matrix -`ci_energy_dressed `_ +`ci_energy_dressed `_ N_states lowest eigenvalues of the dressed CI matrix `delta_ij `_ @@ -62,7 +62,7 @@ Documentation `delta_ij_non_cas `_ Dressing matrix in SD basis -`diagonalize_ci_dressed `_ +`diagonalize_ci_dressed `_ Replace the coefficients of the CI states by the coefficients of the eigenstates of the CI matrix diff --git a/src/MRCC/mrcc_utils.irp.f b/src/MRCC/mrcc_utils.irp.f index 611a70be..872588cf 100644 --- a/src/MRCC/mrcc_utils.irp.f +++ b/src/MRCC/mrcc_utils.irp.f @@ -92,9 +92,10 @@ END_PROVIDER if (diag_algorithm == "Davidson") then - stop 'use Lapack' -! call davidson_diag(psi_det,CI_eigenvectors_dressed,CI_electronic_energy_dressed, & -! size(CI_eigenvectors_dressed,1),N_det,N_states_diag,N_int,output_determinants) + integer :: istate + istate = 1 + call davidson_diag_mrcc(psi_det,CI_eigenvectors_dressed,CI_electronic_energy_dressed, & + size(CI_eigenvectors_dressed,1),N_det,N_states_diag,N_int,output_determinants,istate) else if (diag_algorithm == "Lapack") then diff --git a/src/MRCC/tree_dependancy.png b/src/MRCC/tree_dependancy.png index 9b5920af..ccd5da28 100644 Binary files a/src/MRCC/tree_dependancy.png and b/src/MRCC/tree_dependancy.png differ diff --git a/src/Makefile.config.example b/src/Makefile.config.example index fe6cee8b..410bd0e2 100644 --- a/src/Makefile.config.example +++ b/src/Makefile.config.example @@ -4,6 +4,7 @@ DEBUG = 0 IRPF90_FLAGS+= --align=32 FC = ifort -g +#FC = cache_compile.py ifort -g # Accelerates compilation FCFLAGS= FCFLAGS+= -xHost #FCFLAGS+= -xAVX diff --git a/src/Molden/tree_dependancy.png b/src/Molden/tree_dependancy.png index 887372a5..23007e0f 100644 Binary files a/src/Molden/tree_dependancy.png and b/src/Molden/tree_dependancy.png differ diff --git a/src/NEEDED_MODULES b/src/NEEDED_MODULES index b8ea4e38..6e14bf1d 100644 --- a/src/NEEDED_MODULES +++ b/src/NEEDED_MODULES @@ -1 +1 @@ -AOs Integrals_Bielec Bitmask CAS_SD CID CID_SC2_selected CID_selected CIS CISD CISD_SC2_selected CISD_selected DDCI_selected Determinants Electrons Ezfio_files FCIdump Full_CI Generators_CAS Generators_full Hartree_Fock MOGuess Molden Integrals_Monoelec MOs MP2 MRCC Nuclei Pseudo Selectors_full Utils DensityFit +AOs Integrals_Bielec Bitmask CAS_SD CID CID_SC2_selected CID_selected CIS CISD CISD_SC2_selected CISD_selected DDCI_selected Determinants Electrons Ezfio_files FCIdump Full_CI Generators_CAS Generators_full Hartree_Fock MOGuess Molden Integrals_Monoelec MOs MP2 MRCC Nuclei Pseudo Selectors_full Utils DensityFit QmcChem diff --git a/src/Nuclei/NEEDED_CHILDREN_MODULES b/src/Nuclei/NEEDED_CHILDREN_MODULES index d9eaf57c..dcdb5f86 100644 --- a/src/Nuclei/NEEDED_CHILDREN_MODULES +++ b/src/Nuclei/NEEDED_CHILDREN_MODULES @@ -1 +1 @@ -Ezfio_files +Ezfio_files Utils diff --git a/src/Nuclei/README.rst b/src/Nuclei/README.rst index 7bdf017f..2e99e8fb 100644 --- a/src/Nuclei/README.rst +++ b/src/Nuclei/README.rst @@ -15,6 +15,7 @@ Needed Modules .. image:: tree_dependancy.png * `Ezfio_files `_ +* `Utils `_ Documentation ============= diff --git a/src/Nuclei/tree_dependancy.png b/src/Nuclei/tree_dependancy.png index dfd6f57a..72cfaeee 100644 Binary files a/src/Nuclei/tree_dependancy.png and b/src/Nuclei/tree_dependancy.png differ diff --git a/src/Pseudo/EZFIO.cfg b/src/Pseudo/EZFIO.cfg index 5a401d18..bf3909d1 100644 --- a/src/Pseudo/EZFIO.cfg +++ b/src/Pseudo/EZFIO.cfg @@ -57,10 +57,16 @@ default: False [pseudo_grid_size] type: integer -doc: Size of the QMC grid +doc: Nb of points of the QMC grid interface: input default: 100 +[pseudo_grid_rmax] +type: double precision +doc: R_maxof the QMC grid +interface: input +default: 4.0 + [pseudo_grid] type: double precision doc: QMC grid @@ -71,6 +77,6 @@ size: (pseudo.pseudo_grid_size,ao_basis.ao_num,-pseudo.pseudo_lmax:pseudo.pseudo type: double precision doc: QMC grid interface: output -size: (mo_basis.mo_tot_num,mo_basis.mo_tot_num) +size: (aux_basis.aux_basis_num_sqrt,aux_basis.aux_basis_num_sqrt) diff --git a/src/Pseudo/NEEDED_CHILDREN_MODULES b/src/Pseudo/NEEDED_CHILDREN_MODULES index d9eaf57c..d47a550b 100644 --- a/src/Pseudo/NEEDED_CHILDREN_MODULES +++ b/src/Pseudo/NEEDED_CHILDREN_MODULES @@ -1 +1 @@ -Ezfio_files +Ezfio_files Nuclei diff --git a/src/Pseudo/README.rst b/src/Pseudo/README.rst index 284499d7..f6e7fc6e 100644 --- a/src/Pseudo/README.rst +++ b/src/Pseudo/README.rst @@ -19,4 +19,5 @@ Needed Modules .. image:: tree_dependancy.png * `Ezfio_files `_ +* `Nuclei `_ diff --git a/src/Pseudo/tree_dependancy.png b/src/Pseudo/tree_dependancy.png index 389cb115..98e8be37 100644 Binary files a/src/Pseudo/tree_dependancy.png and b/src/Pseudo/tree_dependancy.png differ diff --git a/src/QmcChem/ASSUMPTIONS.rst b/src/QmcChem/ASSUMPTIONS.rst new file mode 100644 index 00000000..e69de29b diff --git a/src/QmcChem/Makefile b/src/QmcChem/Makefile new file mode 100644 index 00000000..06dc50ff --- /dev/null +++ b/src/QmcChem/Makefile @@ -0,0 +1,6 @@ +# Define here all new external source files and objects.Don't forget to prefix the +# object files with IRPF90_temp/ +SRC= +OBJ= + +include $(QPACKAGE_ROOT)/src/Makefile.common diff --git a/src/QmcChem/NEEDED_CHILDREN_MODULES b/src/QmcChem/NEEDED_CHILDREN_MODULES new file mode 100644 index 00000000..f7ed9913 --- /dev/null +++ b/src/QmcChem/NEEDED_CHILDREN_MODULES @@ -0,0 +1 @@ +Determinants DensityFit diff --git a/src/QmcChem/README.rst b/src/QmcChem/README.rst new file mode 100644 index 00000000..c4b467d7 --- /dev/null +++ b/src/QmcChem/README.rst @@ -0,0 +1,23 @@ +============== +QmcChem Module +============== + +Documentation +============= + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + + + +Needed Modules +============== + +.. Do not edit this section. It was auto-generated from the +.. NEEDED_MODULES file. + +.. image:: tree_dependancy.png + +* `Determinants `_ +* `DensityFit `_ + diff --git a/src/QmcChem/pot_ao_pseudo_ints.irp.f b/src/QmcChem/pot_ao_pseudo_ints.irp.f new file mode 100644 index 00000000..24d75504 --- /dev/null +++ b/src/QmcChem/pot_ao_pseudo_ints.irp.f @@ -0,0 +1,275 @@ +BEGIN_PROVIDER [ double precision, aux_pseudo_integral, (aux_basis_num_sqrt,aux_basis_num_sqrt)] + implicit none + BEGIN_DOC +! Pseudo-potential + END_DOC + if (do_pseudo) then +! aux_pseudo_integral = aux_pseudo_integral_local + aux_pseudo_integral_non_local +! aux_pseudo_integral = aux_pseudo_integral_local + aux_pseudo_integral = aux_pseudo_integral_non_local + else + aux_pseudo_integral = 0.d0 + endif +END_PROVIDER + + BEGIN_PROVIDER [ double precision, aux_pseudo_integral_local, (aux_basis_num_sqrt,aux_basis_num_sqrt)] + implicit none + BEGIN_DOC +! Local pseudo-potential + END_DOC + double precision :: alpha, beta, gama, delta + integer :: num_A,num_B + double precision :: A_center(3),B_center(3),C_center(3) + integer :: power_A(3),power_B(3) + integer :: i,j,k,l,n_pt_in,m + double precision :: Vloc, Vpseudo + + double precision :: cpu_1, cpu_2, wall_1, wall_2, wall_0 + integer :: thread_num + + aux_pseudo_integral_local = 0.d0 + + !! Dump array + integer, allocatable :: n_k_dump(:) + double precision, allocatable :: v_k_dump(:), dz_k_dump(:) + + allocate(n_k_dump(1:pseudo_klocmax), v_k_dump(1:pseudo_klocmax), dz_k_dump(1:pseudo_klocmax)) + + + ! _ + ! / _. | _ | + ! \_ (_| | (_ |_| | + ! + + print*, 'Providing the nuclear electron pseudo integrals ' + + call wall_time(wall_1) + call cpu_time(cpu_1) + + !$OMP PARALLEL & + !$OMP DEFAULT (NONE) & + !$OMP PRIVATE (i,j,k,l,m,alpha,beta,A_center,B_center,C_center,power_A,power_B, & + !$OMP num_A,num_B,Z,c,n_pt_in, & + !$OMP v_k_dump,n_k_dump, dz_k_dump, & + !$OMP wall_0,wall_2,thread_num) & + !$OMP SHARED (aux_basis_num_sqrt,aux_basis_prim_num,aux_basis_expo_transp,aux_basis_power,aux_basis_nucl,nucl_coord,aux_basis_coef_transp, & + !$OMP aux_pseudo_integral_local,nucl_num,nucl_charge, & + !$OMP pseudo_klocmax,pseudo_lmax,pseudo_kmax,pseudo_v_k,pseudo_n_k, pseudo_dz_k, & + !$OMP wall_1) + + !$OMP DO SCHEDULE (guided) + + do j = 1, aux_basis_num_sqrt + + num_A = aux_basis_nucl(j) + power_A(1:3)= aux_basis_power(j,1:3) + A_center(1:3) = nucl_coord(num_A,1:3) + + do i = 1, aux_basis_num_sqrt + + num_B = aux_basis_nucl(i) + power_B(1:3)= aux_basis_power(i,1:3) + B_center(1:3) = nucl_coord(num_B,1:3) + + do l=1,aux_basis_prim_num(j) + alpha = aux_basis_expo_transp(l,j) + + do m=1,aux_basis_prim_num(i) + beta = aux_basis_expo_transp(m,i) + double precision :: c + c = 0.d0 + + do k = 1, nucl_num + double precision :: Z + Z = nucl_charge(k) + + C_center(1:3) = nucl_coord(k,1:3) + + v_k_dump = pseudo_v_k(k,1:pseudo_klocmax) + n_k_dump = pseudo_n_k(k,1:pseudo_klocmax) + dz_k_dump = pseudo_dz_k(k,1:pseudo_klocmax) + + c = c + Vloc(pseudo_klocmax, v_k_dump,n_k_dump, dz_k_dump, & + A_center,power_A,alpha,B_center,power_B,beta,C_center) + + enddo + aux_pseudo_integral_local(i,j) = aux_pseudo_integral_local(i,j) + & + aux_basis_coef_transp(l,j)*aux_basis_coef_transp(m,i)*c + enddo + enddo + enddo + + call wall_time(wall_2) + if (thread_num == 0) then + if (wall_2 - wall_0 > 1.d0) then + wall_0 = wall_2 + print*, 100.*float(j)/float(aux_basis_num_sqrt), '% in ', & + wall_2-wall_1, 's' + endif + endif + enddo + + !$OMP END DO + !$OMP END PARALLEL + + + deallocate(n_k_dump,v_k_dump, dz_k_dump) + + END_PROVIDER + + + BEGIN_PROVIDER [ double precision, aux_pseudo_integral_non_local, (aux_basis_num_sqrt,aux_basis_num_sqrt)] + implicit none + BEGIN_DOC +! Local pseudo-potential + END_DOC + double precision :: alpha, beta, gama, delta + integer :: num_A,num_B + double precision :: A_center(3),B_center(3),C_center(3) + integer :: power_A(3),power_B(3) + integer :: i,j,k,l,n_pt_in,m + double precision :: Vloc, Vpseudo + + double precision :: cpu_1, cpu_2, wall_1, wall_2, wall_0 + integer :: thread_num + + aux_pseudo_integral_non_local = 0.d0 + + !! Dump array + integer, allocatable :: n_kl_dump(:,:) + double precision, allocatable :: v_kl_dump(:,:), dz_kl_dump(:,:) + + allocate(n_kl_dump(pseudo_kmax,0:pseudo_lmax), v_kl_dump(pseudo_kmax,0:pseudo_lmax), dz_kl_dump(pseudo_kmax,0:pseudo_lmax)) + + ! _ + ! / _. | _ | + ! \_ (_| | (_ |_| | + ! + + print*, 'Providing the nuclear electron pseudo integrals ' + + call wall_time(wall_1) + call cpu_time(cpu_1) + + !$OMP PARALLEL & + !$OMP DEFAULT (NONE) & + !$OMP PRIVATE (i,j,k,l,m,alpha,beta,A_center,B_center,C_center,power_A,power_B, & + !$OMP num_A,num_B,Z,c,n_pt_in, & + !$OMP n_kl_dump, v_kl_dump, dz_kl_dump, & + !$OMP wall_0,wall_2,thread_num) & + !$OMP SHARED (aux_basis_num_sqrt,aux_basis_prim_num,aux_basis_expo_transp,aux_basis_power,aux_basis_nucl,nucl_coord,aux_basis_coef_transp, & + !$OMP aux_pseudo_integral_non_local,nucl_num,nucl_charge, & + !$OMP pseudo_klocmax,pseudo_lmax,pseudo_kmax,pseudo_n_kl, pseudo_v_kl, pseudo_dz_kl, & + !$OMP wall_1) + + !$OMP DO SCHEDULE (guided) + + do j = 1, aux_basis_num_sqrt + + num_A = aux_basis_nucl(j) + power_A(1:3)= aux_basis_power(j,1:3) + A_center(1:3) = nucl_coord(num_A,1:3) + + do i = 1, aux_basis_num_sqrt + + num_B = aux_basis_nucl(i) + power_B(1:3)= aux_basis_power(i,1:3) + B_center(1:3) = nucl_coord(num_B,1:3) + + do l=1,aux_basis_prim_num(j) + alpha = aux_basis_expo_transp(l,j) + + do m=1,aux_basis_prim_num(i) + beta = aux_basis_expo_transp(m,i) + double precision :: c + c = 0.d0 + + do k = 1, nucl_num + double precision :: Z + Z = nucl_charge(k) + + C_center(1:3) = nucl_coord(k,1:3) + + n_kl_dump = pseudo_n_kl(k,1:pseudo_kmax,0:pseudo_lmax) + v_kl_dump = pseudo_v_kl(k,1:pseudo_kmax,0:pseudo_lmax) + dz_kl_dump = pseudo_dz_kl(k,1:pseudo_kmax,0:pseudo_lmax) + + c = c + Vpseudo(pseudo_lmax,pseudo_kmax,v_kl_dump,n_kl_dump,dz_kl_dump,A_center,power_A,alpha,B_center,power_B,beta,C_center) + + enddo + aux_pseudo_integral_non_local(i,j) = aux_pseudo_integral_non_local(i,j) + & + aux_basis_coef_transp(l,j)*aux_basis_coef_transp(m,i)*c + enddo + enddo + enddo + + call wall_time(wall_2) + if (thread_num == 0) then + if (wall_2 - wall_0 > 1.d0) then + wall_0 = wall_2 + print*, 100.*float(j)/float(aux_basis_num_sqrt), '% in ', & + wall_2-wall_1, 's' + endif + endif + enddo + + !$OMP END DO + !$OMP END PARALLEL + + + deallocate(n_kl_dump,v_kl_dump, dz_kl_dump) + + + END_PROVIDER + +BEGIN_PROVIDER [ double precision, pseudo_grid, (pseudo_grid_size,ao_num,-pseudo_lmax:pseudo_lmax,0:pseudo_lmax,nucl_num) ] + implicit none + BEGIN_DOC +! Grid points for f(|r-r_A|) = \int Y_{lm}^{C} (|r-r_C|, \Omega_C) \chi_i^{A} (r-r_A) d\Omega_C +! +! + END_DOC + ! l,m : Y(l,m) parameters + ! c(3) : pseudopotential center + ! a(3) : Atomic Orbital center + ! n_a(3) : Powers of x,y,z in the Atomic Orbital + ! g_a : Atomic Orbital exponent + ! r : Distance between the Atomic Orbital center and the considered point + double precision, external :: ylm_orb + integer :: n_a(3) + double precision :: a(3), c(3), g_a + integer :: i,j,k,l,m,n,p + double precision :: r(pseudo_grid_size), dr, Ulc + double precision :: y + + dr = pseudo_grid_rmax/dble(pseudo_grid_size) + r(1) = 0.d0 + do j=2,pseudo_grid_size + r(j) = r(j-1) + dr + enddo + + pseudo_grid = 0.d0 + do k=1,nucl_num + c(1:3) = nucl_coord(k,1:3) + do l=0,pseudo_lmax + do i=1,ao_num + a(1:3) = nucl_coord(ao_nucl(i),1:3) + n_a(1:3) = ao_power(i,1:3) + do j=1,pseudo_grid_size + do p=1,ao_prim_num(i) + g_a = ao_expo_ordered_transp(p,i) + do m=-l,l + y = ylm_orb(l,m,c,a,n_a,g_a,r(j)) + pseudo_grid(j,i,m,l,k) = pseudo_grid(j,i,m,l,k) + & + ao_coef_normalized_ordered_transp(p,i)*y + enddo + enddo + enddo + enddo + enddo + enddo + +END_PROVIDER + + diff --git a/src/QmcChem/pseudo.irp.f b/src/QmcChem/pseudo.irp.f new file mode 100644 index 00000000..b0cd57a5 --- /dev/null +++ b/src/QmcChem/pseudo.irp.f @@ -0,0 +1,170 @@ +subroutine write_pseudopotential + implicit none + BEGIN_DOC +! Write the pseudo_potential into the EZFIO file + END_DOC +! call ezfio_set_pseudo_pseudo_matrix(pseudo_matrix) + call ezfio_set_pseudo_pseudo_grid(pseudo_grid) +end + + +BEGIN_PROVIDER [ double precision, pseudo_matrix, (aux_basis_num_sqrt,aux_basis_num_sqrt) ] + implicit none + BEGIN_DOC + ! Pseudo-potential expressed in the basis of ao products + END_DOC + + integer :: i,j,k,l + integer :: info, m,n, lwork, lda, ldu, ldvt + integer, allocatable :: iwork(:) + character :: jobz + double precision, allocatable :: a(:,:),work(:) + + double precision,allocatable :: U(:,:) + double precision,allocatable :: Vt(:,:) + double precision,allocatable :: S(:), B(:) + + + + jobz = 'A' + m = aux_basis_num + n = aux_basis_num + lda = size(aux_basis_overlap_matrix,1) + ldu = lda + ldvt = lda + lwork = -1 + +! allocate (A(lda,n), U(ldu,n), Vt(ldvt,n), S(n), work(1), b(n), iwork(8*n)) + allocate (A(lda,n), U(ldu,n), Vt(ldvt,n), S(n), work(1), b(n),iwork(1)) + + work(1) = 1 + do i=1,n + do j=1,n + A(i,j) = aux_basis_overlap_matrix(i,j) + enddo + enddo + +! call dgesdd(jobz, m, n, A, lda, s, u, ldu, vt, ldvt, work, lwork, iwork, info) + call dgesvd(jobz, jobz, m, n, a, lda, s, u, ldu, vt, ldvt, work, lwork, info) + lwork = int(work(1)) + deallocate(work) + + print *, 'Fitting pseudo-potentials' + + allocate(work(lwork)) +! call dgesdd(jobz, m, n, A, lda, s, u, ldu, vt, ldvt, work, lwork, iwork, info) + call dgesvd(jobz, jobz, m, n, a, lda, s, u, ldu, vt, ldvt, work, lwork, info) + deallocate(work) + + do i=1,n + print *, i, s(i) + enddo + + do k=1,n + if (s(k) < 1.d-1) then + s(k) = 0.d0 + else + s(k) = 1.d0/s(k) + endif + do m=1,n + Vt(m,k) = S(m) * Vt(m,k) + enddo + enddo + call dgemm('N','N',n,n,n,1.d0,U,lda,Vt,ldvt,0.d0,A,lda) +! do k=1,n +! do l=1,n +! A(k,l) = 0.d0 +! do m=1,n +! A(k,l) = A(k,l) + U(k,m) * Vt(m,l) +! enddo +! enddo + + do k=1,n + i = aux_basis_idx(1,k) + j = aux_basis_idx(2,k) + b(k) = aux_pseudo_integral(i,j) + enddo + + do k=1,n + S(k) = 0.d0 + enddo + do l=1,n + do k=1,n + S(k) = S(k) + A(k,l) * b(l) + enddo + enddo + + do k=1,aux_basis_num + i = aux_basis_idx(1,k) + j = aux_basis_idx(2,k) + pseudo_matrix(i,j) = S(k) + pseudo_matrix(j,i) = S(k) + enddo + deallocate(a,b,s,iwork,u,vt) + +print *, 'Done' + if (info /= 0) then + print *, info + stop 'pseudo fit failed' + endif +END_PROVIDER + + + + + +!BEGIN_PROVIDER [ double precision, pseudo_matrix, (ao_num,ao_num) ] +! implicit none +! BEGIN_DOC +! ! Pseudo-potential expressed in the basis of ao products +! END_DOC +! +! integer :: i,j,k +! integer :: info, n, lwork, lda, ldb, nrhs +! character :: uplo +! integer, allocatable :: ipiv(:) +! double precision, allocatable :: a(:,:),work(:), b(:) +! +! uplo = 'L' +! n = aux_basis_num +! nrhs = 1 +! lda = size(aux_basis_overlap_matrix,1) +! ldb = n +! lwork = -1 +! +! print *, 'Fitting pseudo-potentials' +! allocate(work(1),a(lda,n),ipiv(n),b(n)) +! work(1) = 1 +! do i=1,n +! do j=1,n +! a(i,j) = aux_basis_overlap_matrix(i,j) +! enddo +! enddo +! +! do k=1,n +! i = aux_basis_idx(1,k) +! j = aux_basis_idx(2,k) +! b(k) = ao_pseudo_integral(i,j) +! enddo +! call dsysv( uplo, n, nrhs, a, lda, ipiv, b, ldb, work, lwork, info ) +! lwork = int(work(1)) +! deallocate(work) +! +! allocate(work(lwork)) +! call dsysv( uplo, n, nrhs, a, lda, ipiv, b, ldb, work, lwork, info ) +! deallocate(work,ipiv) +! do k=1,aux_basis_num +! i = aux_basis_idx(1,k) +! j = aux_basis_idx(2,k) +! pseudo_matrix(i,j) = b(k) +! pseudo_matrix(j,i) = b(k) +! enddo +! deallocate(a,b) +! +!print *, 'Done' +! if (info /= 0) then +! print *, info +! stop 'pseudo fit failed' +! endif +!END_PROVIDER + diff --git a/src/QmcChem/save_for_qmcchem.irp.f b/src/QmcChem/save_for_qmcchem.irp.f new file mode 100644 index 00000000..cfc18338 --- /dev/null +++ b/src/QmcChem/save_for_qmcchem.irp.f @@ -0,0 +1,8 @@ +program save_for_qmc + read_wf = .True. + TOUCH read_wf +! call save_dets_qmcchem + call write_spindeterminants +! call write_pseudopotential +! call save_aux_basis +end diff --git a/src/Selectors_full/tree_dependancy.png b/src/Selectors_full/tree_dependancy.png index f88049bd..08ba5eda 100644 Binary files a/src/Selectors_full/tree_dependancy.png and b/src/Selectors_full/tree_dependancy.png differ