From 2aa2c38b8edae700bee95326ec0269c83cf174b2 Mon Sep 17 00:00:00 2001 From: NastaMauger Date: Thu, 27 Mar 2025 12:32:44 -0400 Subject: [PATCH 1/8] Add README_branch --- README_branch.md | 14 ++++++++++++++ 1 file changed, 14 insertions(+) create mode 100644 README_branch.md diff --git a/README_branch.md b/README_branch.md new file mode 100644 index 0000000000..8df7a5fa4d --- /dev/null +++ b/README_branch.md @@ -0,0 +1,14 @@ +# Repository Name + +## Branches Description + +### `develop` branch +This branch contains the latest development work and merges from upstream. + +### `Batched_Converter` branch +This branch includes the work related to the Batched Converter feature. It is based on the forked repository from `anbenali`. + +### `LCAO_NumCenters_Opt` branch +This branch contains optimizations for the `LCAO_NumCenters_Opt` feature. It is based on the forked repository from `anbenali`. + + From d11d00932859dfacf3d53a21aae8e44f6f81f00d Mon Sep 17 00:00:00 2001 From: NastaMauger Date: Mon, 30 Jun 2025 21:49:32 -0400 Subject: [PATCH 2/8] Update PyscfToQmcpack.py: fix multidetector handling --- src/QMCTools/PyscfToQmcpack.py | 14 ++++++++------ 1 file changed, 8 insertions(+), 6 deletions(-) diff --git a/src/QMCTools/PyscfToQmcpack.py b/src/QMCTools/PyscfToQmcpack.py index 45f07413bd..369b3a1d90 100644 --- a/src/QMCTools/PyscfToQmcpack.py +++ b/src/QMCTools/PyscfToQmcpack.py @@ -570,26 +570,28 @@ def get_mo(mo_coeff, cart): #Unsorted Mo_coeffs for Multideterminants order matching QP eigenset_unsorted=GroupDet.create_dataset("eigenset_unsorted_0",(NbMO,NbAO),dtype="f8",data=mo_coeff_unsorted) eigenvalue_unsorted=GroupDet.create_dataset("eigenval_unsorted_0",(1,NbMO),dtype="f8",data=E_g_unsorted) - GroupParameter.create_dataset("numMO",(1,),dtype="i4",data=NbMO) GroupParameter.create_dataset("numAO",(1,),dtype="i4",data=NbAO) - - make_multidet(cell,mf,title, H5_qmcpack) + + is_multidet = isinstance(mf, (mcscf.casci.CASCI, mcscf.mc1step.CASSCF)) + if is_multidet: + make_multidet(cell, mf, title, H5_qmcpack) + + # Close the file before exiting H5_qmcpack.close() print ('Wavefunction successfully saved to QMCPACK HDF5 Format') print ('Use: "convert4qmc -orbitals {}.h5" to generate QMCPACK input files'.format(title)) - # Close the file before exiting + if is_multidet: + print(f'Multideterminant wavefunction saved to {title}_multidet.h5') def make_multidet(cell, mf, title, h5_handle): from pyscf import mcscf import numpy import h5py, re, sys - possible_determinant_source =[mcscf.casci.CASCI, mcscf.mc1step.CASSCF] - if type(mf) in possible_determinant_source: a = mf.fcisolver.large_ci(mf.ci, mf.ncas, mf.nelecas, tol=0.0, return_strs=True) dets_a = [] dets_b = [] From 1f6fc57a1da77f04b8c04de005e80ae2a7923652 Mon Sep 17 00:00:00 2001 From: NastaMauger Date: Mon, 30 Jun 2025 22:14:32 -0400 Subject: [PATCH 3/8] Clang format check --- src/QMCTools/PyscfToQmcpack.py | 191 +++++++++++++++++---------------- 1 file changed, 98 insertions(+), 93 deletions(-) diff --git a/src/QMCTools/PyscfToQmcpack.py b/src/QMCTools/PyscfToQmcpack.py index 369b3a1d90..0e0d5a2948 100644 --- a/src/QMCTools/PyscfToQmcpack.py +++ b/src/QMCTools/PyscfToQmcpack.py @@ -56,7 +56,7 @@ def savetoqmcpack(cell,mf,title="Default",kpts=[],kmesh=[],sp_twist=[],weight=1. def get_supercell(cell,kmesh=[]): latt_vec = cell.lattice_vectors() if len(kmesh)==0: - # Guess kmesh +#Guess kmesh scaled_k = cell.get_scaled_kpts(kpts).round(8) kmesh = (len(numpy.unique(scaled_k[:,0])), len(numpy.unique(scaled_k[:,1])), @@ -69,7 +69,7 @@ def get_supercell(cell,kmesh=[]): R_vec_abs = numpy.einsum('nu, uv -> nv', R_vec_rel, latt_vec) - # R_rel_mesh has to be construct exactly same to the Ts in super_cell function +#R_rel_mesh has to be construct exactly same to the Ts in super_cell function scell = tools.super_cell(cell, kmesh) return scell, kmesh @@ -81,7 +81,7 @@ def get_phase(cell, kpts, kmesh=[]): latt_vec = cell.lattice_vectors() if kmesh is None: - # Guess kmesh +#Guess kmesh scaled_k = cell.get_scaled_kpts(kpts).round(8) kmesh = (len(numpy.unique(scaled_k[:,0])), len(numpy.unique(scaled_k[:,1])), @@ -98,8 +98,8 @@ def get_phase(cell, kpts, kmesh=[]): NR = len(R_vec_abs) phase = numpy.exp(1j*numpy.einsum('Ru, ku -> Rk', R_vec_abs, kpts)) phase /= numpy.sqrt(NR) # normalization in supercell - - # R_rel_mesh has to be construct exactly same to the Ts in super_cell function + +#R_rel_mesh has to be construct exactly same to the Ts in super_cell function scell = tools.super_cell(cell, kmesh) return scell, phase @@ -111,8 +111,8 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): C_k = numpy.asarray(mo_coeff) Nk, Nao, Nmo = C_k.shape NR = phase.shape[0] - - # Transform AO indices + +#Transform AO indices C_gamma = numpy.einsum('Rk, kum -> Rukm', phase, C_k) C_gamma = C_gamma.reshape(Nao*NR, Nk*Nmo) @@ -121,10 +121,10 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): E_g = E_g[E_sort_idx] C_gamma = C_gamma[:,E_sort_idx] s = scell.pbc_intor('int1e_ovlp') - #assert(abs(reduce(numpy.dot, (C_gamma.conj().T, s, C_gamma)) - # - numpy.eye(Nmo*Nk)).max() < 1e-7) - - # Transform MO indices +#assert(abs(reduce(numpy.dot, (C_gamma.conj().T, s, C_gamma)) +#- numpy.eye(Nmo * Nk)).max() <1e-7) + +#Transform MO indices E_k_degen = abs(E_g[1:] - E_g[:-1]).max() < 1e-5 if numpy.any(E_k_degen): print("Entered Strange If statement") @@ -132,24 +132,24 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): shift = min(E_g[degen_mask]) - .1 f = numpy.dot(C_gamma[:,degen_mask] * (E_g[degen_mask] - shift), C_gamma[:,degen_mask].conj().T) - #assert(abs(f.imag).max() < 1e-5) +#assert(abs(f.imag).max() < 1e-5) e, na_orb = la.eigh(f.real, s, type=2) C_gamma[:,degen_mask] = na_orb[:, e>0] - - #assert(abs(C_gamma.imag).max() < 1e-2) - #C_gamma = C_gamma.real - #assert(abs(reduce(numpy.dot, (C_gamma.conj().T, s, C_gamma)) - # - numpy.eye(Nmo*Nk)).max() < 1e-2) + +#assert(abs(C_gamma.imag).max() < 1e-2) +#C_gamma = C_gamma.real +#assert(abs(reduce(numpy.dot, (C_gamma.conj().T, s, C_gamma)) +#- numpy.eye(Nmo * Nk)).max() <1e-2) s_k = cell.pbc_intor('int1e_ovlp', kpts=kpts) - # overlap between k-point unitcell and gamma-point supercell +#overlap between k - point unitcell and gamma - point supercell s_k_g = numpy.einsum('kuv,Rk->kuRv', s_k, phase.conj()).reshape(Nk,Nao,NR*Nao) - # The unitary transformation from k-adapted orbitals to gamma-point orbitals +#The unitary transformation from k - adapted orbitals to gamma - point orbitals mo_phase = lib.einsum('kum,kuv,vi->kmi', C_k.conj(), s_k_g, C_gamma) C_gamma_unsorted=C_gamma[:,E_desort_idx] E_g_unsorted=E_g[E_desort_idx] - #return E_g_unsorted, C_gamma_unsorted, +#return E_g_unsorted, C_gamma_unsorted, return E_g, C_gamma, E_g_unsorted,C_gamma_unsorted @@ -182,14 +182,14 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): natom=loc_cell.natm dt = h5py.special_dtype(vlen=bytes) - #Group Atoms +#Group Atoms groupAtom=H5_qmcpack.create_group("atoms") - #Dataset Number Of Atoms +#Dataset Number Of Atoms groupAtom.create_dataset("number_of_atoms",(1,),dtype="i4",data=natom) - #Dataset Number Of Species - #Species contains (Atom_Name, Atom_Number,Atom_Charge,Atom_Core) +#Dataset Number Of Species +#Species contains(Atom_Name, Atom_Number, Atom_Charge, Atom_Core) l_atoms = [ (loc_cell.atom_pure_symbol(x),IonName[loc_cell.atom_pure_symbol(x)],loc_cell.atom_charge(x),loc_cell.atom_nelec_core(x)) for x in range(natom) ] @@ -213,12 +213,12 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): groupAtom.create_dataset("number_of_species",(1,),dtype="i4",data=NbSpecies) - #Dataset positions +#Dataset positions MyPos=groupAtom.create_dataset("positions",(natom,3),dtype="f8") for x in range(natom): MyPos[x:]=loc_cell.atom_coord(x) - #Group Atoms +#Group Atoms for x in range(NbSpecies): atmname=str(uniq_atoms[x][0]) groupSpecies=groupAtom.create_group("species_"+str(x)) @@ -240,8 +240,7 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): SpeciesID[x:] = idxAtomstoSpecies[x] - - #Parameter Group +#Parameter Group GroupParameter=H5_qmcpack.create_group("parameters") GroupParameter.create_dataset("ECP",(1,),dtype="b1",data=bool(loc_cell.has_ecp())) bohrUnit=True @@ -251,12 +250,12 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): GroupParameter.create_dataset("NbAlpha",(1,),dtype="i4",data=loc_cell.nelec[0]) GroupParameter.create_dataset("NbBeta",(1,),dtype="i4",data=loc_cell.nelec[1]) GroupParameter.create_dataset("NbTotElec",(1,),dtype="i4",data=loc_cell.nelec[0]+loc_cell.nelec[1]) - GroupParameter.create_dataset("spin",(1,),dtype="i4",data=Spin) - + GroupParameter.create_dataset("spin",(1,),dtype="i4",data=Spin) + - #basisset Group +#basisset Group GroupBasisSet=H5_qmcpack.create_group("basisset") - #Dataset Number Of Atoms +#Dataset Number Of Atoms GroupBasisSet.create_dataset("NbElements",(1,),dtype="i4",data=NbSpecies) if Python3: strList=['LCAOBSet'] @@ -266,7 +265,7 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): LCAOName=GroupBasisSet.create_dataset("name",(1,),dtype="S8") LCAOName[0:]="LCAOBSet" - #atomicBasisSets Group +#atomicBasisSets Group for x in range(NbSpecies): MyIdx=idxSpeciestoAtoms[x][0] @@ -385,59 +384,65 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): Val_l=BasisGroup.create_dataset("l",(1,),dtype="i4",data=l) Val_n=BasisGroup.create_dataset("n",(1,),dtype="i4",data=n) RadGroup=BasisGroup.create_group("radfunctions") - #print "" +#print "" IdRad=0 for e,c in zip(contracted_exp,line): DataRadGrp=RadGroup.create_group("DataRad"+str(IdRad)) DataRadGrp.create_dataset("exponent",(1,),dtype="f8",data=e) DataRadGrp.create_dataset("contraction",(1,),dtype="f8",data=c) - #print "" +#print "" IdRad+=1 n+=1 atomicBasisSetGroup.create_dataset("NbBasisGroups",(1,),dtype="i4",data=n) def is_complex(l): - try: - return is_complex(l[0]) - except: - return bool(l.imag) - - - - - - if loc_cell.cart==True: - # Generated from read_order.py in Numerics/codegen - d_gms_order = { - 0:[""], - 1:["x","y","z"], - 2:["xx","yy","zz","xy","xz","yz"], - 3:["xxx","yyy","zzz","xxy","xxz","yyx","yyz","zzx","zzy","xyz"], - 4:["xxxx","yyyy","zzzz","xxxy","xxxz","yyyx","yyyz","zzzx","zzzy","xxyy","xxzz","yyzz","xxyz","yyxz","zzxy"], - 5:["xxxxx","yyyyy","zzzzz","xxxxy","xxxxz","yyyyx","yyyyz","zzzzx","zzzzy","xxxyy","xxxzz","yyyxx","yyyzz","zzzxx","zzzyy","xxxyz","yyyxz","zzzxy","xxyyz","xxzzy","yyzzx"], - 6:["xxxxxx","yyyyyy","zzzzzz","xxxxxy","xxxxxz","yyyyyx","yyyyyz","zzzzzx","zzzzzy","xxxxyy","xxxxzz","yyyyxx","yyyyzz","zzzzxx","zzzzyy","xxxxyz","yyyyxz","zzzzxy","xxxyyy","xxxzzz","yyyzzz","xxxyyz","xxxzzy","yyyxxz","yyyzzx","zzzxxy","zzzyyx","xxyyzz"], - } - - d_l = {'s':0, 'p':1, 'd':2, 'f':3, 'g':4, 'h':5, 'i':6} - - def n_orbital(n): - if n==0: - return 1 - elif n==1: - return 3 - else: - return 2*n_orbital(n-1)-n_orbital(n-2)+1 - - - def compare_gamess_style(item1, item2): - # Warning: - # - d_gms_order is a global variable - n1,n2 = map(len,(item1,item2)) - assert (n1 == n2) - try: - l = d_gms_order[n1] +try: + return is_complex(l[0]) except : return bool(l.imag) + + + if loc_cell.cart == + True : +#Generated from read_order.py in Numerics / codegen + d_gms_order = { + 0 : [""], + 1 : [ "x", "y", "z" ], + 2 : [ "xx", "yy", "zz", "xy", "xz", "yz" ], + 3 : [ "xxx", "yyy", "zzz", "xxy", "xxz", "yyx", "yyz", "zzx", "zzy", "xyz" ], + 4 : [ + "xxxx", "yyyy", "zzzz", "xxxy", "xxxz", "yyyx", "yyyz", "zzzx", "zzzy", "xxyy", "xxzz", "yyzz", "xxyz", "yyxz", + "zzxy" + ], + 5 : [ + "xxxxx", "yyyyy", "zzzzz", "xxxxy", "xxxxz", "yyyyx", "yyyyz", "zzzzx", "zzzzy", "xxxyy", "xxxzz", + "yyyxx", "yyyzz", "zzzxx", "zzzyy", "xxxyz", "yyyxz", "zzzxy", "xxyyz", "xxzzy", "yyzzx" + ], + 6 : [ + "xxxxxx", "yyyyyy", "zzzzzz", "xxxxxy", "xxxxxz", "yyyyyx", "yyyyyz", "zzzzzx", "zzzzzy", "xxxxyy", + "xxxxzz", "yyyyxx", "yyyyzz", "zzzzxx", "zzzzyy", "xxxxyz", "yyyyxz", "zzzzxy", "xxxyyy", "xxxzzz", + "yyyzzz", "xxxyyz", "xxxzzy", "yyyxxz", "yyyzzx", "zzzxxy", "zzzyyx", "xxyyzz" + ], + } + + d_l = {'s' : 0, 'p' : 1, 'd' : 2, 'f' : 3, 'g' : 4, 'h' : 5, 'i' : 6} + + def + n_orbital(n) + : if n == 0 : return 1 elif n == + 1 : return 3 else : return 2 * n_orbital(n - 1) - n_orbital(n - 2) + + 1 + + + def + compare_gamess_style(item1, item2) + : +#Warning: +#- d_gms_order is a global variable + n1, + n2 = map(len, (item1, item2)) assert(n1 == n2) + try : l = d_gms_order[n1] except KeyError: return 0 else: @@ -449,24 +454,24 @@ def compare_python3(item1, item2): return compare_gamess_style(item1[0],item2[0]) ao_label = loc_cell.ao_labels(False) - - - # Create a list of shell + + +#Create a list of shell l_l = [] for label, name, t, l in ao_label: - # Change yyx -> xyy " +#Change yyx->xyy " q = "".join(sorted(l, key=l.count, reverse=True)) l_l.append(q) - - # Pyscf ordering of shell + +#Pyscf ordering of shell l_order = list(range(len(l_l))) - - # Shell ordering indexed + +#Shell ordering indexed n = 1 l_order_new = [] for i,(label, name, t, l) in enumerate(ao_label): r = d_l[t[-1]] - # print r,n_orbital(r) +#print r, n_orbital(r) if n != 1: n-=1 else: @@ -474,16 +479,16 @@ def compare_python3(item1, item2): n = n_orbital(r) unordered_l = l_l[i:i+n] unordered = l_order[i:i+n] - #print i,n,unordered +#print i, n, unordered ordered = [x for _,x in sorted(zip(unordered_l,unordered),key=cmp_to_key(compare_python3))] l_order_new.extend(ordered) def order_mo_coef(ll): - # Order a list of transposed mo_coeff (Ao,Mo) -> (Mo,Ao) ordered - # Warning: - # - l_order_new is used as global variable - # - gamess order +#Order a list of transposed mo_coeff(Ao, Mo)->(Mo, Ao)ordered +#Warning: +#- l_order_new is used as global variable +#- gamess order ll_new= [] for l in zip(*ll): @@ -522,14 +527,14 @@ def order_mo_coef(ll): eigenvalue_dn=GroupDet.create_dataset("eigenval_1",(1,NbMO),dtype="f8",data=mf.mo_energy[1]) else: - #Cell Parameters +#Cell Parameters GroupCell=H5_qmcpack.create_group("Cell") GroupCell.create_dataset("LatticeVectors",(3,3),dtype="f8",data=loc_cell.lattice_vectors()) def get_mo(mo_coeff, cart): return order_mo_coef(mo_coeff) if cart else list(zip(*mo_coeff)) - #Supertwist Coordinate +#Supertwist Coordinate GroupDet.create_dataset("Coord",(1,3),dtype="f8",data=sp_twist) if Gamma: @@ -566,8 +571,8 @@ def get_mo(mo_coeff, cart): eigenset=GroupDet.create_dataset("eigenset_0",(NbMO,NbAO),dtype="f8",data=mo_coeff_) eigenvalue=GroupDet.create_dataset("eigenval_0",(1,NbMO),dtype="f8",data=E_g) - - #Unsorted Mo_coeffs for Multideterminants order matching QP + +#Unsorted Mo_coeffs for Multideterminants order matching QP eigenset_unsorted=GroupDet.create_dataset("eigenset_unsorted_0",(NbMO,NbAO),dtype="f8",data=mo_coeff_unsorted) eigenvalue_unsorted=GroupDet.create_dataset("eigenval_unsorted_0",(1,NbMO),dtype="f8",data=E_g_unsorted) @@ -579,7 +584,7 @@ def get_mo(mo_coeff, cart): if is_multidet: make_multidet(cell, mf, title, H5_qmcpack) - # Close the file before exiting +#Close the file before exiting H5_qmcpack.close() print ('Wavefunction successfully saved to QMCPACK HDF5 Format') From aada0e3d8d6312770633aec55a2d3a486fa41cae Mon Sep 17 00:00:00 2001 From: NastaMauger Date: Mon, 30 Jun 2025 22:35:24 -0400 Subject: [PATCH 4/8] Improve output messages for multideterminant --- src/QMCTools/PyscfToQmcpack.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/QMCTools/PyscfToQmcpack.py b/src/QMCTools/PyscfToQmcpack.py index 0e0d5a2948..18c810f995 100644 --- a/src/QMCTools/PyscfToQmcpack.py +++ b/src/QMCTools/PyscfToQmcpack.py @@ -591,7 +591,7 @@ def get_mo(mo_coeff, cart): print ('Use: "convert4qmc -orbitals {}.h5" to generate QMCPACK input files'.format(title)) if is_multidet: - print(f'Multideterminant wavefunction saved to {title}_multidet.h5') + print(f'Multideterminant part saved to {title}_multidet.h5') def make_multidet(cell, mf, title, h5_handle): from pyscf import mcscf From 1f4b7f1cc18e295529cb09f9ee268a4afe589380 Mon Sep 17 00:00:00 2001 From: NastaMauger Date: Tue, 1 Jul 2025 10:48:11 -0400 Subject: [PATCH 5/8] Remove README_branch.md from current branch --- README_branch.md | 14 -------------- 1 file changed, 14 deletions(-) delete mode 100644 README_branch.md diff --git a/README_branch.md b/README_branch.md deleted file mode 100644 index 8df7a5fa4d..0000000000 --- a/README_branch.md +++ /dev/null @@ -1,14 +0,0 @@ -# Repository Name - -## Branches Description - -### `develop` branch -This branch contains the latest development work and merges from upstream. - -### `Batched_Converter` branch -This branch includes the work related to the Batched Converter feature. It is based on the forked repository from `anbenali`. - -### `LCAO_NumCenters_Opt` branch -This branch contains optimizations for the `LCAO_NumCenters_Opt` feature. It is based on the forked repository from `anbenali`. - - From 833315c4ae4992f51846313cd247b6e374ca02e2 Mon Sep 17 00:00:00 2001 From: NastaMauger Date: Tue, 1 Jul 2025 10:52:37 -0400 Subject: [PATCH 6/8] Remove formating --- src/QMCTools/PyscfToQmcpack.py | 265 ++++++++++++++++----------------- 1 file changed, 131 insertions(+), 134 deletions(-) diff --git a/src/QMCTools/PyscfToQmcpack.py b/src/QMCTools/PyscfToQmcpack.py index 18c810f995..5a9f30696f 100644 --- a/src/QMCTools/PyscfToQmcpack.py +++ b/src/QMCTools/PyscfToQmcpack.py @@ -16,7 +16,7 @@ def savetoqmcpack(cell,mf,title="Default",kpts=[],kmesh=[],sp_twist=[],weight=1. import h5py, re, sys from collections import defaultdict from pyscf.pbc import gto, scf, df, dft - from pyscf import lib + from pyscf import lib, mcscf, fci from pyscf.pbc import tools from numpy import empty import numpy @@ -56,7 +56,7 @@ def savetoqmcpack(cell,mf,title="Default",kpts=[],kmesh=[],sp_twist=[],weight=1. def get_supercell(cell,kmesh=[]): latt_vec = cell.lattice_vectors() if len(kmesh)==0: -#Guess kmesh + # Guess kmesh scaled_k = cell.get_scaled_kpts(kpts).round(8) kmesh = (len(numpy.unique(scaled_k[:,0])), len(numpy.unique(scaled_k[:,1])), @@ -69,7 +69,7 @@ def get_supercell(cell,kmesh=[]): R_vec_abs = numpy.einsum('nu, uv -> nv', R_vec_rel, latt_vec) -#R_rel_mesh has to be construct exactly same to the Ts in super_cell function + # R_rel_mesh has to be construct exactly same to the Ts in super_cell function scell = tools.super_cell(cell, kmesh) return scell, kmesh @@ -81,7 +81,7 @@ def get_phase(cell, kpts, kmesh=[]): latt_vec = cell.lattice_vectors() if kmesh is None: -#Guess kmesh + # Guess kmesh scaled_k = cell.get_scaled_kpts(kpts).round(8) kmesh = (len(numpy.unique(scaled_k[:,0])), len(numpy.unique(scaled_k[:,1])), @@ -98,8 +98,8 @@ def get_phase(cell, kpts, kmesh=[]): NR = len(R_vec_abs) phase = numpy.exp(1j*numpy.einsum('Ru, ku -> Rk', R_vec_abs, kpts)) phase /= numpy.sqrt(NR) # normalization in supercell - -#R_rel_mesh has to be construct exactly same to the Ts in super_cell function + + # R_rel_mesh has to be construct exactly same to the Ts in super_cell function scell = tools.super_cell(cell, kmesh) return scell, phase @@ -111,8 +111,8 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): C_k = numpy.asarray(mo_coeff) Nk, Nao, Nmo = C_k.shape NR = phase.shape[0] - -#Transform AO indices + + # Transform AO indices C_gamma = numpy.einsum('Rk, kum -> Rukm', phase, C_k) C_gamma = C_gamma.reshape(Nao*NR, Nk*Nmo) @@ -121,10 +121,10 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): E_g = E_g[E_sort_idx] C_gamma = C_gamma[:,E_sort_idx] s = scell.pbc_intor('int1e_ovlp') -#assert(abs(reduce(numpy.dot, (C_gamma.conj().T, s, C_gamma)) -#- numpy.eye(Nmo * Nk)).max() <1e-7) - -#Transform MO indices + #assert(abs(reduce(numpy.dot, (C_gamma.conj().T, s, C_gamma)) + # - numpy.eye(Nmo*Nk)).max() < 1e-7) + + # Transform MO indices E_k_degen = abs(E_g[1:] - E_g[:-1]).max() < 1e-5 if numpy.any(E_k_degen): print("Entered Strange If statement") @@ -132,24 +132,24 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): shift = min(E_g[degen_mask]) - .1 f = numpy.dot(C_gamma[:,degen_mask] * (E_g[degen_mask] - shift), C_gamma[:,degen_mask].conj().T) -#assert(abs(f.imag).max() < 1e-5) + #assert(abs(f.imag).max() < 1e-5) e, na_orb = la.eigh(f.real, s, type=2) C_gamma[:,degen_mask] = na_orb[:, e>0] - -#assert(abs(C_gamma.imag).max() < 1e-2) -#C_gamma = C_gamma.real -#assert(abs(reduce(numpy.dot, (C_gamma.conj().T, s, C_gamma)) -#- numpy.eye(Nmo * Nk)).max() <1e-2) + + #assert(abs(C_gamma.imag).max() < 1e-2) + #C_gamma = C_gamma.real + #assert(abs(reduce(numpy.dot, (C_gamma.conj().T, s, C_gamma)) + # - numpy.eye(Nmo*Nk)).max() < 1e-2) s_k = cell.pbc_intor('int1e_ovlp', kpts=kpts) -#overlap between k - point unitcell and gamma - point supercell + # overlap between k-point unitcell and gamma-point supercell s_k_g = numpy.einsum('kuv,Rk->kuRv', s_k, phase.conj()).reshape(Nk,Nao,NR*Nao) -#The unitary transformation from k - adapted orbitals to gamma - point orbitals + # The unitary transformation from k-adapted orbitals to gamma-point orbitals mo_phase = lib.einsum('kum,kuv,vi->kmi', C_k.conj(), s_k_g, C_gamma) C_gamma_unsorted=C_gamma[:,E_desort_idx] E_g_unsorted=E_g[E_desort_idx] -#return E_g_unsorted, C_gamma_unsorted, + #return E_g_unsorted, C_gamma_unsorted, return E_g, C_gamma, E_g_unsorted,C_gamma_unsorted @@ -182,14 +182,14 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): natom=loc_cell.natm dt = h5py.special_dtype(vlen=bytes) -#Group Atoms + #Group Atoms groupAtom=H5_qmcpack.create_group("atoms") -#Dataset Number Of Atoms + #Dataset Number Of Atoms groupAtom.create_dataset("number_of_atoms",(1,),dtype="i4",data=natom) -#Dataset Number Of Species -#Species contains(Atom_Name, Atom_Number, Atom_Charge, Atom_Core) + #Dataset Number Of Species + #Species contains (Atom_Name, Atom_Number,Atom_Charge,Atom_Core) l_atoms = [ (loc_cell.atom_pure_symbol(x),IonName[loc_cell.atom_pure_symbol(x)],loc_cell.atom_charge(x),loc_cell.atom_nelec_core(x)) for x in range(natom) ] @@ -213,12 +213,12 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): groupAtom.create_dataset("number_of_species",(1,),dtype="i4",data=NbSpecies) -#Dataset positions + #Dataset positions MyPos=groupAtom.create_dataset("positions",(natom,3),dtype="f8") for x in range(natom): MyPos[x:]=loc_cell.atom_coord(x) -#Group Atoms + #Group Atoms for x in range(NbSpecies): atmname=str(uniq_atoms[x][0]) groupSpecies=groupAtom.create_group("species_"+str(x)) @@ -240,7 +240,8 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): SpeciesID[x:] = idxAtomstoSpecies[x] -#Parameter Group + + #Parameter Group GroupParameter=H5_qmcpack.create_group("parameters") GroupParameter.create_dataset("ECP",(1,),dtype="b1",data=bool(loc_cell.has_ecp())) bohrUnit=True @@ -250,12 +251,12 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): GroupParameter.create_dataset("NbAlpha",(1,),dtype="i4",data=loc_cell.nelec[0]) GroupParameter.create_dataset("NbBeta",(1,),dtype="i4",data=loc_cell.nelec[1]) GroupParameter.create_dataset("NbTotElec",(1,),dtype="i4",data=loc_cell.nelec[0]+loc_cell.nelec[1]) - GroupParameter.create_dataset("spin",(1,),dtype="i4",data=Spin) + GroupParameter.create_dataset("spin",(1,),dtype="i4",data=Spin) + - -#basisset Group + #basisset Group GroupBasisSet=H5_qmcpack.create_group("basisset") -#Dataset Number Of Atoms + #Dataset Number Of Atoms GroupBasisSet.create_dataset("NbElements",(1,),dtype="i4",data=NbSpecies) if Python3: strList=['LCAOBSet'] @@ -265,7 +266,7 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): LCAOName=GroupBasisSet.create_dataset("name",(1,),dtype="S8") LCAOName[0:]="LCAOBSet" -#atomicBasisSets Group + #atomicBasisSets Group for x in range(NbSpecies): MyIdx=idxSpeciestoAtoms[x][0] @@ -384,65 +385,59 @@ def mo_k2gamma(cell, mo_energy, mo_coeff, kpts, kmesh=None): Val_l=BasisGroup.create_dataset("l",(1,),dtype="i4",data=l) Val_n=BasisGroup.create_dataset("n",(1,),dtype="i4",data=n) RadGroup=BasisGroup.create_group("radfunctions") -#print "" + #print "" IdRad=0 for e,c in zip(contracted_exp,line): DataRadGrp=RadGroup.create_group("DataRad"+str(IdRad)) DataRadGrp.create_dataset("exponent",(1,),dtype="f8",data=e) DataRadGrp.create_dataset("contraction",(1,),dtype="f8",data=c) -#print "" + #print "" IdRad+=1 n+=1 atomicBasisSetGroup.create_dataset("NbBasisGroups",(1,),dtype="i4",data=n) def is_complex(l): -try: - return is_complex(l[0]) except : return bool(l.imag) - - - if loc_cell.cart == - True : -#Generated from read_order.py in Numerics / codegen - d_gms_order = { - 0 : [""], - 1 : [ "x", "y", "z" ], - 2 : [ "xx", "yy", "zz", "xy", "xz", "yz" ], - 3 : [ "xxx", "yyy", "zzz", "xxy", "xxz", "yyx", "yyz", "zzx", "zzy", "xyz" ], - 4 : [ - "xxxx", "yyyy", "zzzz", "xxxy", "xxxz", "yyyx", "yyyz", "zzzx", "zzzy", "xxyy", "xxzz", "yyzz", "xxyz", "yyxz", - "zzxy" - ], - 5 : [ - "xxxxx", "yyyyy", "zzzzz", "xxxxy", "xxxxz", "yyyyx", "yyyyz", "zzzzx", "zzzzy", "xxxyy", "xxxzz", - "yyyxx", "yyyzz", "zzzxx", "zzzyy", "xxxyz", "yyyxz", "zzzxy", "xxyyz", "xxzzy", "yyzzx" - ], - 6 : [ - "xxxxxx", "yyyyyy", "zzzzzz", "xxxxxy", "xxxxxz", "yyyyyx", "yyyyyz", "zzzzzx", "zzzzzy", "xxxxyy", - "xxxxzz", "yyyyxx", "yyyyzz", "zzzzxx", "zzzzyy", "xxxxyz", "yyyyxz", "zzzzxy", "xxxyyy", "xxxzzz", - "yyyzzz", "xxxyyz", "xxxzzy", "yyyxxz", "yyyzzx", "zzzxxy", "zzzyyx", "xxyyzz" - ], - } - - d_l = {'s' : 0, 'p' : 1, 'd' : 2, 'f' : 3, 'g' : 4, 'h' : 5, 'i' : 6} - - def - n_orbital(n) - : if n == 0 : return 1 elif n == - 1 : return 3 else : return 2 * n_orbital(n - 1) - n_orbital(n - 2) + - 1 - - - def - compare_gamess_style(item1, item2) - : -#Warning: -#- d_gms_order is a global variable - n1, - n2 = map(len, (item1, item2)) assert(n1 == n2) - try : l = d_gms_order[n1] + try: + return is_complex(l[0]) + except: + return bool(l.imag) + + + + + + if loc_cell.cart==True: + # Generated from read_order.py in Numerics/codegen + d_gms_order = { + 0:[""], + 1:["x","y","z"], + 2:["xx","yy","zz","xy","xz","yz"], + 3:["xxx","yyy","zzz","xxy","xxz","yyx","yyz","zzx","zzy","xyz"], + 4:["xxxx","yyyy","zzzz","xxxy","xxxz","yyyx","yyyz","zzzx","zzzy","xxyy","xxzz","yyzz","xxyz","yyxz","zzxy"], + 5:["xxxxx","yyyyy","zzzzz","xxxxy","xxxxz","yyyyx","yyyyz","zzzzx","zzzzy","xxxyy","xxxzz","yyyxx","yyyzz","zzzxx","zzzyy","xxxyz","yyyxz","zzzxy","xxyyz","xxzzy","yyzzx"], + 6:["xxxxxx","yyyyyy","zzzzzz","xxxxxy","xxxxxz","yyyyyx","yyyyyz","zzzzzx","zzzzzy","xxxxyy","xxxxzz","yyyyxx","yyyyzz","zzzzxx","zzzzyy","xxxxyz","yyyyxz","zzzzxy","xxxyyy","xxxzzz","yyyzzz","xxxyyz","xxxzzy","yyyxxz","yyyzzx","zzzxxy","zzzyyx","xxyyzz"], + } + + d_l = {'s':0, 'p':1, 'd':2, 'f':3, 'g':4, 'h':5, 'i':6} + + def n_orbital(n): + if n==0: + return 1 + elif n==1: + return 3 + else: + return 2*n_orbital(n-1)-n_orbital(n-2)+1 + + + def compare_gamess_style(item1, item2): + # Warning: + # - d_gms_order is a global variable + n1,n2 = map(len,(item1,item2)) + assert (n1 == n2) + try: + l = d_gms_order[n1] except KeyError: return 0 else: @@ -454,24 +449,24 @@ def compare_python3(item1, item2): return compare_gamess_style(item1[0],item2[0]) ao_label = loc_cell.ao_labels(False) - - -#Create a list of shell + + + # Create a list of shell l_l = [] for label, name, t, l in ao_label: -#Change yyx->xyy " + # Change yyx -> xyy " q = "".join(sorted(l, key=l.count, reverse=True)) l_l.append(q) - -#Pyscf ordering of shell + + # Pyscf ordering of shell l_order = list(range(len(l_l))) - -#Shell ordering indexed + + # Shell ordering indexed n = 1 l_order_new = [] for i,(label, name, t, l) in enumerate(ao_label): r = d_l[t[-1]] -#print r, n_orbital(r) + # print r,n_orbital(r) if n != 1: n-=1 else: @@ -479,16 +474,16 @@ def compare_python3(item1, item2): n = n_orbital(r) unordered_l = l_l[i:i+n] unordered = l_order[i:i+n] -#print i, n, unordered + #print i,n,unordered ordered = [x for _,x in sorted(zip(unordered_l,unordered),key=cmp_to_key(compare_python3))] l_order_new.extend(ordered) def order_mo_coef(ll): -#Order a list of transposed mo_coeff(Ao, Mo)->(Mo, Ao)ordered -#Warning: -#- l_order_new is used as global variable -#- gamess order + # Order a list of transposed mo_coeff (Ao,Mo) -> (Mo,Ao) ordered + # Warning: + # - l_order_new is used as global variable + # - gamess order ll_new= [] for l in zip(*ll): @@ -527,14 +522,14 @@ def order_mo_coef(ll): eigenvalue_dn=GroupDet.create_dataset("eigenval_1",(1,NbMO),dtype="f8",data=mf.mo_energy[1]) else: -#Cell Parameters + #Cell Parameters GroupCell=H5_qmcpack.create_group("Cell") GroupCell.create_dataset("LatticeVectors",(3,3),dtype="f8",data=loc_cell.lattice_vectors()) def get_mo(mo_coeff, cart): return order_mo_coef(mo_coeff) if cart else list(zip(*mo_coeff)) -#Supertwist Coordinate + #Supertwist Coordinate GroupDet.create_dataset("Coord",(1,3),dtype="f8",data=sp_twist) if Gamma: @@ -571,61 +566,63 @@ def get_mo(mo_coeff, cart): eigenset=GroupDet.create_dataset("eigenset_0",(NbMO,NbAO),dtype="f8",data=mo_coeff_) eigenvalue=GroupDet.create_dataset("eigenval_0",(1,NbMO),dtype="f8",data=E_g) - -#Unsorted Mo_coeffs for Multideterminants order matching QP + + #Unsorted Mo_coeffs for Multideterminants order matching QP eigenset_unsorted=GroupDet.create_dataset("eigenset_unsorted_0",(NbMO,NbAO),dtype="f8",data=mo_coeff_unsorted) eigenvalue_unsorted=GroupDet.create_dataset("eigenval_unsorted_0",(1,NbMO),dtype="f8",data=E_g_unsorted) + GroupParameter.create_dataset("numMO",(1,),dtype="i4",data=NbMO) GroupParameter.create_dataset("numAO",(1,),dtype="i4",data=NbAO) - + is_multidet = isinstance(mf, (mcscf.casci.CASCI, mcscf.mc1step.CASSCF)) + if is_multidet: make_multidet(cell, mf, title, H5_qmcpack) -#Close the file before exiting + # Close the file before exiting H5_qmcpack.close() print ('Wavefunction successfully saved to QMCPACK HDF5 Format') print ('Use: "convert4qmc -orbitals {}.h5" to generate QMCPACK input files'.format(title)) if is_multidet: - print(f'Multideterminant part saved to {title}_multidet.h5') + print(f'Multideterminant wavefunction saved to {title}_multidet.h5') + def make_multidet(cell, mf, title, h5_handle): - from pyscf import mcscf import numpy import h5py, re, sys - a = mf.fcisolver.large_ci(mf.ci, mf.ncas, mf.nelecas, tol=0.0, return_strs=True) - dets_a = [] - dets_b = [] - coeffs = [] - cas_mo_start_a = cell.nelec[0]-mf.nelecas[0] - cas_mo_start_b = cell.nelec[1]-mf.nelecas[1] - n = 64 # chunk length - for idx,i in enumerate(a): - occ_a = numpy.array(list(i[1][2:]),dtype=int) - occ_b = numpy.array(list(i[2][2:]),dtype=int) - string_a = '0'*(len(mf.mo_coeff) - cas_mo_start_a - len(occ_a)) + i[1][2:] + '1'*cas_mo_start_a - string_b = '0'*(len(mf.mo_coeff) - cas_mo_start_b - len(occ_b)) + i[2][2:] + '1'*cas_mo_start_b - chunks_a = [int(string_a[j:j+n],2) for j in range(0, len(string_a), n)] - chunks_b = [int(string_b[j:j+n],2) for j in range(0, len(string_b), n)] - dets_a.append(chunks_a) - dets_b.append(chunks_b) - coeffs.append(i[0]) - H5_qmcpack_multidet = h5py.File(title+'_multidet.h5','w') - groupApp=H5_qmcpack_multidet.create_group("MultiDet") - dets_a = numpy.array(dets_a) - dets_b = numpy.array(dets_b) - - dt = numpy.dtype(numpy.uint64) - groupApp.create_dataset('CI_Alpha',dets_a.shape,dtype=dt, data = dets_a) - groupApp.create_dataset('CI_Beta',dets_b.shape,dtype=dt, data = dets_b) - groupApp.create_dataset('Coeff', (len(coeffs),),dtype = float,data = coeffs) - groupApp.create_dataset('NbDet', (1,),dtype = "i4",data = len(coeffs)) - groupApp.create_dataset('Nbits', (1,),dtype = "i4",data = len(dets_a[0])) - groupApp.create_dataset('nstate', (1,),dtype = "i4",data = mf.mo_coeff.shape[0]) - groupApp.create_dataset('nexcitedstate', (1,),dtype = "i4",data = 2) - - H5_qmcpack_multidet.close() + a = mf.fcisolver.large_ci(mf.ci, mf.ncas, mf.nelecas, tol=0.0, return_strs=True) + dets_a = [] + dets_b = [] + coeffs = [] + cas_mo_start_a = cell.nelec[0]-mf.nelecas[0] + cas_mo_start_b = cell.nelec[1]-mf.nelecas[1] + n = 64 # chunk length + for idx,i in enumerate(a): + occ_a = numpy.array(list(i[1][2:]),dtype=int) + occ_b = numpy.array(list(i[2][2:]),dtype=int) + string_a = '0'*(len(mf.mo_coeff) - cas_mo_start_a - len(occ_a)) + i[1][2:] + '1'*cas_mo_start_a + string_b = '0'*(len(mf.mo_coeff) - cas_mo_start_b - len(occ_b)) + i[2][2:] + '1'*cas_mo_start_b + chunks_a = [int(string_a[j:j+n],2) for j in range(0, len(string_a), n)] + chunks_b = [int(string_b[j:j+n],2) for j in range(0, len(string_b), n)] + dets_a.append(chunks_a) + dets_b.append(chunks_b) + coeffs.append(i[0]) + H5_qmcpack_multidet = h5py.File(title+'_multidet.h5','w') + groupApp=H5_qmcpack_multidet.create_group("MultiDet") + dets_a = numpy.array(dets_a) + dets_b = numpy.array(dets_b) + + dt = numpy.dtype(numpy.uint64) + groupApp.create_dataset('CI_Alpha',dets_a.shape,dtype=dt, data = dets_a) + groupApp.create_dataset('CI_Beta',dets_b.shape,dtype=dt, data = dets_b) + groupApp.create_dataset('Coeff', (len(coeffs),),dtype = float,data = coeffs) + groupApp.create_dataset('NbDet', (1,),dtype = "i4",data = len(coeffs)) + groupApp.create_dataset('Nbits', (1,),dtype = "i4",data = len(dets_a[0])) + groupApp.create_dataset('nstate', (1,),dtype = "i4",data = mf.mo_coeff.shape[0]) + groupApp.create_dataset('nexcitedstate', (1,),dtype = "i4",data = 2) + + H5_qmcpack_multidet.close() From f20137577d9652fe13c9fef38154e4e5ae999ddd Mon Sep 17 00:00:00 2001 From: Ye Luo Date: Mon, 7 Jul 2025 13:35:29 -0500 Subject: [PATCH 7/8] Move print closer. --- src/QMCTools/PyscfToQmcpack.py | 4 +--- 1 file changed, 1 insertion(+), 3 deletions(-) diff --git a/src/QMCTools/PyscfToQmcpack.py b/src/QMCTools/PyscfToQmcpack.py index 5a9f30696f..b77796677e 100644 --- a/src/QMCTools/PyscfToQmcpack.py +++ b/src/QMCTools/PyscfToQmcpack.py @@ -580,6 +580,7 @@ def get_mo(mo_coeff, cart): if is_multidet: make_multidet(cell, mf, title, H5_qmcpack) + print(f'Multideterminant wavefunction saved to {title}_multidet.h5') # Close the file before exiting H5_qmcpack.close() @@ -587,9 +588,6 @@ def get_mo(mo_coeff, cart): print ('Wavefunction successfully saved to QMCPACK HDF5 Format') print ('Use: "convert4qmc -orbitals {}.h5" to generate QMCPACK input files'.format(title)) - if is_multidet: - print(f'Multideterminant wavefunction saved to {title}_multidet.h5') - def make_multidet(cell, mf, title, h5_handle): import numpy From d4b70afa38ebec44972f7108f4b5642157fb2292 Mon Sep 17 00:00:00 2001 From: NastaMauger Date: Tue, 8 Jul 2025 14:38:27 -0400 Subject: [PATCH 8/8] Correct formatting in mol.py --- utils/afqmctools/afqmctools/hamiltonian/mol.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/utils/afqmctools/afqmctools/hamiltonian/mol.py b/utils/afqmctools/afqmctools/hamiltonian/mol.py index 32d94add73..aabdfc1187 100644 --- a/utils/afqmctools/afqmctools/hamiltonian/mol.py +++ b/utils/afqmctools/afqmctools/hamiltonian/mol.py @@ -196,7 +196,7 @@ def chunked_cholesky(mol, max_error=1e-6, verbose=False, cmax=10): nu = numpy.argmax(diag) delta_max = diag[nu] if verbose: - print(" # Generating Cholesky decomposition of ERIs."%nchol_max) + print(" # Generating Cholesky decomposition of ERIs. nchol_max = %d "%nchol_max) print(" # max number of cholesky vectors = %d"%nchol_max) header = ['iteration', 'max_residual', 'time'] print(format_fixed_width_strings(header))