From 7a6450b6fa9038fa88e326d252fd2fde57dfd0ef Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Malte=20Sch=C3=BCler?= Date: Mon, 16 Sep 2019 10:59:46 +0200 Subject: [PATCH] implemented multiple ncsf VASP cycles --- doc/tutorials/images_scripts/nio_csc.py | 2 +- doc/tutorials/nio_csc.rst | 36 ++++++++------- python/converters/plovasp/elstruct.py | 8 +++- python/converters/plovasp/sc_dmft.py | 60 ++++++++++++++++++++----- shells/vasp_dmft.bash.in | 40 ++++++++++++++--- 5 files changed, 111 insertions(+), 35 deletions(-) diff --git a/doc/tutorials/images_scripts/nio_csc.py b/doc/tutorials/images_scripts/nio_csc.py index 993293da..75f2c2fa 100644 --- a/doc/tutorials/images_scripts/nio_csc.py +++ b/doc/tutorials/images_scripts/nio_csc.py @@ -105,7 +105,7 @@ def dmft_cycle(): ar['DMFT_input']['code_versions']["cthyb_version"] = cthyb_version.version ar['DMFT_input']['code_versions']["cthyb_git"] = cthyb_version.cthyb_hash ar['DMFT_input']['code_versions']["dft_tools_version"] = dft_tools_version.version - ar['DMFT_input']['code_versions']["dft_tools_version"] = dft_tools_version.dft_tools_hash + ar['DMFT_input']['code_versions']["dft_tools_git"] = dft_tools_version.dft_tools_hash if 'iteration_count' in ar['DMFT_results']: iteration_offset = ar['DMFT_results']['iteration_count']+1 S.Sigma_iw = ar['DMFT_results']['Iterations']['Sigma_it'+str(iteration_offset-1)] diff --git a/doc/tutorials/nio_csc.rst b/doc/tutorials/nio_csc.rst index 23f60f43..ab12cbdf 100644 --- a/doc/tutorials/nio_csc.rst +++ b/doc/tutorials/nio_csc.rst @@ -53,7 +53,7 @@ DMFT dmft script ------------------------------ -Since the python script for performing the dmft loop pretty much resembles that presented in the tutorial on :ref:`srvo3`, we will not go into detail here but simply provide the script :ref:`nio.py`. Following Kunes et al. in `PRB 75 165115 (2007) `_ we use :math:`U=8` and :math:`J=1`. We se;ect :math:`\beta=5` instead of :math:`\beta=10` to ease the problem slightly. For simplicity we fix the double-counting potential to :math:`\mu_{DC}=59` eV by:: +Since the python script for performing the dmft loop pretty much resembles that presented in the tutorial on :ref:`SrVO3 `, we will not go into detail here but simply provide the script :ref:`nio.py`. Following Kunes et al. in `PRB 75 165115 (2007) `_ we use :math:`U=8` and :math:`J=1`. We se;ect :math:`\beta=5` instead of :math:`\beta=10` to ease the problem slightly. For simplicity we fix the double-counting potential to :math:`\mu_{DC}=59` eV by:: DC_value = 59.0 SK.calc_dc(dm, U_interact=U, J_hund=J, orb=0, use_dc_value=DC_value) @@ -62,11 +62,11 @@ For sensible results run this script in parallel on at least 20 cores. As a quic Local lattice Green's function for all projected orbitals ----------------------- +--------------------------------------------------------- We calculate the local lattice Green's function - now also for the uncorrelated orbitals, i.e., the O p states, for what we use the script :ref:`NiO_local_lattice_GF.py`. The result is saved in the h5 file as `G_latt_orb_it`, where `` is the number of the last DMFT iteration. Spectral function on real axis: MaxEnt ----------------------- +-------------------------------------- To compare to results from literature we make use of the `maxent triqs application `_ and calculate the spectral function on real axis. Use this script to perform a crude but quick calculation: :ref:`maxent.py` using a linear real axis and a line-fit analyzer to determine the optimal :math:`\alpha`. The output is saved in the h5 file in `DMFT_results/Iterations/G_latt_orb_w_o_it`, where `` is the number of the orbital and `n_it` is again the number of the last iteration. The real axis information is stored in `DMFT_results/Iterations/w_it`. @@ -74,26 +74,28 @@ To compare to results from literature we make use of the `maxent triqs applicati :width: 400 :align: center -.. Charge self-consistent DMFT -.. ================================================== +Charge self-consistent DMFT +================================================== -.. In this part we will perform charge self-consistent DMFT calculations. To do so we have to adapt the VASP `INCAR` such that :program:`VASP` reads the updated charge density after each step. We add the lines +In this part we will perform charge self-consistent DMFT calculations. To do so we have to adapt the VASP `INCAR` such that :program:`VASP` reads the updated charge density after each step. We add the lines:: -.. ICHARG = 5 -.. NELM = 1000 -.. NELMIN = 1000 + ICHARG = 5 + NELM = 1000 + NELMIN = 1000 -.. which makes VASP wait after each step of its iterative diagonalization until the file vasp.lock is created. It then reads the update of the charge density in the file `GAMMA`. It is terminated by an external script after a desired amount of steps, such that we deactivate all automatic stoping criterion by setting the number of steps to a very high number. +which makes VASP wait after each step of its iterative diagonalization until the file vasp.lock is created. It then reads the update of the charge density in the file `GAMMA`. It is terminated by an external script after a desired amount of steps, such that we deactivate all automatic stoping criterion by setting the number of steps to a very high number. -.. We take the respective converged DFT and DMFT calculations from before as a starting point. I.e., we copy the `CHGCAR` and `nio.h5` together with the other :program:`VASP` input files and :file:`plo.cfg` in a new directory. We use a script called :program:`vasp_dmft` to invoke :program:`VASP` in the background and start the DMFT calculation together with :program:`plovasp` and the converter. This script assumes that the dmft sript contains a function `dmft_cycle()` and also the conversion from text files to the h5 file. Additionally we have to add a few lines to calculate the density correction and calculate the correlation energy. We adapt the script straightforardly and call it :ref:`nio_csc.py`. The most important new lines are +We take the respective converged DFT and DMFT calculations from before as a starting point. I.e., we copy the `CHGCAR` and `nio.h5` together with the other :program:`VASP` input files and :file:`plo.cfg` in a new directory. We use a script called :program:`vasp_dmft` to invoke :program:`VASP` in the background and start the DMFT calculation together with :program:`plovasp` and the converter. This script assumes that the dmft sript contains a function `dmft_cycle()` and also the conversion from text files to the h5 file. Additionally we have to add a few lines to calculate the density correction and calculate the correlation energy. We adapt the script straightforardly (for a working example see :ref:`nio_csc.py`). The most important new lines are:: -.. SK.chemical_potential = SK.calc_mu( precision = 0.000001 ) -.. SK.calc_density_correction(dm_type='vasp') -.. correnerg = 0.5 * (S.G_iw * S.Sigma_iw).total_density() + SK.chemical_potential = SK.calc_mu( precision = 0.000001 ) + SK.calc_density_correction(dm_type='vasp') + correnerg = 0.5 * (S.G_iw * S.Sigma_iw).total_density() -.. where the chemical potential is determined to a greater precision than before, the correction to the dft density matrix is calculated and output to the file :file:`GAMMA`. The correlation energy is calculated via Migdal-Galitzki formula. We also slightly increase the tolerance for the detection of blocks since the DFT calculation now includes some QMC noise. +where the chemical potential is determined to a greater precision than before, the correction to the dft density matrix is calculated and output to the file :file:`GAMMA`. The correlation energy is calculated via Migdal-Galitzki formula. We also slightly increase the tolerance for the detection of blocks since the DFT calculation now includes some QMC noise. -.. We can start the whole machinery by excectuing +To help convergence, we keep the density (i.e., the GAMMA file) fixed for a few DFT iterations. We do so since VASP uses an iterative diagonalization. -.. vasp_dmft -n -i -p nio_csc.py +We can start the whole machinery by excecuting:: + + vasp_dmft -n -i -j -p nio_csc.py diff --git a/python/converters/plovasp/elstruct.py b/python/converters/plovasp/elstruct.py index c30b88e3..9421e768 100644 --- a/python/converters/plovasp/elstruct.py +++ b/python/converters/plovasp/elstruct.py @@ -140,6 +140,7 @@ class ElectronicStructure: ions = sorted(list(set([param['isite'] for param in self.proj_params]))) nions = len(ions) norb = nproj / nions + print nproj, ns, nk, nb, ions, nions, norb # Spin factor sp_fac = 2.0 if ns == 1 and not self.nc_flag else 1.0 @@ -162,12 +163,15 @@ class ElectronicStructure: for ispin in xrange(ns): print print " Spin:", ispin + 1 + print self.proj_params for io, ion in enumerate(ions): print " Site:", ion iorb_inds = [(ip, param['m']) for ip, param in enumerate(self.proj_params) if param['isite'] == ion] + print iorb_inds norb = len(iorb_inds) - dm = np.zeros((norb, norb)) - ov = np.zeros((norb, norb)) + norb2 = self.proj_params[0]['l']*2+1 + dm = np.zeros((norb2, norb2)) + ov = np.zeros((norb2, norb2)) for ind, iorb in iorb_inds: for ind2, iorb2 in iorb_inds: dm[iorb, iorb2] = den_mat[ispin, ind, ind2] diff --git a/python/converters/plovasp/sc_dmft.py b/python/converters/plovasp/sc_dmft.py index 7ae107fc..11fdfa2a 100644 --- a/python/converters/plovasp/sc_dmft.py +++ b/python/converters/plovasp/sc_dmft.py @@ -32,6 +32,7 @@ import signal import sys import pytriqs.utility.mpi as mpi import converter +from shutil import copyfile xch = sys.excepthook def excepthook(typ, value, traceback): @@ -40,7 +41,7 @@ def excepthook(typ, value, traceback): mpi.MPI.COMM_WORLD.Abort(1) sys.excepthook = excepthook -debug = True +debug = False # # Helper functions # @@ -99,7 +100,7 @@ class bcolors: # # Main self-consistent cycle # -def run_all(vasp_pid, dmft_cycle, cfg_file, n_iter): +def run_all(vasp_pid, dmft_cycle, cfg_file, n_iter, n_iter_dft, vasp_version): """ """ mpi.report(" Waiting for VASP lock to appear...") @@ -120,6 +121,7 @@ def run_all(vasp_pid, dmft_cycle, cfg_file, n_iter): mpi.report(" VASP stopped") vasp_running = False break + # Tell VASP to stop if the maximum number of iterations is reached iter += 1 if iter == n_iter: @@ -154,10 +156,31 @@ def run_all(vasp_pid, dmft_cycle, cfg_file, n_iter): print " DFT DC: ", dft_dc print "="*80 print - - if mpi.is_master_node() and vasp_running: - open('./vasp.lock', 'a').close() - + + # check if we should do additional VASP calculations + # in the standard VASP version, VASP writes out GAMMA itself + # so that if we want to keep GAMMA fixed we have to copy it to + # GAMMA_recent and copy it back after VASP has completed an iteration + # if we are using a hacked Version of VASP where the write out is + # disabled we can skip this step. + # the hack consists of removing the call of LPRJ_LDApU in VASP src file + # electron.F around line 644 + iter_dft = 0 + if vasp_version == 'standard': + copyfile(src='GAMMA',dst='GAMMA_recent') + while iter_dft < n_iter_dft: + if mpi.is_master_node(): + open('./vasp.lock', 'a').close() + while is_vasp_lock_present(): + time.sleep(1) + if not is_vasp_running(vasp_pid): + mpi.report(" VASP stopped") + vasp_running = False + break + iter_dft += 1 + if vasp_version == 'standard': + copyfile(src='GAMMA_recent',dst='GAMMA') + if mpi.is_master_node(): total_energy = dft_energy + corr_energy - dft_dc with open('TOTENERGY', 'w') as f: @@ -170,7 +193,9 @@ def run_all(vasp_pid, dmft_cycle, cfg_file, n_iter): mpi.report("***Done") def main(): + import importlib + try: vasp_pid = int(sys.argv[1]) except (ValueError, KeyError): @@ -184,19 +209,34 @@ def main(): if mpi.is_master_node(): print "ERROR: Number of iterations must be provided as the second argument" raise + + try: + n_iter_dft = int(sys.argv[3]) + except (ValueError, KeyError): + if mpi.is_master_node(): + print "ERROR: Number of VASP iterations with fixed charge density must be provided as the third argument" + raise try: - dmft_script = re.sub("\.py$", "", sys.argv[3]) + dmft_script = re.sub("\.py$", "", sys.argv[4]) except: if mpi.is_master_node(): - print "ERROR: User-defined DMFT script must be provided as the third argument" + print "ERROR: User-defined DMFT script must be provided as the fourth argument" raise # Optional parameter: config-file name try: - cfg_file = sys.argv[4] + cfg_file = sys.argv[5] except KeyError: cfg_file = 'plo.cfg' + + try: + vasp_version = sys.argv[6] + except KeyError: + vasp_version = 'standard' + + if vasp_version != 'standard' and vasp_version != 'no_gamma_write': + raise Exception('vasp_version has to be standard or no_gamma_write') # if len(sys.argv) > 1: # vasp_path = sys.argv[1] @@ -211,7 +251,7 @@ def main(): dmft_mod = importlib.import_module(dmft_script) - run_all(vasp_pid, dmft_mod.dmft_cycle, cfg_file, n_iter) + run_all(vasp_pid, dmft_mod.dmft_cycle, cfg_file, n_iter, n_iter_dft, vasp_version) if __name__ == '__main__': main() diff --git a/shells/vasp_dmft.bash.in b/shells/vasp_dmft.bash.in index 838b2e39..884e2ec6 100755 --- a/shells/vasp_dmft.bash.in +++ b/shells/vasp_dmft.bash.in @@ -5,9 +5,9 @@ MPIRUN_CMD=mpirun show_help() { echo " -Usage: vasp_dmft [-n ] -i [-p ] [] +Usage: vasp_dmft [-n ] -i -j [-v ] [-p ] [] - If the number of cores is note specified it is set to 1 by default. + If the number of cores is not specified it is set to 1 by default. must provide an importable function 'dmft_cycle()' which is invoked once per DFT+DMFT iteration. If the script name is omitted the default name @@ -18,7 +18,7 @@ Usage: vasp_dmft [-n ] -i [-p &2 show_help @@ -55,6 +71,7 @@ if [ -z "$NITER" ]; then exit 1 fi + if [ -z "$VASP_DIR" ]; then echo " Error: A path to VASP directory must be given either with option -p or by setting variable VASP_DIR" >&2 show_help @@ -66,6 +83,16 @@ if [ -z "$NPROC" ]; then NPROC=1 fi +if [ -z "$NDFTITER" ]; then + echo " Number of VASP iterations without updating density not specified, setting to 1" + NDFTITER=1 +fi + +if [ -z "$VASP_VERSION" ]; then + echo " VASP version not specified, setting to standard" + VASP_VERSION="standard" +fi + shift $((OPTIND-1)) if [ -z "$1" ]; then @@ -76,11 +103,14 @@ fi echo " Number of cores: $NPROC" echo " Number of iterations: $NITER" +echo " Number of iterations with fixed density: $NDFTITER" +echo " VASP version: $VASP_VERSION" echo " Script name: $DMFT_SCRIPT" rm -f vasp.lock -stdbuf -o 0 $MPIRUN_CMD -np $NPROC "$VASP_DIR" & +# stdbuf -o 0 $MPIRUN_CMD -np $NPROC "$VASP_DIR" & +vasp_intel.sh & -$MPIRUN_CMD -np $NPROC python -m triqs_dft_tools.converters.plovasp.sc_dmft $(jobs -p) $NITER $DMFT_SCRIPT 'plo.cfg' || kill %1 +$MPIRUN_CMD -np $NPROC python -m triqs_dft_tools.converters.plovasp.sc_dmft $(jobs -p) $NITER $NDFTITER $DMFT_SCRIPT 'plo.cfg' $VASP_VERSION || kill %1