From 7a031b3356dcd1c50be1af774044d4ab28b83868 Mon Sep 17 00:00:00 2001 From: Abdallah Ammar Date: Sat, 31 Jul 2021 14:07:03 +0200 Subject: [PATCH] post QMC py --- devel/dmc_dress/dress_SVD.py | 174 +++++++++++++ devel/svdwf/SVD.slurm | 24 ++ devel/svdwf/cr2_work/LOBD_VMC.py | 261 ++++++++++++++++++++ devel/svdwf/cr2_work/LOBD_compactVMC.py | 225 +++++++++++++++++ devel/svdwf/cr2_work/LOBD_postVMC.py | 175 +++++++++++++ devel/svdwf/cr2_work/LOBD_postcompactSVD.py | 173 +++++++++++++ devel/svdwf/qmcchemRun.slurm | 31 +++ devel/svdwf/qpRun.slurm | 35 +++ devel/svdwf/run_py.slurm | 26 ++ 9 files changed, 1124 insertions(+) create mode 100755 devel/dmc_dress/dress_SVD.py create mode 100644 devel/svdwf/SVD.slurm create mode 100644 devel/svdwf/cr2_work/LOBD_VMC.py create mode 100644 devel/svdwf/cr2_work/LOBD_compactVMC.py create mode 100644 devel/svdwf/cr2_work/LOBD_postVMC.py create mode 100644 devel/svdwf/cr2_work/LOBD_postcompactSVD.py create mode 100644 devel/svdwf/qmcchemRun.slurm create mode 100644 devel/svdwf/qpRun.slurm create mode 100644 devel/svdwf/run_py.slurm diff --git a/devel/dmc_dress/dress_SVD.py b/devel/dmc_dress/dress_SVD.py new file mode 100755 index 0000000..0e52d38 --- /dev/null +++ b/devel/dmc_dress/dress_SVD.py @@ -0,0 +1,174 @@ +#!/usr/bin/env python + +import numpy as np +import subprocess +import sys, os +import time +from datetime import datetime + +QP_PATH=os.environ["QP_ROOT"]+"/external/ezfio/Python/" +sys.path.insert(0,QP_PATH) + +from ezfio import ezfio + + +#____________________________________________________________________________________ +# +def get_energy(): + energy_buffer = subprocess.check_output( + ['qmcchem', 'result', '-e', 'e_loc', EZFIO_file] + , encoding='UTF-8') + + if energy_buffer.strip() != "": + energy_buffer = energy_buffer.splitlines()[-1] + _, energy, error = [float(x) for x in energy_buffer.split()] + return energy, error + else: + return None, None +#____________________________________________________________________________________ + + + +#____________________________________________________________________________________ +# +def read_scalar(scalar): + + scal_txt = scalar + ' : ' + beg_read = results.find(scal_txt) + len(scal_txt) + scal_buf = results[beg_read:beg_read+100].split('\n')[0].split() + + scal_val , scal_err = float( scal_buf[0] ) , float( scal_buf[2] ) + + return(scal_val,scal_err) +#____________________________________________________________________________________ + + + + +#____________________________________________________________________________________ +# +def read_ci_dress_svd(): + + Ci_dress_svd_val = np.zeros((n_svd)) + Ci_dress_svd_err = np.zeros((n_svd)) + + beg_read = results.find('Ci_dress_svd : [ ') + len( 'Ci_dress_svd : [ ' ) + end_read = len(results) + vect_buf = results[beg_read:end_read] + vect_buf = vect_buf.split( '\n' ) + + for iline in range(1, n_svd+1): + line = vect_buf[iline].split() + + # check + indc = int( line[0] ) + if( indc != iline ): + print(' Error in reading Ci_dress_svd') + print(' indc = {}'.format(indc) ) + print(' iline = {}'.format(iline) ) + break + + else: + Ci_dress_svd_val[indc-1] = float( line[2] ) + Ci_dress_svd_err[indc-1] = float( line[4] ) + # !!! + # !!! + # !!! + + return(Ci_dress_svd_val,Ci_dress_svd_err) +#____________________________________________________________________________________ + + +#____________________________________________________________________________________ +# +def get_ci_dress(): + + U_svd = np.array(ezfio.get_spindeterminants_psi_svd_alpha()) + V_svd = np.array(ezfio.get_spindeterminants_psi_svd_beta()) + U_svd = U_svd[0,:,:].T + V_svd = V_svd[0,:,:].T + + ci_dress_mat_val = np.dot(U_svd,np.dot(np.diagflat(Ci_dress_svd_val),V_svd.T)) + ci_dress_mat_err = np.dot(U_svd,np.dot(np.diagflat(Ci_dress_svd_err),V_svd.T)) + print(ci_dress_mat_val.shape) + + n_det = ezfio.get_spindeterminants_n_det() + i_row = np.array(ezfio.get_spindeterminants_psi_coef_matrix_rows()) + j_col = np.array(ezfio.get_spindeterminants_psi_coef_matrix_columns()) + + Ci_dress_val = np.zeros((n_det)) + Ci_dress_err = np.zeros((n_det)) + for k in range(n_det): + i = i_row[k] - 1 + j = j_col[k] - 1 + Ci_dress_val[k] = ci_dress_mat_val[i,j] + Ci_dress_err[k] = ci_dress_mat_err[i,j] + + Ci_dress_val = Ci_dress_val / psi_norm + Ci_dress_err = Ci_dress_err / psi_norm + + return(Ci_dress_val,Ci_dress_err) +#____________________________________________________________________________________ + + + + +#____________________________________________________________________________________ +# +def save_to_file(v_val,v_err, name_f): + + with open(name_f, 'a') as ff: + for i in range(len(v_val)): + ff.write(' {} {} {}\n'.format(i, v_val[i], v_err[i])) +#____________________________________________________________________________________ + + + +if __name__ == '__main__': + t0 = time.time() + + # _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ + # + #EZFIO_file = sys.argv[1] + EZFIO_file = "/home/aammar/qp2/src/svdwf/cr2_work/cr2_1p6_dress" + # _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ + + ezfio.set_file(EZFIO_file) + + n_svd = ezfio.get_spindeterminants_n_svd_coefs() + + print(" ") + print(" Today's date:", datetime.now() ) + print(" EZFIO file = {}".format(EZFIO_file)) + print(" n_svd = {}\n".format(n_svd) ) + + t_read = time.time() + print(' getting QMCCHEM results from {}'.format(EZFIO_file) ) + results = subprocess.check_output(['qmcchem', 'result', EZFIO_file], encoding='UTF-8') + print(' getting results after {} minutes \n'.format( (time.time()-t_read)/60. )) + + # read QMC results + E_loc, ErrE_loc = get_energy() + print(' E_loc = {} +/- {}'.format(E_loc, ErrE_loc)) + print(' {} <= E_loc <= {}\n'.format(E_loc-ErrE_loc,E_loc+ErrE_loc)) + + E_loc_qmcvar, ErrE_loc_qmcvar = read_scalar('E_loc_qmcvar') + print(' E_loc_qmcvar = {} +/- {} \n'.format(E_loc_qmcvar,ErrE_loc_qmcvar)) + + psi_norm, Errpsi_norm = read_scalar('Psi_norm') + print(' psi_norm = {} +/- {}\n'.format(psi_norm,Errpsi_norm)) + + psi_norm_qmcvar, Errpsi_norm_qmcvar = read_scalar('Psi_norm_qmcvar') + print(' Psi_norm_qmcvar = {} +/- {} \n'.format(psi_norm_qmcvar, Errpsi_norm_qmcvar)) + + Ci_dress_svd_val, Ci_dress_svd_err = read_ci_dress_svd() + save_to_file(Ci_dress_svd_val,Ci_dress_svd_err, 'ci_dress_svd.txt') + + Ci_dress_val, Ci_dress_err = get_ci_dress() + save_to_file(Ci_dress_val,Ci_dress_err, 'ci_dress.txt') + #ezfio.set_dmc_dress_dmc_delta_h(Ci_dress_val) + + + print("end after {:.3f} minutes".format((time.time()-t0)/60.) ) + # !!! +# !!! diff --git a/devel/svdwf/SVD.slurm b/devel/svdwf/SVD.slurm new file mode 100644 index 0000000..18eaf05 --- /dev/null +++ b/devel/svdwf/SVD.slurm @@ -0,0 +1,24 @@ +#!/bin/bash +# +#SBATCH --job-name=SLURM_h2o_631g_SVD +#SBATCH --output=SLURM_h2o_631g_SVD.out +# +#SBATCH --mail-type=ALL +#SBATCH --mail-user=aammar@irsamc.ups-tlse.fr +# +#SBATCH -N 1 +#SBATCH -n 1 +#SBATCH -c 4 +#SBATCH -t 0:10:00 +##SBATCH --exclusive +# +# +srun hostname +echo "Running on `hostname`" +# +echo "Kicking off $SLURM_NPROCS compute servers on $SLURM_NNODES nodes" +echo "Nodes are\n `scontrol show hostname`" +# +~/qp2/bin/qpsh +qp_run Evar_TruncSVD h2o_631g >> h2o_631g_SVD.out +# diff --git a/devel/svdwf/cr2_work/LOBD_VMC.py b/devel/svdwf/cr2_work/LOBD_VMC.py new file mode 100644 index 0000000..335b3c2 --- /dev/null +++ b/devel/svdwf/cr2_work/LOBD_VMC.py @@ -0,0 +1,261 @@ +# !!! +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 + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +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) ) + SErr_H = 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( indc != iline ): + print(' Error in reading Ci_h_matrix_svd') + print(' indc = {}'.format(indc) ) + print(' iline = {}'.format(iline) ) + break + else: + irow = indc % n_svd + icol = indc // n_svd + if( irow!=0 ): + Ci_h_matrix_svd[irow-1][icol] = float( line[2] ) + SErr_H [irow-1][icol] = float( line[4] ) + + #ci_tmp = Ci_h_matrix_svd[irow-1][icol] + #er_tmp = SErr_H [irow-1][icol] + #if(ci_tmp>0): + # Ci_h_matrix_svd[irow-1][icol] = max(0 , ci_tmp - er_tmp) + #else: + # Ci_h_matrix_svd[irow-1][icol] = min(0 , ci_tmp + er_tmp) + + else: + Ci_h_matrix_svd[n_svd-1][icol-1] = float( line[2] ) + SErr_H [n_svd-1][icol-1] = float( line[4] ) + + #ci_tmp = Ci_h_matrix_svd[n_svd-1][icol-1] + #er_tmp = SErr_H [n_svd-1][icol-1] + #if(ci_tmp>0): + # Ci_h_matrix_svd[n_svd-1][icol-1] = max(0 , ci_tmp - er_tmp) + #else: + # Ci_h_matrix_svd[n_svd-1][icol-1] = min(0 , ci_tmp + er_tmp) + + return(Ci_h_matrix_svd,SErr_H) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +def get_Ci_overlap_matrix_svd(): + + Ci_overlap_matrix_svd = np.zeros( (n_svd,n_svd) ) + SErr_O = 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( indc != iline ): + print(' Error in reading Ci_overlap_matrix_svd') + print(' indc = {}'.format(indc)) + print(' iline = {}'.format(iline)) + break + else: + irow = indc % n_svd + icol = indc // n_svd + if( irow!=0 ): + Ci_overlap_matrix_svd[irow-1][icol] = float( line[2] ) + SErr_O [irow-1][icol] = float( line[4] ) + + #ci_tmp = Ci_overlap_matrix_svd[irow-1][icol] + #er_tmp = SErr_O [irow-1][icol] + #if(ci_tmp>0): + # Ci_overlap_matrix_svd[irow-1][icol] = max(0 , ci_tmp - er_tmp) + #else: + # Ci_overlap_matrix_svd[irow-1][icol] = min(0 , ci_tmp + er_tmp) + + else: + Ci_overlap_matrix_svd[n_svd-1][icol-1] = float( line[2] ) + SErr_O [n_svd-1][icol-1] = float( line[4] ) + + #ci_tmp = Ci_overlap_matrix_svd[n_svd-1][icol-1] + #er_tmp = SErr_O [n_svd-1][icol-1] + #if(ci_tmp>0): + # Ci_overlap_matrix_svd[n_svd-1][icol-1] = max(0 , ci_tmp - er_tmp) + #else: + # Ci_overlap_matrix_svd[n_svd-1][icol-1] = min(0 , ci_tmp + er_tmp) + + # !!! + # !!! + # !!! + + return(Ci_overlap_matrix_svd,SErr_O) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +def check_symmetric(a, tol=1e-3): + return np.all(np.abs(a-a.T) < tol) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +def get_Esvd(): + + # read CI_SVD matrices + Ci_h_matrix_svd , SErr_H = get_Ci_h_matrix_svd() + Ci_overlap_matrix_svd, SErr_O = 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)) ) + + if( save_SVD_QMCmat ): + with open('SVD_QMCmat.txt', 'a') as ff: + for i in range(n_svd): + for j in range(n_svd): + ff.write('{} {} {} {} {} {} \n'.format(i,j + , Ci_h_matrix_svd[i,j] , SErr_H[i,j] + , Ci_overlap_matrix_svd[i,j] , SErr_O[i,j] )) + + # 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] ) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + + +if __name__ == '__main__': + + t0 = time.time() + + # ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ # + #EZFIO_file = "Cr2_TSVD30" + #EZFIO_file = "Cr2_TSVD100" + #EZFIO_file = "Cr2_TSVD200" + #EZFIO_file = "Cr2_TSVD924_0" + #EZFIO_file = "Cr2_TSVD924_1" + #EZFIO_file = "Cr2_TSVD924_2" + + #EZFIO_file = "Cr2_TSVD30_it1" + #EZFIO_file = "Cr2_TSVD100_it1" + #EZFIO_file = "Cr2_TSVD200_it1" + #EZFIO_file = "Cr2_TSVD924_0_it1" + #EZFIO_file = "Cr2_TSVD924_1_it1" + EZFIO_file = "Cr2_TSVD924_2_it1" + + #EZFIO_file = "Cr2_TSVD200_minmaxit1" + + #EZFIO_file = "Cr2_TSVD30_it2" + + E_toadd = 64.8242130309350 + # ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ # + + ezfio.set_file(EZFIO_file) + n_svd = ezfio.get_spindeterminants_n_svd_coefs() + psi_svd_coeff = np.array(ezfio.get_spindeterminants_psi_svd_coefs()) + + print(" ") + print(" Today's date:", datetime.now() ) + print(" EZFIO file = {}".format(EZFIO_file)) + print(" n_svd = {}\n".format(n_svd) ) + + t_read = time.time() + print(' getting QMCCHEM results from {}'.format(EZFIO_file) ) + results = subprocess.check_output(['qmcchem', 'result', EZFIO_file], encoding='UTF-8') + print(' getting results after {} minutes \n'.format( (time.time()-t_read)/60. )) + + if( len(results)<100): + print(' result is empty \n') + exit() + + Eloc, ErrEloc = get_energy() + print(' Eloc = {} +/- {}'.format(Eloc, ErrEloc)) + print(' {} <= Eloc <= {}'.format(Eloc-ErrEloc,Eloc+ErrEloc)) + #print(' initial svd coeff = {}\n'.format(psi_svd_coeff)) + + save_SVD_QMCmat = True + #save_SVD_QMCmat = False + + print('') + E_svd, sigma_svd = get_Esvd() + print(' svd enegry = {}'.format(E_svd) ) + sigma_svd = sigma_svd * np.sign(sigma_svd[0]) + #print(' svd coeff = {}'.format(sigma_svd)) + print('') + + # ___________________________________________________________________________ + # ___________________________________________________________________________ + # + save_to_EZFIO = 'svd_QMC' + #save_to_EZFIO = '' + + 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) ) + else: + print(' exit with no modifications in {} EZFIO'.format(EZFIO_file)) + # ___________________________________________________________________________ + # ___________________________________________________________________________ + # + + print("end after {:.3f} minutes".format((time.time()-t0)/60.) ) + + # !!! +# !!! + +~ diff --git a/devel/svdwf/cr2_work/LOBD_compactVMC.py b/devel/svdwf/cr2_work/LOBD_compactVMC.py new file mode 100644 index 0000000..dea6941 --- /dev/null +++ b/devel/svdwf/cr2_work/LOBD_compactVMC.py @@ -0,0 +1,225 @@ +# !!! +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 + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +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_modifsvd(): + + Ci_h_matrix_modifsvd = np.zeros( (n_modifsvd,n_modifsvd) ) + SErr_H = np.zeros( (n_modifsvd,n_modifsvd) ) + + beg_Ci_h_matrix_modifsvd = results.find('Ci_h_matrix_modifsvd : [ ') + len( 'Ci_h_matrix_modifsvd : [ ' ) + end_Ci_h_matrix_modifsvd = len(results) + Ci_h_matrix_modifsvd_buf = results[beg_Ci_h_matrix_modifsvd:end_Ci_h_matrix_modifsvd] + Ci_h_matrix_modifsvd_buf = Ci_h_matrix_modifsvd_buf.split( '\n' ) + + for iline in range(1, n_modifsvd**2+1): + + line = Ci_h_matrix_modifsvd_buf[iline].split() + indc = int( line[0] ) + errS = float( line[4] ) + if( indc != iline ): + print(' Error in reading Ci_h_matrix_modifsvd') + print(' indc = {}'.format(indc) ) + print(' iline = {}'.format(iline) ) + break + else: + irow = indc % n_modifsvd + icol = indc // n_modifsvd + if( irow!=0 ): + Ci_h_matrix_modifsvd[irow-1][icol] = float( line[2] ) + SErr_H [irow-1][icol] = float( line[4] ) + else: + Ci_h_matrix_modifsvd[n_modifsvd-1][icol-1] = float( line[2] ) + SErr_H [n_modifsvd-1][icol-1] = float( line[4] ) + # !!! + # !!! + + return(Ci_h_matrix_modifsvd,SErr_H) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +def get_Ci_overlap_matrix_modifsvd(): + + Ci_overlap_matrix_modifsvd = np.zeros( (n_modifsvd,n_modifsvd) ) + SErr_O = np.zeros( (n_modifsvd,n_modifsvd) ) + + beg_Ci_overlap_matrix_modifsvd = results.find('Ci_overlap_matrix_modifsvd : [ ') + len( 'Ci_overlap_matrix_modifsvd : [ ' ) + end_Ci_overlap_matrix_modifsvd = len(results) + Ci_overlap_matrix_modifsvd_buf = results[beg_Ci_overlap_matrix_modifsvd:end_Ci_overlap_matrix_modifsvd] + Ci_overlap_matrix_modifsvd_buf = Ci_overlap_matrix_modifsvd_buf.split( '\n' ) + + for iline in range(1, n_modifsvd**2+1): + + line = Ci_overlap_matrix_modifsvd_buf[iline].split() + indc = int( line[0] ) + errS = float( line[4] ) + if( indc != iline ): + print(' Error in reading Ci_overlap_matrix_modifsvd') + print(' indc = {}'.format(indc)) + print(' iline = {}'.format(iline)) + break + else: + irow = indc % n_modifsvd + icol = indc // n_modifsvd + if( irow!=0 ): + Ci_overlap_matrix_modifsvd[irow-1][icol] = float( line[2] ) + SErr_O [irow-1][icol] = float( line[4] ) + else: + Ci_overlap_matrix_modifsvd[n_modifsvd-1][icol-1] = float( line[2] ) + SErr_O [n_modifsvd-1][icol-1] = float( line[4] ) + + return(Ci_overlap_matrix_modifsvd,SErr_O) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +def check_symmetric(a, tol=1e-3): + return np.all(np.abs(a-a.T) < tol) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +def get_Esvd(): + + Ci_h_matrix_modifsvd , SErr_H = get_Ci_h_matrix_modifsvd() + Ci_overlap_matrix_modifsvd, SErr_O = get_Ci_overlap_matrix_modifsvd() + + if( save_compactSVD_QMCmat ): + with open('compactSVD_QMCmat.txt', 'a') as ff: + for i in range(n_modifsvd): + for j in range(n_modifsvd): + ff.write('{} {} {} {} {} {} \n'.format(i,j + , Ci_h_matrix_modifsvd[i,j] , SErr_H[i,j] + , Ci_overlap_matrix_modifsvd[i,j] , SErr_O[i,j] )) + + # symmetrise and diagonalise + aa = Ci_h_matrix_modifsvd + aa = 0.5*( aa + aa.T ) + bb = Ci_overlap_matrix_modifsvd + eigvals_modifsvd, vr = eig(aa, bb, left=False, right=True, overwrite_a=True, overwrite_b=True, + check_finite=True, homogeneous_eigvals=False) + + psi_tmp = psi_svd_coefs[0,0:n_modifsvd] + recouvre_modifsvd = np.abs(psi_tmp @ vr) + ind_gssvd = np.argmax(recouvre_modifsvd) + + E_modifsvd = eigvals_modifsvd[ind_gssvd] + E_toadd + + return( E_modifsvd, vr[:,ind_gssvd] ) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + + +if __name__ == '__main__': + + t0 = time.time() + + # ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ # + #EZFIO_file = "Cr2_compactSVD30" + + + + EZFIO_file = "Cr2_compactSVD30_it1" + + E_toadd = 64.8242130309350 + # ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ # + + ezfio.set_file(EZFIO_file) + n_svd = ezfio.get_spindeterminants_n_svd_coefs() + n_modifsvd = ezfio.get_spindeterminants_n_modifsvd_coefs() + psi_svd_coefs = np.array(ezfio.get_spindeterminants_psi_svd_coefs()) + + print(" ") + print(" Today's date:", datetime.now() ) + print(" EZFIO file = {}".format(EZFIO_file)) + print(" n_svd = {}".format(n_svd) ) + print(" n_modifsvd = {}\n".format(n_modifsvd) ) + + t_read = time.time() + print(' getting QMCCHEM results from {}'.format(EZFIO_file) ) + results = subprocess.check_output(['qmcchem', 'result', EZFIO_file], encoding='UTF-8') + print(' getting results after {} minutes \n'.format( (time.time()-t_read)/60. )) + + if( len(results)<100): + print(' result is empty \n') + exit() + + Eloc, ErrEloc = get_energy() + print(' Eloc = {} +/- {}'.format(Eloc, ErrEloc)) + print(' {} <= Eloc <= {}'.format(Eloc-ErrEloc,Eloc+ErrEloc)) + #print(' initial svd coeff = {}\n'.format(psi_svd_coefs)) + + save_compactSVD_QMCmat = True + + print('') + E_compactsvd, sigma_compactsvd = get_Esvd() + print(' compact svd enegry = {}'.format(E_compactsvd) ) + sigma_compactsvd = sigma_compactsvd * np.sign(sigma_compactsvd[0]) + print(' compact svd coeff = {}\n'.format(sigma_compactsvd)) + + #print(' compact svd factor = {}\n'.format(sigma_compactsvd[n_modifsvd-1])) + + # ___________________________________________________________________________ + # ___________________________________________________________________________ + # + #save_to_EZFIO = 'compactsvd_QMC' + save_to_EZFIO = '' + + if( save_to_EZFIO == 'compactsvd_QMC' ): + ezfio_sigma_svd = np.zeros((n_svd)) + + for i in range(n_modifsvd-1): + ezfio_sigma_svd[i] = sigma_compactsvd[i] + for i in range(n_modifsvd-1, n_svd): + ezfio_sigma_svd[i] = psi_svd_coefs[0,i] * sigma_compactsvd[n_modifsvd-1] + + #print(' new svd coeff = {}\n'.format(ezfio_sigma_svd)) + + ezfio.set_spindeterminants_psi_svd_coefs( ezfio_sigma_svd ) + direc_svdcoeff = EZFIO_file + '/spindeterminants/psi_svd_coefs.gz' + print(' {} is modified'.format(direc_svdcoeff) ) + else: + print(' no modifications has been done in {}'.format(EZFIO_file)) + # ___________________________________________________________________________ + # ___________________________________________________________________________ + + print("end after {:.3f} minutes".format((time.time()-t0)/60.) ) + + # !!! +# !!! + +~ +~ + diff --git a/devel/svdwf/cr2_work/LOBD_postVMC.py b/devel/svdwf/cr2_work/LOBD_postVMC.py new file mode 100644 index 0000000..0d2f933 --- /dev/null +++ b/devel/svdwf/cr2_work/LOBD_postVMC.py @@ -0,0 +1,175 @@ +# !!! +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 + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +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) ) + SErr_H = np.zeros( (n_svd,n_svd) ) + + f = open(VMCmat_file, 'r') + Lines = f.readlines() + for line in Lines: + line = line.split() + ii = int ( line[0] ) + jj = int ( line[1] ) + Ci_h_matrix_svd[ii][jj] = float( line[2] ) + SErr_H [ii][jj] = float( line[3] ) + + return(Ci_h_matrix_svd,SErr_H) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +def get_Ci_overlap_matrix_svd(): + + Ci_overlap_matrix_svd = np.zeros( (n_svd,n_svd) ) + SErr_O = np.zeros( (n_svd,n_svd) ) + + f = open(VMCmat_file, 'r') + Lines = f.readlines() + for line in Lines: + line = line.split() + ii = int ( line[0] ) + jj = int ( line[1] ) + Ci_overlap_matrix_svd[ii][jj] = float( line[4] ) + SErr_O [ii][jj] = float( line[5] ) + + return(Ci_overlap_matrix_svd,SErr_O) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +def check_symmetric(a, tol=1e-3): + return np.all(np.abs(a-a.T) < tol) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +def get_Esvd(): + + # read CI_SVD matrices + Ci_h_matrix_svd , SErr_H = get_Ci_h_matrix_svd() + Ci_overlap_matrix_svd, SErr_O = 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] ) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + + +if __name__ == '__main__': + + t0 = time.time() + + # ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ # + #EZFIO_file = "Cr2_TSVD30" + #EZFIO_file = "Cr2_TSVD100" + #EZFIO_file = "Cr2_TSVD200" + #EZFIO_file = "Cr2_TSVD924" + #EZFIO_file = "Cr2_TSVD924_1" + + #EZFIO_file = "Cr2_TSVD30_it1" + #EZFIO_file = "Cr2_TSVD100_it1" + #EZFIO_file = "Cr2_TSVD200_it1" + EZFIO_file = "Cr2_TSVD924_1_it1" + + #EZFIO_file = "Cr2_TSVD200_minmaxit1" + + #EZFIO_file = "Cr2_TSVD30_it2" + + E_toadd = 64.8242130309350 + # ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ # + + ezfio.set_file(EZFIO_file) + n_svd = ezfio.get_spindeterminants_n_svd_coefs() + psi_svd_coeff = np.array(ezfio.get_spindeterminants_psi_svd_coefs()) + + print(" ") + print(" Today's date:", datetime.now() ) + print(" EZFIO file = {}".format(EZFIO_file)) + print(" n_svd = {}\n".format(n_svd) ) + + #VMCmat_list = [ 'SVD_QMCmat_FSVD_2b.txt' , 'SVD_QMCmat_FSVD_3b.txt' , 'SVD_QMCmat_FSVD_0.txt' ] + #VMCmat_list = [ 'SVD_QMCmat_FSVD_2b_1.txt', 'SVD_QMCmat_FSVD_3b_1.txt' + # , 'SVD_QMCmat_FSVD_4b_1.txt', 'SVD_QMCmat_FSVD_5b_1.txt' + # , 'SVD_QMCmat_FSVD_1.txt' ] + #for VMCmat_file in VMCmat_list: + # print(" SVD matrices file = {}".format(VMCmat_file)) + # E_svd, sigma_svd = get_Esvd() + # print(' svd enegry = {} \n'.format(E_svd) ) + + + VMCmat_file = 'SVD_QMCmat_FSVD_4b_1.txt' + print(" SVD matrices file = {}".format(VMCmat_file)) + E_svd, sigma_svd = get_Esvd() + sigma_svd = sigma_svd * np.sign(sigma_svd[0]) + print(' svd coeff = {}'.format(sigma_svd)) + + # ___________________________________________________________________________ + # ___________________________________________________________________________ + # + #save_to_EZFIO = 'svd_QMC' + save_to_EZFIO = '' + + 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) ) + else: + pass + # ___________________________________________________________________________ + # ___________________________________________________________________________ + + print("end after {:.3f} minutes".format((time.time()-t0)/60.) ) + + # !!! +# !!! + +~ diff --git a/devel/svdwf/cr2_work/LOBD_postcompactSVD.py b/devel/svdwf/cr2_work/LOBD_postcompactSVD.py new file mode 100644 index 0000000..a4c30b2 --- /dev/null +++ b/devel/svdwf/cr2_work/LOBD_postcompactSVD.py @@ -0,0 +1,173 @@ +# !!! +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 + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +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_modifsvd(): + + Ci_h_matrix_modifsvd = np.zeros( (n_modifsvd,n_modifsvd) ) + SErr_H = np.zeros( (n_modifsvd,n_modifsvd) ) + + f = open(compactVMCmat_file, 'r') + Lines = f.readlines() + for line in Lines: + line = line.split() + ii = int ( line[0] ) + jj = int ( line[1] ) + Ci_h_matrix_modifsvd[ii][jj] = float( line[2] ) + SErr_H [ii][jj] = float( line[3] ) + + return(Ci_h_matrix_modifsvd,SErr_H) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +def get_Ci_overlap_matrix_modifsvd(): + + Ci_overlap_matrix_modifsvd = np.zeros( (n_modifsvd,n_modifsvd) ) + SErr_O = np.zeros( (n_modifsvd,n_modifsvd) ) + + f = open(compactVMCmat_file, 'r') + Lines = f.readlines() + for line in Lines: + line = line.split() + ii = int ( line[0] ) + jj = int ( line[1] ) + Ci_overlap_matrix_modifsvd[ii][jj] = float( line[4] ) + SErr_O [ii][jj] = float( line[5] ) + + return(Ci_overlap_matrix_modifsvd,SErr_O) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +def check_symmetric(a, tol=1e-3): + return np.all(np.abs(a-a.T) < tol) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! +def get_Esvd(): + + Ci_h_matrix_modifsvd , SErr_H = get_Ci_h_matrix_modifsvd() + Ci_overlap_matrix_modifsvd, SErr_O = get_Ci_overlap_matrix_modifsvd() + + # symmetrise and diagonalise + aa = Ci_h_matrix_modifsvd + aa = 0.5*( aa + aa.T ) + bb = Ci_overlap_matrix_modifsvd + eigvals_modifsvd, vr = eig(aa, bb, left=False, right=True, overwrite_a=True, overwrite_b=True, + check_finite=True, homogeneous_eigvals=False) + + psi_tmp = psi_svd_coefs[0,0:n_modifsvd] + recouvre_modifsvd = np.abs(psi_tmp @ vr) + ind_gssvd = np.argmax(recouvre_modifsvd) + + E_modifsvd = eigvals_modifsvd[ind_gssvd] + E_toadd + + return( E_modifsvd, vr[:,ind_gssvd] ) +# ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! ~ ! + + + + +if __name__ == '__main__': + + t0 = time.time() + + # ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ # + #EZFIO_file = "Cr2_compactSVD30" + + #EZFIO_file = "Cr2_compactSVD30_it1" + EZFIO_file = "Cr2_compactSVD30_it1_1" + + E_toadd = 64.8242130309350 + # ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ # + + ezfio.set_file(EZFIO_file) + n_svd = ezfio.get_spindeterminants_n_svd_coefs() + n_modifsvd = ezfio.get_spindeterminants_n_modifsvd_coefs() + psi_svd_coefs = np.array(ezfio.get_spindeterminants_psi_svd_coefs()) + + print(" ") + print(" Today's date:", datetime.now() ) + print(" EZFIO file = {}".format(EZFIO_file)) + print(" n_svd = {}".format(n_svd) ) + print(" n_modifsvd = {}\n".format(n_modifsvd) ) + #print(' initial svd coeff = {}\n'.format(psi_svd_coefs)) + + compactVMCmat_list = [ 'compactSVD_QMCmat_1.txt' , 'compactSVD_QMCmat_2.txt' ] + for compactVMCmat_file in compactVMCmat_list: + print(" compactSVD matrices file = {}".format(compactVMCmat_file)) + E_compactsvd, sigma_compactsvd = get_Esvd() + print(' compact svd enegry = {}'.format(E_compactsvd) ) + + print('') + compactVMCmat_file = 'compactSVD_QMCmat_1.txt' + print(" compactSVD matrices file = {}".format(compactVMCmat_file)) + E_compactsvd, sigma_compactsvd = get_Esvd() + sigma_compactsvd = sigma_compactsvd * np.sign(sigma_compactsvd[0]) + print(' compact svd enegry = {}'.format(E_compactsvd) ) + print(' compact svd coeff = {}\n'.format(sigma_compactsvd)) + + print(' compact svd factor = {}\n'.format(sigma_compactsvd[n_modifsvd-1])) + + # ___________________________________________________________________________ + # ___________________________________________________________________________ + # + #save_to_EZFIO = 'compactsvd_QMC' + save_to_EZFIO = '' + + if( save_to_EZFIO == 'compactsvd_QMC' ): + ezfio_sigma_svd = np.zeros((n_svd)) + + for i in range(n_modifsvd-1): + ezfio_sigma_svd[i] = sigma_compactsvd[i] + for i in range(n_modifsvd-1, n_svd): + ezfio_sigma_svd[i] = psi_svd_coefs[0,i] * sigma_compactsvd[n_modifsvd-1] + + #print(' new svd coeff = {}\n'.format(ezfio_sigma_svd)) + + ezfio.set_spindeterminants_psi_svd_coefs( ezfio_sigma_svd ) + direc_svdcoeff = EZFIO_file + '/spindeterminants/psi_svd_coefs.gz' + print(' {} is modified'.format(direc_svdcoeff) ) + else: + print(' no modifications has been done in {}'.format(EZFIO_file)) + # ___________________________________________________________________________ + # ___________________________________________________________________________ + + print("end after {:.3f} minutes".format((time.time()-t0)/60.) ) + + # !!! +# !!! + +~ diff --git a/devel/svdwf/qmcchemRun.slurm b/devel/svdwf/qmcchemRun.slurm new file mode 100644 index 0000000..58e0ad6 --- /dev/null +++ b/devel/svdwf/qmcchemRun.slurm @@ -0,0 +1,31 @@ +#!/bin/bash +# +#SBATCH --job-name=QMC_n20_2h2 +#SBATCH --output=%x_%j.out +# +#SBATCH --mail-type=ALL +#SBATCH --mail-user=aammar@irsamc.ups-tlse.fr +# +#SBATCH -N 1 +#SBATCH -n 32 +#SBATCH -c 1 +#SBATCH -t 24:00:00 +#SBATCH -p xeonv4 +#SBATCH --exclusive +# +# +# +srun hostname +echo "Running on `hostname`" +# +echo "Kicking off $SLURM_NPROCS compute servers on $SLURM_NNODES nodes" +echo "Nodes are\n `scontrol show hostname`" +# +source ~/qmcchem/qmcchemrc +# +##qmcchem edit -c -j None -l 10 -t 120 h2o_631g >> qmcchem_h2o_631g.out +#qmcchem edit -c -j None -l 600 -t 10000 h2o_631g >> qmcchem_h2o_631g.out +#qmcchem run h2o_631g_frez_Jopt_nsvd10_iter2postsvd >> qmcchem_h2o_631g_frez_Jopt_nsvd10_iter2postsvd.out +# +qmcchem run 2h2_work/2h2_cisd_nsvd20 >> QMC_2h2_nsvd20.out +# diff --git a/devel/svdwf/qpRun.slurm b/devel/svdwf/qpRun.slurm new file mode 100644 index 0000000..f627081 --- /dev/null +++ b/devel/svdwf/qpRun.slurm @@ -0,0 +1,35 @@ +#!/bin/bash +# +##SBATCH --job-name=buildpsi_h2o_pv +#SBATCH --job-name=diag_f2 +#SBATCH --output=%x_%j.out +#SBATCH --mail-type=ALL +#SBATCH --mail-user=aammar@irsamc.ups-tlse.fr +# +#SBATCH -N 1 +#SBATCH -n 1 +#SBATCH -c 32 +##SBATCH -c 1 +#SBATCH -t 72:00:00 +#SBATCH -p xeonv2 +#SBATCH --exclusive +# +srun hostname +echo "Running on `hostname`" +echo "Kicking off $SLURM_NPROCS compute servers on $SLURM_NNODES nodes" +echo "Nodes are\n `scontrol show hostname`" +# +#qp create_ezfio -b 6-31g h2o.xyz -o h2o_631g +#qp create_ezfio -p bfd -b cc-pvdz_ecp_bfd h2o.xyz -o h2o_pseudo +# +~/qp2/bin/qpsh +# +#qp_run scf h2o_631g >> scf_h2o_631g.out +#qp_run cisd h2o_631g >> cisd_h2o_631g.out +# +#qp_run buildpsi_postsvd_pv h2o_work/h2o_631g_frez_cisd >> buildpsi_postsvd_h2o_631g_frez_cisd.dat +# +#qp_run diag_H_npostsvd h2o_work/h2o_631g_frez_fci >> diag_H_npostsvd.dat +qp_run diag_H_npostsvd_pv f2_work/f2_cisd >> diag_H_npostsvd_pv_f2.dat +#qp_run diag_H_npostsvd_pv h2o_work/h2o_631g_frez_cisd >> diag_H_npostsvd_pv_h2o.dat +# diff --git a/devel/svdwf/run_py.slurm b/devel/svdwf/run_py.slurm new file mode 100644 index 0000000..6157724 --- /dev/null +++ b/devel/svdwf/run_py.slurm @@ -0,0 +1,26 @@ +#!/bin/bash +#SBATCH --job-name=RSVD +#SBATCH --output=%x_%j.out +# +#SBATCH --mail-type=ALL +#SBATCH --mail-user=aammar@irsamc.ups-tlse.fr +# +#SBATCH -N 1 +#SBATCH -n 1 +#SBATCH -c 1 +#SBATCH -t 250:00:00 +#SBATCH -p xeonv2 +###SBATCH --exclusive +# +# +# +srun hostname +echo "Running on `hostname`" +# +echo "Kicking off $SLURM_NPROCS compute servers on $SLURM_NNODES nodes" +echo "Nodes are\n `scontrol show hostname`" +# +source ~/qmcchem/qmcchemrc +# +python perform_RSVD.py >> RSVD_py.out +