LCPQ/qmcchem
10
1
Fork 1
mirror of https://gitlab.com/scemama/qmcchem.git synced 2024-06-02 11:25:18 +02:00
Code Releases Activity
12d8e84ddf
qmcchem/src/test_SVD/h2o/QMCCHEM_withSVD.py
Abdallah Ammar 12d8e84ddf build_psi_SVD
2021-05-26 09:34:22 +02:00

481 lines
18 KiB
Python
Executable File
Raw Blame History

# !!!
import sys, os
QMCCHEM_PATH=os.environ["QMCCHEM_PATH"]
sys.path.insert(0,QMCCHEM_PATH+"/EZFIO/Python/")
# !!!
from ezfio import ezfio
from math import sqrt
from datetime import datetime
import time
import numpy as np
import subprocess
from scipy.linalg import eig, eigh
from RSVD import powit_RSVD
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
def get_energy():
buffer = subprocess.check_output(['qmcchem', 'result', '-e', 'e_loc', EZFIO_file], encoding='UTF-8')
if buffer.strip() != "":
buffer = buffer.splitlines()[-1]
_, energy, error = [float(x) for x in buffer.split()]
return energy, error
else:
return None, None
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
def get_Ci_h_matrix_svd():
# !!!
Ci_h_matrix_svd = np.zeros( (n_svd,n_svd) )
# !!!
beg_Ci_h_matrix_svd = results.find('Ci_h_matrix_svd : [ ') + len( 'Ci_h_matrix_svd : [ ' )
end_Ci_h_matrix_svd = len(results)
Ci_h_matrix_svd_buf = results[beg_Ci_h_matrix_svd:end_Ci_h_matrix_svd]
Ci_h_matrix_svd_buf = Ci_h_matrix_svd_buf.split( '\n' )
# !!!
for iline in range(1, n_svd**2+1):
# !!!
line = Ci_h_matrix_svd_buf[iline].split()
indc = int( line[0] )
errS = float( line[4] )
#if( errS>eps ):
#print( line )
if( indc != iline ):
print('Error in reading Ci_h_matrix_svd')
break
else:
#Ci_h_matrix_svd[indc-1] = float( line[2] )
irow = indc % n_svd
icol = indc // n_svd
if( irow!=0 ):
Ci_h_matrix_svd[irow-1][icol] = float( line[2] )
else:
Ci_h_matrix_svd[n_svd-1][icol-1] = float( line[2] )
# !!!
# !!!
# Ci_h_matrix_svd = np.reshape(Ci_h_matrix_svd, (n_svd, n_svd), order='F')
# !!!
return(Ci_h_matrix_svd)
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
def get_Ci_overlap_matrix_svd():
# !!!
Ci_overlap_matrix_svd = np.zeros( (n_svd,n_svd) )
# !!!
beg_Ci_overlap_matrix_svd = results.find('Ci_overlap_matrix_svd : [ ') + len( 'Ci_overlap_matrix_svd : [ ' )
end_Ci_overlap_matrix_svd = len(results)
Ci_overlap_matrix_svd_buf = results[beg_Ci_overlap_matrix_svd:end_Ci_overlap_matrix_svd]
Ci_overlap_matrix_svd_buf = Ci_overlap_matrix_svd_buf.split( '\n' )
# !!!
for iline in range(1, n_svd**2+1):
# !!!
line = Ci_overlap_matrix_svd_buf[iline].split()
indc = int( line[0] )
# !!!
errS = float( line[4] )
#if( errS>eps ):
#print( line )
if( indc != iline ):
print('Error in reading Ci_overlap_matrix_svd')
break
# !!!
# !!!
else:
#Ci_overlap_matrix_svd[indc-1] = float( line[2] )
irow = indc % n_svd
icol = indc // n_svd
if( irow!=0 ):
Ci_overlap_matrix_svd[irow-1][icol] = float( line[2] )
else:
Ci_overlap_matrix_svd[n_svd-1][icol-1] = float( line[2] )
# !!!
# !!!
# !!!
#Ci_overlap_matrix_svd = np.reshape(Ci_overlap_matrix_svd, (n_svd, n_svd), order='F')
# !!!
return(Ci_overlap_matrix_svd)
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
def get_Ci_h_matrix_postsvd():
#file = open('verif_order.txt','w')
# !!!
Ci_h_matrix_postsvd = np.zeros( (n_svd*n_svd , n_svd*n_svd) )
# !!!
beg_Ci_h_matrix_postsvd = results.find('Ci_h_matrix_postsvd : [ ') + len( 'Ci_h_matrix_postsvd : [ ' )
end_Ci_h_matrix_postsvd = len(results)
Ci_h_matrix_postsvd_buf = results[beg_Ci_h_matrix_postsvd:end_Ci_h_matrix_postsvd]
Ci_h_matrix_postsvd_buf = Ci_h_matrix_postsvd_buf.split( '\n' )
# !!!
for iline in range(1, n_svd**4+1):
# !!!
line = Ci_h_matrix_postsvd_buf[iline].split()
indc = int( line[0] )
errS = float( line[4] )
#if( errS>eps ):
#print( line )
if( indc != iline ):
print('Error in reading Ci_h_matrix_postsvd')
break
else:
# !!!
kp = indc % n_svd
if( ( indc % n_svd ) !=0 ):
kp = indc % n_svd
else:
kp = n_svd
indc1 = int( ( indc - kp ) / n_svd )
k = indc1 % n_svd + 1
indc2 = int( ( indc1 - (k-1) ) / n_svd )
lp = indc2 % n_svd + 1
l = int( ( indc2 - (lp-1) ) / n_svd ) + 1
# !!!
#indcrep = kp + (k-1)*n_svd + (lp-1)*n_svd**2 + (l-1)*n_svd**3
#file.write( '{:5} {:5} {:5} {:5} {:5} {:5} \n'.format(indc, indc-indcrep, kp, k, lp, l ) )
# !!!
irow = kp + (k-1)*n_svd - 1
icol = lp + (l-1)*n_svd - 1
Ci_h_matrix_postsvd[irow][icol] = float( line[2] )
#Ci_h_matrix_postsvd[indc-1] = float( line[2] )
# !!!
#Ci_h_matrix_postsvd = np.reshape(Ci_h_matrix_postsvd, (n_svd*n_svd, n_svd*n_svd), order='F')
# !!!
#file.close()
return(Ci_h_matrix_postsvd)
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
def get_Ci_overlap_matrix_postsvd():
# !!!
Ci_overlap_matrix_postsvd = np.zeros( (n_svd*n_svd , n_svd*n_svd) )
# !!!
beg_Ci_overlap_matrix_postsvd = results.find('Ci_overlap_matrix_postsvd : [ ') + len( 'Ci_overlap_matrix_postsvd : [ ' )
end_Ci_overlap_matrix_postsvd = len(results)
Ci_overlap_matrix_postsvd_buf = results[beg_Ci_overlap_matrix_postsvd:end_Ci_overlap_matrix_postsvd]
Ci_overlap_matrix_postsvd_buf = Ci_overlap_matrix_postsvd_buf.split( '\n' )
# !!!
for iline in range(1, n_svd**4+1):
# !!!
line = Ci_overlap_matrix_postsvd_buf[iline].split()
indc = int( line[0] )
errS = float( line[4] )
#if( errS>eps ):
#print( line )
if( indc != iline ):
print('Error in reading Ci_overlap_matrix_postsvd')
break
else:
# !!!
#kp = indc % n_svd
if( ( indc % n_svd ) !=0 ):
kp = indc % n_svd
else:
kp = n_svd
indc1 = int( ( indc - kp ) / n_svd )
k = indc1 % n_svd + 1
indc2 = int( ( indc1 - (k-1) ) / n_svd )
lp = indc2 % n_svd + 1
l = int( ( indc2 - (lp-1) ) / n_svd ) + 1
# !!!
irow = kp + (k-1)*n_svd - 1
icol = lp + (l-1)*n_svd - 1
Ci_overlap_matrix_postsvd[irow][icol] = float( line[2] )
#Ci_overlap_matrix_postsvd[indc-1] = float( line[2] )
# !!!
# !!!
#Ci_overlap_matrix_postsvd = np.reshape(Ci_overlap_matrix_postsvd, (n_svd*n_svd, n_svd*n_svd), order='F')
# !!!
return(Ci_overlap_matrix_postsvd)
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
def check_symmetric(a, tol=1e-3):
return np.all(np.abs(a-a.T) < tol)
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
def save_results_to_resultsQMC():
file = open('resultsQMC.txt','a')
file.write('\n \n \n')
file.write('- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -\n \n')
file.write("Today's date: {}\n".format(datetime.now()))
file.write("EZFIO file = {}\n".format(EZFIO_file))
file.write('\n')
file.write( results )
file.write('\n')
file.write('+ + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + + +\n \n')
file.close()
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
def get_Esvd():
# !!!
# read CI_SVD matrices
Ci_h_matrix_svd = get_Ci_h_matrix_svd()
Ci_overlap_matrix_svd = get_Ci_overlap_matrix_svd()
#print( 'Ci_h_matrix_svd is symmetric ? {}' .format(check_symmetric(Ci_h_matrix_svd)) )
#print( 'Ci_overlap_matrix_svd is symmetric ? {}' .format(check_symmetric(Ci_overlap_matrix_svd)) )
# !!!
# symmetrise and diagonalise
aa = Ci_h_matrix_svd
aa = 0.5*( aa + aa.T )
bb = Ci_overlap_matrix_svd
eigvals_svd, vr = eig(aa, bb, left=False, right=True, overwrite_a=True, overwrite_b=True,
check_finite=True, homogeneous_eigvals=False)
#print( eigvals_svd + E_toadd )
recouvre_svd = np.abs(psi_svd_coeff @ vr)
ind_gssvd = np.argmax(recouvre_svd)
# !!!
E_svd = eigvals_svd[ind_gssvd] + E_toadd
return( E_svd, vr[:,ind_gssvd] )
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
def get_Epostsvd():
# !!!
# read CI_postSVD matrices
Ci_h_matrix_postsvd = get_Ci_h_matrix_postsvd()
Ci_overlap_matrix_postsvd = get_Ci_overlap_matrix_postsvd()
#print( 'Ci_h_matrix_postsvd is symmetric ? {}' .format(check_symmetric(Ci_h_matrix_postsvd)) )
#print( 'Ci_overlap_matrix_postsvd is symmetric ? {}' .format(check_symmetric(Ci_overlap_matrix_postsvd)) )
# !!!
# symmetrise and diagonalise
aa = Ci_h_matrix_postsvd
aa = 0.5*( aa + aa.T )
bb = Ci_overlap_matrix_postsvd
eigvals_postsvd, vr = eig(aa, bb, left=False, right=True, overwrite_a=True, overwrite_b=True,
check_finite=True, homogeneous_eigvals=False)
#print( eigvals_postsvd + E_toadd )
d_postsvd = np.diagflat(psi_svd_coeff)
d_postsvd = d_postsvd.reshape( (1,n_svd*n_svd) )
recouvre_postsvd = np.abs(d_postsvd @ vr)
ind_gspostsvd = np.argmax(recouvre_postsvd)
#print(recouvre_postsvd, ind_gspostsvd)
# !!!
#print( eigvals_postsvd.shape )
E_postsvd = eigvals_postsvd[ind_gspostsvd] + E_toadd
# !!!
return( E_postsvd, vr[:,ind_gspostsvd] )
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
def SVD_postsvd(sigma_postsvd):
# !!!
print( 'performing new SVD for the post SVD eigenvector:' )
# !!!
#sigma_postsvd = sigma_postsvd.reshape( (n_svd,n_svd) )
#sigma_postsvd = sigma_postsvd.reshape( n_svd, n_svd, order='F' )
#print( 'sigma_postsvd is symmetric ? {}' .format(check_symmetric(sigma_postsvd)) )
# !!!
sigma_postsvd_mat = np.zeros( (n_svd,n_svd) )
for indc in range(1, n_svd**2+1):
if( ( indc % n_svd ) !=0 ):
kp = indc % n_svd
else:
kp = n_svd
indc1 = int( ( indc - kp ) / n_svd )
k = indc1 % n_svd + 1
irow = kp + (k-1)*n_svd - 1
sigma_postsvd_mat[kp-1][k-1] = sigma_postsvd[irow]
sigma_postsvd = sigma_postsvd_mat
# !!!
# construct the new matrix Y
Y = U_svd @ sigma_postsvd @ V_svd.T
normY = np.linalg.norm(Y, ord='fro')
# !!!
# parameters of RSVD
rank = n_svd
npow = 10
nb_oversamp = 10
# !!!
# call RSV
U_postSVD, sigma_postsvd_diag, VT_postsvd = powit_RSVD(Y, rank, npow, nb_oversamp)
# !!!
# check precision
Y_SVD = np.dot( U_postSVD , np.dot( np.diag(sigma_postsvd_diag) , VT_postsvd ) )
energy = np.sum( np.square(sigma_postsvd_diag) ) / normY**2
err_SVD = 100. * np.linalg.norm( Y - Y_SVD, ord="fro") / normY
print('energy = {}, error = {}\n'.format(energy, err_SVD))
# !!!
return(U_postSVD, sigma_postsvd_diag, VT_postsvd)
# !!!
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
def get_Hsvd_QP(Hsvd_qp_txt):
Hsvd_qp = np.zeros( (n_svd,n_svd) )
Hsvd_qp_file = open(Hsvd_qp_txt, 'r')
for line in Hsvd_qp_file:
line = line.split()
i = int(line[0]) - 1
j = int(line[1]) - 1
Hsvd_qp[i,j] = float(line[2])
return(Hsvd_qp)
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
def get_Esvd_QP(Hsvd_qp):
# !!!
# symmetrise and diagonalise
aa = Hsvd_qp
aa = 0.5*( aa + aa.T )
bb = np.identity(n_svd)
eigvals_svd, vr = eig(aa, bb, left=False, right=True, overwrite_a=True, overwrite_b=True,
check_finite=True, homogeneous_eigvals=False)
recouvre_svd = np.abs(psi_svd_coeff @ vr)
ind_gssvd = np.argmax(recouvre_svd)
E_svd = eigvals_svd[ind_gssvd] + E_toadd
return( E_svd, vr[:,ind_gssvd] )
# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ !
if __name__ == '__main__':
t0 = time.time()
# !!!
# ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ #
EZFIO_file = "/home/aammar/qp2/src/svdwf/h2o_631g_frez_nsvd10"
E_toadd = 9.194966082434476 #6.983610961797779
# ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ #
# !!!
ezfio.set_file(EZFIO_file)
n_svd = ezfio.get_spindeterminants_n_svd_coefs()
psi_svd_coeff = np.array(ezfio.get_spindeterminants_psi_svd_coefs())
U_svd = np.array(ezfio.get_spindeterminants_psi_svd_alpha())
V_svd = np.array(ezfio.get_spindeterminants_psi_svd_beta())
# !!!
#print( U_svd.shape )
U_svd = U_svd[0,:,:].T
#print( U_svd.shape )
V_svd = V_svd[0,:,:].T
# !!!
print("Today's date:", datetime.now() )
print("EZFIO file = {}".format(EZFIO_file))
print("nuclear energy = {}".format(E_toadd) )
print("n_svd = {} \n".format(n_svd) )
# !!!
#set_vmc_params()
#run_qmc()
#print("start QMC:")
#stop_qmc()
# !!!
print( 'getting QMCCHEM results from {}'.format(EZFIO_file) )
results = subprocess.check_output(['qmcchem', 'result', EZFIO_file], encoding='UTF-8')
# !!!
Eloc, ErrEloc = get_energy()
print('Eloc = {} +/- {}'.format(Eloc, ErrEloc))
print('{} <= Eloc <= {}\n'.format(Eloc-ErrEloc,Eloc+ErrEloc))
# !!!
save_resultsQMC = input( 'save QMC results from {}? (y/n) '.format(EZFIO_file) )
if( save_resultsQMC == 'y' ):
print('saving in resultsQMC.txt')
save_results_to_resultsQMC()
print('\n')
# !!!
read_QPsvd = input( 'read QP Hsvd matrix ? (y/n) ')
if( read_QPsvd == 'y' ):
Hsvd_qp_txt = input('name of file with QM H_svd matrix:')
Hsvd_qp = get_Hsvd_QP(Hsvd_qp_txt)
E_svd_QP, sigma_svd_QP = get_Esvd_QP(Hsvd_qp)
print('QP SVD enegry = {} \n'.format(E_svd_QP) )
print('\n')
# !!!
E_svd, sigma_svd = get_Esvd()
print('svd enegry = {} \n '.format(E_svd) )
# !!!
E_postsvd, sigma_postsvd = get_Epostsvd()
print('post svd energy = {} \n'.format(E_postsvd) )
# !!!
save_to_EZFIO = input( "modify EZFIO (with svd_QP, svd_QMC or postsvd)? " )
# !!!
if( save_to_EZFIO == 'svd_QMC' ):
ezfio.set_spindeterminants_psi_svd_coefs( sigma_svd )
direc_svdcoeff = EZFIO_file + '/spindeterminants/psi_svd_coefs.gz'
print( '{} is modified'.format(direc_svdcoeff) )
print( 'precedent svd coeff:' )
print(psi_svd_coeff)
print( 'current svd coeff:' )
print(sigma_svd)
# !!!
elif( save_to_EZFIO == 'svd_QP'):
ezfio.set_spindeterminants_psi_svd_coefs( sigma_svd_QP )
direc_svdcoeff = EZFIO_file + '/spindeterminants/psi_svd_coefs.gz'
print( '{} is modified'.format(direc_svdcoeff) )
print( 'precedent svd coeff:' )
print(psi_svd_coeff)
print( 'current svd coeff:' )
print(sigma_svd_QP)
# !!!
elif( save_to_EZFIO == 'postsvd' ):
# SVD first
U_postSVD, sigma_postsvd_diag, V_postSVD = SVD_postsvd(sigma_postsvd)
V_postSVD = V_postSVD.T
# save next in ezfio file
U_postSVD_toEZFIO = np.zeros( ( U_postSVD.shape[0], U_postSVD.shape[1], 1) )
V_postSVD_toEZFIO = np.zeros( ( V_postSVD.shape[0], V_postSVD.shape[1], 1) )
U_postSVD_toEZFIO[:,:,0] = U_postSVD
V_postSVD_toEZFIO[:,:,0] = V_postSVD
#
ezfio.set_spindeterminants_psi_svd_alpha( U_postSVD_toEZFIO )
ezfio.set_spindeterminants_psi_svd_coefs( sigma_postsvd_diag )
ezfio.set_spindeterminants_psi_svd_beta( V_postSVD_toEZFIO )
# !!!
else:
print("end after {:.3f} minutes".format((time.time()-t0)/60.) )
exit()
# !!!
print("end after {:.3f} minutes".format((time.time()-t0)/60.) )
# !!!
# !!!
View Git Blame Copy Permalink