PTEROSOR/QuantumPackage
10
0
Fork 0
mirror of https://github.com/QuantumPackage/qp2.git synced 2024-11-19 04:22:32 +01:00
Code Releases Activity

working on pw integral converter

Browse Source
This commit is contained in:
Kevin Gasperich 2022-08-24 14:19:55 -05:00
parent 7b3b10e143
commit 1a9a6b240d
1 changed files with 103 additions and 32 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

135
bin/qp_convert_h5_to_ezfio
View File

@ -22,12 +22,29 @@ import gzip
#fname = sys.argv[1] #fname = sys.argv[1]
#qph5name = sys.argv[2] #qph5name = sys.argv[2]
def kconserv_p_from_qkk2_mk(qkk2,mk):
nk, nk2 = qkk2.shape
assert(nk == nk2)
kcon_p = np.zeros((nk,nk,nk),dtype=int)
for i in range(nk):
for j in range(nk):
for k in range(nk):
kcon_p[i,j,k] = qkk2[mk[j],qkk2[k,i]]
return kcon_p
def get_full_path(file_path): def get_full_path(file_path):
file_path = os.path.expanduser(file_path) file_path = os.path.expanduser(file_path)
file_path = os.path.expandvars(file_path) file_path = os.path.expandvars(file_path)
# file_path = os.path.abspath(file_path) # file_path = os.path.abspath(file_path)
return file_path return file_path
def make_reim_identity_kblocks(nk,nm,na=None):
if na is None:
na = nm
single_block = np.eye(nm, na, dtype=np.complex128).view(dtype=np.float64).reshape((nm, na, 2))
kblocks = np.tile(single_block,[nk, 1, 1, 1])
return kblocks
def flatten(l): def flatten(l):
res = [] res = []
@ -174,7 +191,7 @@ def convert_mol(filename,qph5path):
return return
def convert_kpts_cd(filename,qph5path,qmcpack=True): def convert_kpts_cd(filename,qph5path,qmcpack=True,is_ao=True):
import json import json
ezfio.set_file(filename) ezfio.set_file(filename)
@ -186,20 +203,36 @@ def convert_kpts_cd(filename,qph5path,qmcpack=True):
ezfio.set_nuclei_nucl_coord( [ [0.], [0.], [0.] ]*nucl_num ) ezfio.set_nuclei_nucl_coord( [ [0.], [0.], [0.] ]*nucl_num )
ezfio.set_nuclei_nucl_label( ['He'] * nucl_num ) ezfio.set_nuclei_nucl_label( ['He'] * nucl_num )
with h5py.File(qph5path,'r') as qph5: with h5py.File(qph5path,'r') as qph5:
kpt_num = qph5['Hamiltonian']['KPoints'].shape[0] kpt_num = qph5['Hamiltonian/KPoints'][()].shape[0]
elec_alpha_num = qph5['Hamiltonian']['dims'][4] ham_dims = qph5['Hamiltonian/dims'][()]
elec_beta_num = qph5['Hamiltonian']['dims'][5] NMOPerKP = qph5['Hamiltonian/NMOPerKP'][()]
orb_num = qph5['Hamiltonian']['dims'][3] _, _, kpt_num, orb_num, elec_alpha_num_tot, elec_beta_num_tot, _, nchol_maybe = ham_dims
try:
is_ao = json.loads(qph5['metadata'][()].decode("utf-8").replace("'",'"'))['ortho_ao'] #for now, all kpts must have same number of MOs
if is_ao: for nmoi in NMOPerKP:
ao_num = orb_num if nmoi != NMOPerKP[0]:
elif is_ao ==False: print("ERROR: all KPs must have same number of MOs")
mo_num = orb_num raise ValueError
else:
raise ValueError('Problem with ortho_ao key in metadata') #TODO: fix na, nb in rmg
except: assert(elec_alpha_num_tot % kpt_num == 0)
raise UnicodeDecodeError('metadata not correctly parsed from HDF5 file') assert(elec_beta_num_tot % kpt_num == 0)
elec_alpha_num_per_kpt = elec_alpha_num_tot // kpt_num
elec_beta_num_per_kpt = elec_beta_num_tot // kpt_num
#elec_alpha_num_per_kpt = qph5['Hamiltonian']['dims'][4]
#elec_beta_num_per_kpt = qph5['Hamiltonian']['dims'][5]
#orb_num = qph5['Hamiltonian']['dims'][3]
#try:
# is_ao = json.loads(qph5['metadata'][()].decode("utf-8").replace("'",'"'))['ortho_ao']
# if is_ao:
# ao_num = orb_num
# elif is_ao ==False:
# mo_num = orb_num
# else:
# raise ValueError('Problem with ortho_ao key in metadata')
#except:
# raise UnicodeDecodeError('metadata not correctly parsed from HDF5 file')
ezfio.set_nuclei_kpt_num(kpt_num) ezfio.set_nuclei_kpt_num(kpt_num)
kpt_pair_num = (kpt_num*kpt_num + kpt_num)//2 kpt_pair_num = (kpt_num*kpt_num + kpt_num)//2
@ -210,12 +243,20 @@ def convert_kpts_cd(filename,qph5path,qmcpack=True):
ezfio.set_nuclei_nucl_num(nucl_num_per_kpt) ezfio.set_nuclei_nucl_num(nucl_num_per_kpt)
# these are totals (kpt_num * num_per_kpt) # these are totals (kpt_num * num_per_kpt)
# need to change if we want to truncate orbital space within pyscf # need to change if we want to truncate orbital space within pyscf
ezfio.set_ao_basis_ao_num(ao_num) #if is_ao:
ezfio.set_mo_basis_mo_num(mo_num) # ao_num = orb_num*kpt_num
ezfio.set_ao_basis_ao_num_per_kpt(ao_num//kpt_num) #TODO: fix this?
ezfio.set_mo_basis_mo_num_per_kpt(mo_num//kpt_num) ao_num_tot = orb_num
ezfio.electrons_elec_alpha_num = elec_alpha_num ao_num_per_kpt = NMOPerKP[0]
ezfio.electrons_elec_beta_num = elec_beta_num mo_num_tot = orb_num
mo_num_per_kpt = NMOPerKP[0]
#mo_num_per_kpt = ao_num_per_kpt
ezfio.set_ao_basis_ao_num(ao_num_per_kpt * kpt_num)
ezfio.set_mo_basis_mo_num(mo_num_per_kpt * kpt_num)
ezfio.set_ao_basis_ao_num_per_kpt(ao_num_per_kpt)
ezfio.set_mo_basis_mo_num_per_kpt(mo_num_per_kpt)
ezfio.electrons_elec_alpha_num = elec_alpha_num_per_kpt * kpt_num
ezfio.electrons_elec_beta_num = elec_beta_num_per_kpt * kpt_num
########################################## ##########################################
# # # #
# Basis # # Basis #
@ -228,13 +269,38 @@ def convert_kpts_cd(filename,qph5path,qmcpack=True):
# # # #
########################################## ##########################################
#TODO #TODO
#coef_per_kpt = np.eye(mo_num_per_kpt, ao_num_per_kpt, dtype=np.complex128).view(dtype=np.float64).reshape((mo_num_per_kpt, ao_num_per_kpt, 2))
#mo_coef_kpts = np.tile(coef_per_kpt,[kpt_num, 1, 1, 1])
#qph5.create_dataset('mo_basis/mo_coef_kpts',data=make_reim_identity_kblocks(kpt_num, mo_num_per_kpt, ao_num_per_kpt))
ezfio.set_mo_basis_mo_coef_kpts(make_reim_identity_kblocks(kpt_num, mo_num_per_kpt, ao_num_per_kpt))
########################################## ##########################################
# # # #
# Integrals Mono # # Integrals Mono #
# # # #
########################################## ##########################################
with h5py.File(qph5path,'r') as qph5:
# we don't have separate kinetic, nuc-elec, pseudo 1e ints, so just combine in nuc-elec and set rest to zero
mono_ints_tot = np.zeros((kpt_num,ao_num_per_kpt,ao_num_per_kpt,2),dtype=np.float64)
for i in range(kpt_num):
mono_ints_tot[i] = qph5[f'Hamiltonian/H1_kp{i}'][()]
ovlp_ao_reim = make_reim_identity_kblocks(kpt_num,ao_num_per_kpt,ao_num_per_kpt)
kin_ao_reim = np.zeros((kpt_num,ao_num_per_kpt,ao_num_per_kpt,2),dtype=np.float64)
ne_ao_reim = mono_ints_tot
ezfio.set_ao_one_e_ints_ao_integrals_kinetic_kpts(kin_ao_reim)
ezfio.set_ao_one_e_ints_ao_integrals_overlap_kpts(ovlp_ao_reim)
ezfio.set_ao_one_e_ints_ao_integrals_n_e_kpts(ne_ao_reim)
ezfio.set_ao_one_e_ints_io_ao_integrals_kinetic('Read')
ezfio.set_ao_one_e_ints_io_ao_integrals_overlap('Read')
ezfio.set_ao_one_e_ints_io_ao_integrals_n_e('Read')
""" """
with h5py.File(qph5path,'r') as qph5: with h5py.File(qph5path,'r') as qph5:
if is_ao: if is_ao:
kin_ao_reim= kin_ao_reim=
ovlp_ao_reim= ovlp_ao_reim=
@ -267,16 +333,19 @@ def convert_kpts_cd(filename,qph5path,qmcpack=True):
########################################## ##########################################
#TODO #TODO
with h5py.File(qph5path,'r') as qph5: with h5py.File(qph5path,'r') as qph5:
kconserv = qph5['nuclei/kconserv'][()].tolist() #kconserv = qph5['nuclei/kconserv'][()].tolist()
minusk = qph5['Hamiltonian']['MinusK'][:]+1 #TODO: change this after rmg is fixed
QKTok2 = qph5['Hamiltonian']['QPTok2'][:]+1 #minusk = qph5['Hamiltonian']['MinusK'][:]+1
QKTok2 = qph5['Hamiltonian']['QKTok2'][:]+1
minusk = QKTok2[:,0]
kconserv = kconserv_p_from_qkk2_mk(QKTok2-1,minusk-1)+1
unique_kpt_num = len(qph5['Hamiltonian']['KPFactorized']) unique_kpt_num = len(qph5['Hamiltonian']['KPFactorized'])
unique_k_idx = [] unique_k_idx = []
for i in qph5['Hamiltonian']['KPFactorized'].keys(): for i in qph5['Hamiltonian']['KPFactorized'].keys():
unique_k_idx.append(int(i[1:])+1) unique_k_idx.append(int(i[1:])+1)
kpt_sparse_map = np.zeros(kpt_num) kpt_sparse_map = np.zeros(kpt_num)
for i in range(kpt_num): for i in range(kpt_num):
if i+1 in uniq_k_idx: if i+1 in unique_k_idx:
kpt_sparse_map[i] = i+1 kpt_sparse_map[i] = i+1
else: else:
kpt_sparse_map[i] = -minusk[i] kpt_sparse_map[i] = -minusk[i]
@ -298,16 +367,18 @@ def convert_kpts_cd(filename,qph5path,qmcpack=True):
# should this be in ao_basis? ao_two_e_ints? # should this be in ao_basis? ao_two_e_ints?
with h5py.File(qph5path,'r') as qph5: with h5py.File(qph5path,'r') as qph5:
nchol_per_kpt = qph5['Hamiltonian']['NCholPerKP'][:] nchol_per_kpt = qph5['Hamiltonian']['NCholPerKP'][:]
print(nchol_per_kpt)
nchol_per_kpt = nchol_per_kpt[nchol_per_kpt != 0] nchol_per_kpt = nchol_per_kpt[nchol_per_kpt != 0]
print(nchol_per_kpt)
nchol_per_kpt_max = max(nchol_per_kpt) nchol_per_kpt_max = max(nchol_per_kpt)
ezfio.set_chol_num(nchol_per_kpt) ezfio.set_ao_two_e_ints_chol_num(nchol_per_kpt)
ezfio.set_chol_num_max(nchol_per_kpt_max) ezfio.set_ao_two_e_ints_chol_num_max(nchol_per_kpt_max)
if is_ao: if is_ao:
ao_num_per_kpt = ao_num//kpt_num #ao_num_per_kpt = ao_num//kpt_num
ezfio.set_io_chol_ao_integrals('Read') ezfio.set_ao_two_e_ints_io_chol_ao_integrals('Read')
#ao_chol_two_e_ints = np.zeros((2, ao_num_per_kpt, ao_num_per_kpt, nchol_per_kpt_max, kpt_num, len(nchol_per_kpt))) #ao_chol_two_e_ints = np.zeros((2, ao_num_per_kpt, ao_num_per_kpt, nchol_per_kpt_max, kpt_num, len(nchol_per_kpt)))
L_list = [] L_list = []
for i in len(nchol_per_kpt): for i in range(len(nchol_per_kpt)):
L = qph5['Hamiltonian']['KPFactorized'][f'L{i}'][:] L = qph5['Hamiltonian']['KPFactorized'][f'L{i}'][:]
L.reshape(kpt_num, ao_num_per_kpt, ao_num_per_kpt, nchol_per_kpt[i], 2) L.reshape(kpt_num, ao_num_per_kpt, ao_num_per_kpt, nchol_per_kpt[i], 2)
L = np.einsum("ijklm->ilkjm", A, B) L = np.einsum("ijklm->ilkjm", A, B)
@ -344,7 +415,7 @@ def convert_kpts_cd(filename,qph5path,qmcpack=True):
ezfio.set_mo_two_e_ints_df_mo_integrals_complex(dfmo_reim) ezfio.set_mo_two_e_ints_df_mo_integrals_complex(dfmo_reim)
ezfio.set_mo_two_e_ints_io_df_mo_integrals('Read') ezfio.set_mo_two_e_ints_io_df_mo_integrals('Read')
""" """
mo_num_per_kpt = ao_num//kpt_num #mo_num_per_kpt = ao_num//kpt_num
ezfio.set_io_chol_mo_integrals('Read') ezfio.set_io_chol_mo_integrals('Read')
#ao_chol_two_e_ints = np.zeros((2, ao_num_per_kpt, ao_num_per_kpt, nchol_per_kpt_max, kpt_num, len(nchol_per_kpt))) #ao_chol_two_e_ints = np.zeros((2, ao_num_per_kpt, ao_num_per_kpt, nchol_per_kpt_max, kpt_num, len(nchol_per_kpt)))
L_list = [] L_list = []