x init_lapw
Enter reduction in %: 0 Use old or new scheme (o/N) : N Do you want to accept these radii; .... (a/d/r) : a specify nn-bondlength factor: 2 Ctrl-X-C continue with sgroup : c Ctrl-X-C continue with symmetry : c Ctrl-X-C continue with lstart : c specify switches for instgen_lapw (or press ENTER): ENTER SELECT XCPOT: 5 SELECT ENERGY : -6 Ctrl-X-C continue with kgen : c Ctrl-X-C Ctrl-X-C NUMBER OF K-POINTS : 500 Shift of k-mesh allowed. Do you want to shift: (0=no, 1=shift): 1 Ctrl-X-C continue with dstart or execute kgen again or exit (c/e/x) : c Ctrl-X-C do you want to perform a spinpolarized calculation : n
run_lapw
x kgen
init_dmft.py
There are 4 atoms in the unit cell: 1 Fe 2+ 2 Fe 2+ 3 Se 2- 4 Se 2- Specify correlated atoms (ex: 1-4,7,8): 1-2 You have chosen the following atoms to be correlated: 1 Fe 2+ 2 Fe 2+ Do you want to continue; or edit again? (c/e): c For each atom, specify correlated orbital(s) (ex: d,f): 1 Fe 2+: d 2 Fe 2+: d You have chosen to apply correlations to the following orbitals: 1 Fe 2+-1 d 2 Fe 2+-2 d Do you want to continue; or edit again? (c/e): c Specify qsplit for each correlated orbital (default = 0): Qsplit Description ------ ------------------------------------------------------------ 2 real harmonics basis, no symmetry, except spin (up=dn) ------ ------------------------------------------------------------ 1 Fe 2+-1 d: 2 2 Fe 2+-2 d: 2 Do you want to continue; or edit again? (c/e): c Specify projector type (default = 2): ------ ------------------------------------------------------------ Projector Drscription 5 fixed projector, which is written to projectorw.dat. You can generate projectorw.dat with the tool wavef.py ------ ------------------------------------------------------------ > 5 Do you want to continue; or edit again? (c/e): c Do you want to group any of these orbitals into cluster-DMFT problems? (y/n): n Enter the correlated problems forming each unique correlated problem, separated by spaces (ex: 1,3 2,4 5-8): 1-2 Each set of equivalent correlated problems are listed below: 1 (Fe 2+1 d) (Fe 2+2 d) are equivalent. Do you want to continue; or edit again? (c/e): c Range (in eV) of hybridization taken into account in impurity problems; default -10.0, 10.0: <ENTER> Perform calculation on real; or imaginary axis? (r/i): i Ctrl-X-C Is this a spin-orbit run? (y/n): n Ctrl-X-C
19 44 1 5 # hybridization band index nemin and nemax, renormalize for interstitials, projection type 1 0.025 0.025 200 -3.000000 1.000000 # matsubara, broadening-corr, broadening-noncorr, nomega, omega_min, omega_max (in eV) 2 # number of correlated atoms 1 1 0 # iatom, nL, locrot 2 2 1 # L, qsplit, cix 2 1 0 # iatom, nL, locrot 2 2 2 # L, qsplit, cix
#================ # Siginds and crystal-field transformations for correlated orbitals ================ 2 5 5 # Number of independent kcix blocks, max dimension, max num-independent-components 1 5 5 # cix-num, dimension, num-independent-components #---------------- # Independent components are -------------- 'z^2' 'x^2-y^2' 'xz' 'yz' 'xy' #---------------- # Sigind follows -------------------------- 1 0 0 0 0 0 2 0 0 0 0 0 3 0 0 0 0 0 4 0 0 0 0 0 5 #---------------- # Transformation matrix follows ----------- 0.00000000 0.00000000 0.00000000 0.00000000 1.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.70710679 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.70710679 0.00000000 0.00000000 0.00000000 0.70710679 0.00000000 0.00000000 0.00000000 -0.70710679 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.70710679 0.00000000 0.00000000 0.00000000 0.70710679 0.00000000 0.00000000 -0.70710679 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.00000000 0.70710679 0.00000000
1 # number of sigind blocks 5 # dimension of this sigind block 1 0 0 0 0 0 2 0 0 0 0 0 3 0 0 0 0 0 4 0 0 0 0 0 5
dmft_copy.py <dft_results_directory>
solver = 'CTQMC' # impurity solver DCs = 'nominal' # double counting scheme max_dmft_iterations = 1 # number of iteration of the dmft-loop only max_lda_iterations = 100 # number of iteration of the LDA-loop only finish = 50 # number of iterations of full charge loop (1 = no charge self-consistency) ntail = 300 # on imaginary axis, number of points in the tail of the logarithmic mesh cc = 2e-6 # the charge density precision to stop the LDA+DMFT run ec = 2e-6 # the energy precision to stop the LDA+DMFT run recomputeEF = 1 # Recompute EF in dmft2 step. If recomputeEF = 2, it tries to find an insulating gap. wbroad = 0.0 # broadening of sigma on the imaginary axis kbroad = 0.0 # broadening of sigma on the imaginary axis # Impurity problem number 0 iparams0={"exe" : ["ctqmc" , "# Name of the executable"], "U" : [5.0 , "# Coulomb repulsion (F0)"], "J" : [0.8 , "# Coulomb repulsion (F0)"], "CoulombF" : ["'Ising'" , "# Can be set to 'Full'"], "beta" : [50 , "# Inverse temperature"], "svd_lmax" : [25 , "# We will use SVD basis to expand G, with this cutoff"], "M" : [5e6 , "# Total number of Monte Carlo steps"], "mode" : ["SH" , "# We will use self-energy sampling, and Hubbard I tail"], "nom" : [200 , "# Number of Matsubara frequency points sampled"], "tsample" : [400 , "# How often to record measurements"], "GlobalFlip" : [1000000 , "# How often to try a global flip"], "warmup" : [3e5 , "# Warmup number of QMC steps"], "nf0" : [6.0 , "# Nominal occupancy nd for double-counting"], }
solver = 'CTQMC'
DCs = 'nominal' # nominal double-counting
max_dmft_iterations = 1 max_lda_iterations = 100 finish = 50
ntail = 300
cc = 2e-6 # the charge density precision to stop the LDA+DMFT run ec = 2e-6 # the energy precision to stop the LDA+DMFT run
recomputeEF = 1
wbroad = 0.0 # broadening of sigma on the imaginary axis kbroad = 0.0 # broadening of sigma on the imaginary axis
"iparams0=...."
"U" : [5.0]
"J" : [0.8]
"CoulombF" : ["'Ising'"]
"beta" : [50. ]
"svd_lmax" : [25 ],
"M" : [10e6]
"mode" : ["SH" ]
"nom" : [200]
"nf0" : [6.0 ]
"Ncout" : [1000000]
"GlobalFlip" : [1000000 , "# How often to try a global flip"] "warmup" : [3e5 , "# Warmup number of QMC steps"], "tsample" : [400 , "# How often to record measurements"],
szero.py
run_dmft.py
but the precise form of the command depends on the supercomputer and its environment.
If many more CPU-cores are available than the number of k-points, we can
optimize the execution further. Namely, a single k-point can be
executed on many cores using open_mp instructions. To use this
feature, one needs to specify "mpi_prefix.dat2" in addition to
"mpi_prefix.dat" file. The former is than used in the dmft1, and dmft2 part of the
loop, while the latter is used by the impurity solver.
Clearly, "mpi_prefix.dat2" should specify the mpi command, which is started
on a subset of machines. Therefore "-np" should be smaller than
"$NSLOTS" and "OMP_NUM_THREADS" should be larger than 1.
We provide a script that attempts to do that, hence you might want to
insert the following line between "....> mpi_prefix.dat" and "run_dmft.py":
$WIEN_DMFT_ROOT/createDMFprefix.py $JOBNAME.klist mpi_prefix.dat > mpi_prefix.dat2
#!/bin/bash ######################################################################## # SUN Grid Engine job wrapper ######################################################################## #$ -N FeSe #$ -pe mpi2 250 #$ -q <your_queue> #$ -j y #$ -M <email@your_institution> #$ -m e ######################################################################## source $TMPDIR/sge_init.sh ######################################################################## export WIENROOT=<your_w2k_root> export WIEN_DMFT_ROOT=<your_dmft_root> export PATH=$PATH:$WIENROOT:$WIEN_DMFT_ROOT export PYTHONPATH=$PYTHONPATH:.:$WIEN_DMFT_ROOT export SCRATCH='.' JOBNAME="FeSe" # must be equal to case name mpi_prefix="mpiexec -np $NSLOTS -machinefile $TMPDIR/machines -port $port -env OMP_NUM_THREADS 1 -envlist SCRATCH,WIEN_DMFT_ROOT,PYTHONPATH,WIENROOT" echo $mpi_prefix > mpi_prefix.dat $WIEN_DMFT_ROOT/createDMFprefix.py $JOBNAME.klist mpi_prefix.dat > mpi_prefix.dat2 $WIEN_DMFT_ROOT/run_dmft.py > nohup.dat 2>&1 rm -f $JOBNAME.vector* find . -size 0|xargs rm
# FeSe.in2, line 1 : FOR (TOT,FOR,QTL,EFG,FERMI)
TOTAL FORCE IN mRy/a.u. = |F| Fx Fy Fz with/without FOR in case.in2 :FOR001: 1.ATOM 0.000000 0.000000 0.000000 0.000000 total forces :FOR002: 2.ATOM 19.938644 0.000000 0.000000 -19.938644 total forces TOTAL FORCE WITH RESPECT TO GLOBAL CARTESIAN COORDINATES: :FCA001: 1.ATOM 0.000000 0.000000 0.000000 total forces :FCA002: 2.ATOM 0.000000 0.000000 -19.938644 total forces TOTAL FORCE WITH RESPECT TO THE GLOBAL COORDINATE SYSTEM: :FGL001: 1.ATOM 0.000000000 0.000000000 0.000000000 total forces :FGL002: 2.ATOM 0.000000000 0.000000000 -19.938644357 total forces
grep ':FGL002' FeSe.scf
# FeSe.inm, line 1: MSR1a 0.0 YES (BROYD/PRATT, extra charge (+1 for additional e), norm)
:APOS001 X=0.25000000 Y=0.75000000 Z=0.00000000 : DX= 0.00 DY= 0.00 DZ= 0.00 TOT= 0.00 mau :APOS002 X=0.75000000 Y=0.75000000 Z=0.72948443 : DX= 0.00 DY= 0.00 DZ= -21.28 TOT= 21.28 mau :FRMSA: (mRyd/au) 73.175 103.485 RMS (au) 1.50E-02 MAX 2.13E-02 :F-condition (mRyd/au) 2.000 F
force_Rho_diff = 1e-4
grep ':FRM' FeSe.scf
:FRMSA: (mRyd/au) 658.782 931.658 RMS (au) 2.47E-03 MAX 3.50E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 622.141 879.841 RMS (au) 2.41E-03 MAX 3.41E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 599.243 847.458 RMS (au) 5.76E-02 MAX 8.15E-02 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 73.175 103.485 RMS (au) 1.50E-02 MAX 2.13E-02 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 22.986 32.507 RMS (au) 3.19E-03 MAX 4.50E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 14.839 20.985 RMS (au) 2.63E-03 MAX 3.72E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 4.120 5.826 RMS (au) 2.07E-03 MAX 2.92E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 1.355 1.916 RMS (au) 9.79E-03 MAX 1.38E-02 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 9.552 13.509 RMS (au) 1.60E-02 MAX 2.26E-02 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 5.092 7.201 RMS (au) 1.24E-02 MAX 1.75E-02 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 7.584 10.725 RMS (au) 2.91E-03 MAX 4.11E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 3.919 5.542 RMS (au) 4.41E-03 MAX 6.23E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 0.782 1.105 RMS (au) 4.38E-03 MAX 6.20E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 0.759 1.073 RMS (au) 6.89E-03 MAX 9.75E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 0.727 1.028 RMS (au) 6.85E-03 MAX 9.68E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 0.089 0.126 RMS (au) 5.45E-03 MAX 7.71E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 0.649 0.918 RMS (au) 1.87E-03 MAX 2.65E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 1.281 1.811 RMS (au) 2.39E-03 MAX 3.38E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 0.475 0.671 RMS (au) 2.51E-03 MAX 3.55E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 0.991 1.402 RMS (au) 1.37E-03 MAX 1.94E-03 :F-condition (mRyd/au) 2.000 T :FRMSA: (mRyd/au) 0.056 0.079 RMS (au) 1.14E-04 MAX 1.61E-04 :F-condition (mRyd/au) 2.000 T :FRMSA: (mRyd/au) 0.046 0.066 RMS (au) 9.53E-05 MAX 1.35E-04 :F-condition (mRyd/au) 2.000 T :FRMSA: (mRyd/au) 0.005 0.007 RMS (au) 1.51E-05 MAX 2.14E-05 :F-condition (mRyd/au) 2.000 T :FRMSA: (mRyd/au) 0.425 0.600 RMS (au) 1.06E-05 MAX 1.50E-05 :F-condition (mRyd/au) 2.000 T :FRMSA: (mRyd/au) 0.153 0.216 RMS (au) 1.14E-05 MAX 1.61E-05 :F-condition (mRyd/au) 2.000 T :FRMSA: (mRyd/au) 3.312 4.684 RMS (au) 2.50E-04 MAX 3.53E-04 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 13.764 19.465 RMS (au) 1.90E-04 MAX 2.68E-04 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 7.423 10.498 RMS (au) 2.81E-03 MAX 3.98E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 0.956 1.352 RMS (au) 3.48E-05 MAX 4.92E-05 :F-condition (mRyd/au) 2.000 T :FRMSA: (mRyd/au) 1.352 1.912 RMS (au) 2.41E-03 MAX 3.41E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 1.693 2.394 RMS (au) 3.51E-03 MAX 4.96E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 2.032 2.874 RMS (au) 4.66E-03 MAX 6.60E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 0.223 0.315 RMS (au) 7.32E-03 MAX 1.04E-02 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 1.224 1.731 RMS (au) 8.78E-03 MAX 1.24E-02 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 0.098 0.139 RMS (au) 2.79E-03 MAX 3.94E-03 :F-condition (mRyd/au) 2.000 F :FRMSA: (mRyd/au) 0.804 1.136 RMS (au) 7.93E-04 MAX 1.12E-03 :F-condition (mRyd/au) 2.000 T :FRMSA: (mRyd/au) 0.001 0.001 RMS (au) 1.27E-04 MAX 1.79E-04 :F-condition (mRyd/au) 2.000 T :FRMSA: (mRyd/au) 0.149 0.211 RMS (au) 1.80E-04 MAX 2.55E-04 :F-condition (mRyd/au) 2.000 T :FRMSA: (mRyd/au) 0.001 0.001 RMS (au) 4.60E-05 MAX 6.50E-05 :F-condition (mRyd/au) 2.000 T :FRMSA: (mRyd/au) 0.010 0.014 RMS (au) 1.59E-05 MAX 2.25E-05 :F-condition (mRyd/au) 2.000 T
grep ':APOS002' FeSe.scf
grep ':APOS002' FeSe.scf| sed 's/Z=/ /'|strans '$j++' '1-#4'
saverage.py saverage.py sig.inp.2[5-9].1
params={'statistics': 'fermi', # fermi/bose 'Ntau' : 300, # Number of time points 'L' : 20.0, # cutoff frequency on real axis 'x0' : 0.005, # low energy cut-off 'bwdth' : 0.004, # smoothing width 'Nw' : 300, # number of frequency points on real axis 'gwidth' : 2*15.0, # width of gaussian 'idg' : 1, # error scheme: idg=1 -> sigma=deltag ; idg=0 -> sigma=deltag*G(tau) 'deltag' : 0.005, # error 'Asteps' : 4000, # anealing steps 'alpha0' : 1000, # starting alpha 'min_ratio' : 0.001, # condition to finish, what should be the ratio 'iflat' : 1, # iflat=0 : constant model, iflat=1 : gaussian of width gwidth, iflat=2 : input using file model$ 'Nitt' : 1000, # maximum number of outside iterations 'Nr' : 0, # number of smoothing runs 'Nf' : 40, # to perform inverse Fourier, high frequency limit is computed from the last Nf points }
maxent_run.py sig.inpx
dmft_copy.py ../
cp ../maxent/Sig.out sig.inp
19 42 1 5 # hybridization band index nemin and nemax, renormalize for interstitials, projection type 0 0.025 0.025 200 -3.000000 1.000000 # matsubara, broadening-corr, broadening-noncorr, nomega, omega_min, omega_max (in eV) 2 # number of correlated atoms
x lapw0 -f FeSe x_dmft.py lapw1 x_dmft.py dmft1
19 42 1 5 # hybridization band index nemin and nemax, renormalize for interstitials, projection type 0 0.025 0.025 200 -1.000000 1.000000 # matsubara, broadening-corr, broadening-noncorr, nomega, omega_min, omega_max (in eV)
x_dmft.py lapw1 --band x_dmft.py dmftp
wakplot.py 1.0
To proceed, we will create a new directory, and while in that directory, copy results from previous dmft run into this new directory, i.e.,
dmft_copy.py <dmft_previous_run>
"CoulombF" : ["'Full'" , "# Can be set to 'Full'"], "beta" : [100 , "# Inverse temperature"], "nom" : [400 , "# Number of Matsubara frequency points sampled"], "M" : [20e6 , "# Total number of Monte Carlo steps"], "mode" : ["SH" , "# We will use self-energy sampling, and Hubbard I tail"], "PChangeOrder" : [0.97 , "# How often to add/rm versus move a kink"], "warmup" : [1e6 , "# Warmup number of QMC steps"],
We will also change the mixing method back to MSR1 in
FeSe.inm file, because the structure will not change
appreciably anymore. Once we converge the solution with the full
Coulomb repulsion, using previously determined optimized structure,
the force should still be appreciably smaller than our previous cutoff
(0.5mRy/a.u.).
The first iteration is not very good, because we are starting without
status.xxx files, therefore the impurity problem starts from
atomic configuration. But nevertheless we achieve rapid convergence of
the self-energy. After a few dmft steps, we can also
increase slightly the number of MC steps, to have better statistics:
"M" : [30e6 , "# Total number of Monte Carlo steps"]
saverage.py sig.inp.2?.1 sig.inp.30.1
maxent_run.py sig.inpx
dmft_copy.py ../
cp ../maxent/Sig.out sig.inp
0 0.025 0.025 200 -1.000000 1.000000 # matsubara, broadening-corr, broadening-noncorr, nomega, omega_min, omega_max (in eV)
x lapw0 -f FeSe x_dmft.py lapw1 x_dmft.py dmft1
x_dmft.py lapw1 --band
x_dmft.py dmftp
wakplot.py 0.1
x_dmft.py dmftu -g --band
e # mode [g/e]: we use mode to compute eigenvalues and eigenvectors BasicArrays.dat # filename for projector 0 # matsubara FeSe.energy # LDA-energy-file, case.energy(so)(updn) FeSe.klist_band # k-list FeSe.rotlm # for reciprocal vectors Udmft.0 # filename for projector 0.0025 # gamma for non-correlated 0.0025 # gammac sig.inp1_band sig.inp2_band # self-energy name, sig.inp(x) eigenvalues.dat # eigenvalues UR.dat # right eigenvector in orbital basis UL.dat # left eigenvector in orbital basis -2. # emin for printed eigenvalues 2. # emax for printed eigenvalues
dmftgk dmftgke.in
wakplot_sophisticated.py
if __name__ == '__main__': aspect_scale=0.3 small = 1e-5 DY = 0.0 # intensity = [3.,3.,3.,1.] COHERENCE=True orb_plot=[1,1,2,2,1,1,0,0,0,0,0,0,0,0] #2,2,2,2,2,2] # should have integers 0,1,2 for red,green,blue #orb_plot=[] # plotting unfolded band structure, but all orbits with hot map
intensity = [0.1,3.,3.,1.]
aspect_scale=1.
./wakplot_sophisticated.py eigenvalues.dat
ln -s UL.dat UL.dat_ ln -s UR.dat UR.dat_
# z2, x2, xz, yz, xy, z2, x2, xz, yz, xy orb_plot=[ 2, 2, 1, 1, 0, 2, 2, 1, 1, 0]
intensity = [3.,3.,3.,1.]
# red,green,blue,transparent intensity = [3., 2., 2., 3.]