Spectral functions for LSCO
This tutorial explains how to generate densitiesofstates and kresolved spectral functions $A(\mathbf{k},\omega)$ from a local dynamical selfenergies $\Sigma(\omega)$ generated by DFMT.
It has close similarities to generating spectral functions from dynamical selfenergies $\Sigma(\mathbf{k},\omega)$ generated by the GW code; see for example this tutorial.
Table of Contents
Command summary
For any of the sequences below, set up conditions with
touch ctrl.lscoq rm f *.lscoq cp ~/lm/dmft/test/lscoq/{{ctrl,basp,sigm,rst,syml,syml2,indmfl&125;.lscoq,sig.inp%125; .
The following sequence generates “interacting energy bands” along symmetry lines as defined in syml.lscoq.
lmfdmft lscoq rs=1,0 ldadc=82.2 job=1 pade~nw=121~window=4/13.606,4/13.606~icut=30,75
lmchk lscoq  grep A 3 Qlat  tail 3  awk '{print $4, $5, $6}' > qlat
~/lm/startup/symlinepoints qlat=qlat nk=101 syml.lscoq > syml2.lscoq
lmfdmft lscoq vnkabc=8 ldadc=82.2 job=1 gprt~band,fn=syml2~rdsigr=sig2~mode=19 pr45
lmfgws lscoq 'sfuned~units eV~readsek@useef@irrmesh@minmax@ib=30:38~eps 0.02~se band@fn=syml2 nw=2 isp=1 range=4,4'
plbnds sp~atop=10~window=4,4 spq.lscoq ; gnuplot gnu.plt ; open spf.ps
The following sequence generates the noninteracting DOS and the spectrum DOS (the k integrated spectral function)
lmf lscoq vnkabc=12 dos@npts=721@ef0@rdm@window=6/13.606,6/13.606 quit=dos
lmfdmft lscoq vnkabc=8 ldadc=82.2 job=1 gprt~rdsigr=sig2~mode=19
lmfgws lscoq vnkgw=8 'sfuned~units eV~readsek@useef@ib=30:38~eps 0.02~dos isp=1 range=4,4 getev=12 nw=4~savesea~q'
Note: the main part of both sequences are done automatically by invoking ~/lm/dmft/test/test.dmft lscoq 2
Preliminaries
This tutorial works with paramagnetic La_{2}CuO_{4}. QSGW code has no CPAlike capability (in contrast to DMFT, which may be thought of as a dynamical generalization of the Coherent Potential Approximation). At the GW level the material must be treated either as an ordered antiferromagnet, in which case it predicts La_{2}CuO_{4} to be an insulator with a gap of ~4 eV, or nonmagnetically, in which case it predicts La_{2}CuO_{4} to be a metal.
Experimentally, La_{2}CuO_{4} is an insulator with an optical gap slightly less than 2 eV. Both LDA+DMFT and QSGW+DMFT yield an insulator in the paramagnetic state: the $x^2{}y^2$ orbital (which in the nonmagnetic case forms a continuous band that crosses the Fermi level) gets split off.
This tutorial starts from a static nonlocal QSGW selfenergy $\Sigma^0(\mathbf{k})$ and a local, dynamical $\Sigma(i\omega_n)$, projected onto the Cu d states, generated by DMFT for frequencies on the Matsubara axis. The materials system is La_{2}CuO_{4}; files are taken from a validation test taken from the Questaal installation.
It tutorial assumes your build directory is ~/lm, and the executables are your path. The test itself can be run by the script ~/lm/dmft/test/test.dmft.
It requires the following:
 lmfdmft^{†} is used for Pade extrapolation of Σ to the real axis, and also to embed Σ into the QSGW bath.
 lmchk^{†} is used to create file qlat, needed by symlinepoints
 ~/lm/startup/symlinepoints maps kpoint lines in syml.lscoq to syml2.lscoq.
 lmfgws^{†} Makes $A(\mathbf{k},\omega)$ from the embedded $\Sigma(\mathbf{k},\omega)$. lmfgws uses the selfenergy editor, described here
 plbnds^{†} band plot utility, operating in the spectral function mode
 gnuplot^{†}
 open^{†} [use any utility that views postscript files]
^{†}assumed to be in your path.
Interacting Energy bands
This section parallels the tutorial for the interacting band structure generated by GW. In the GW case, $\Sigma(\mathbf{k},\omega)$ is generated from QSGW on the real axis. Here $\Sigma(\mathbf{k},\omega)$ is generated from a local $\Sigma_d(\omega)$ generated by DMFT on the Matsubara axis.
To generate physical observables, we need to analytically continue $\Sigma(i\omega_n)$, given by DMFT for a set of discrete Matsubara frequencies, to the real axis. This is an important difference between Questaal’s QSGW and DMFT implementations.
In this tutorial we use lmfdmft’s built in Pade extrapolator to carry out the analytic continuation.
Note: analytic continuation is a tricky business. Pade works reliably only when the polynomial order is not too large, which means in can operate reliably only over a small energy window. Also the selfenergy should be reasonably smooth.
Continue the local DMFT $\Sigma(i\omega)$ to the real axis in an energy window [−4,4] eV, using a Pade extrapolation:
$ lmfdmft lscoq rs=1,0 ldadc=82.2 job=1 pade~nw=121~window=4/13.606,4/13.606~icut=30,75
The switches to pade are delimited by ~ which is the first character following pade:
 ~nw=121
Create an energy mesh of 121 points  ~window=−4/13.606,4/13.606
Use an energy window of [−4,4] eV (specify in Ry units)  ~icut=30,75
Specifies which Matsubara points to use in the Pade approximant. Use the first 30 points; thereafter increment the spacing by 1,2,4 … points.
Continue until the counter passes point 75. In this case the Pade approximant includes Matsubara points 1:30,32,36,44,60
Next, we need a list of qpoints where the spectral function $A(\mathbf{k},\omega)$ will be made (band plot). Here we use generic symmetry lines from syml.lscoq, which contain highsymmetry points in units of the reciprocal lattice vectors. The symlinepoints utility will convert these points to Cartesian coordinates and at the same time specify a number of points on each line to make the spacing between points approximately uniform. Do the following:
$ lmchk lscoq  grep A 3 Qlat  tail 3  awk '{print $4, $5, $6}' > qlat
$ ~/lm/startup/symlinepoints qlat=qlat nk=101 syml.lscoq > syml2.lscoq
symlinepoints writes points in Questaal’s symmetry line format to file syml2.lscoq.
Next, we need to embed the local DMFT selfenergy $\Sigma_d(\omega)$ (now analytically continued to the real axis) in the bath defined by the oneparticle QSGW hamiltonian. The next step will:
 Determine the Fermi level, so the oneparticle states can be referenced relative to it

Embed $\Sigma_d(\omega)$ into $\Sigma(\mathbf{k},\omega)$, in the basis of QSGW eigenfunctions, for each of the k points supplied through syml2.lscoq.
$ lmfdmft lscoq vnkabc=8 ldadc=82.2 job=1 gprt~band,fn=syml2~rdsigr=sig2~mode=19 pr45
−gprt tells lmfdmft to operate in a special purpose mode which creates and save some property of the Green’s function or self energy. Switches to –gprt are documented here.
lmfdmft first determines the chemical potential μ from $\Sigma_d$ on the Matsubara axis, avoiding errors in the Pade approximation. Next it makes another band pass for the points defined through syml2.lscoq. For each point it embeds $\Sigma_d(\omega)$ (now continued to the real axis) to construct the crystal $\Sigma(\mathbf{k},\omega)$, for the set of bands defined by the DMFT hybridization function. It writes $\Sigma(\mathbf{k},\omega)$ to se.lscoq.
Note: lmfdmft writes only the diagonal part of Σ to se.lscoq. It is assumed that the offdiagonal elements do not significantly modify the spectral function.
We are now in a position to make spectral functions. It is done by the dynamical selfenergy editor lmfgws. $\Sigma(\mathbf{k},\omega)$ constructed either by DMFT and GW can use this utility to analyze spectral functions. For the GW case see this Fe tutorial.
lmfgws can work interactively by adding sfuned unadorned by modifiers:
$ lmfgws lscoq sfuned
or it can be invoked in batch mode. We do the latter here. The following instruction will read se.lscoq, interpolate it to a finer frequency mesh, and write $A(\mathbf{k},\omega)$ to file spq.lscoq.
$ lmfgws lscoq 'sfuned~units eV~readsek@useef@irrmesh@minmax@ib=30:38~eps 0.02~se band@fn=syml2 nw=2 isp=1 range=4,4'
The editor’s instructions are delimited by the first character following sfuned, which is ~ in this case. They perform the following functions:

units eV
Select data is to be output in electron volt units 
readsek@flags
Reads the selfenergy from se.lscoq. Optional flags (separated by a delimiter) modify what is read.
Note: you cannot use ~ for the flags delimiter because it was already used to separate instructions to sfuned (batch mode). The flags perform the following functions: useef
Tells the editor to use the chemical potential in se.lscoq. The QSGW and DMFT Fermi levels are different, and this switch can have some effect if lmfgws is used to interpolate between k points. However, in this context it is not because lmfgws already has the oneparticle levels and selfenergy relative to the Fermi level for every k point; there is no interpolation (in contradistinction to the DOS calculation  irrmesh
Tells the editor that data is not being read on a regular mesh of points suitable for k integration. Operations that require k integration or k interpolation cannot be performed.  minmax
Instructs the editor to print out the minimum and maximum values of QP levels among all the q points it reads. Values are printed for each band. This is for informational purposes only.  ib=30:38
Restricts the bands entering into the spectral function to bands 30 through 38. These are the bands in the energy window [−4,4] eV, as the minmax printout will show.
 useef

eps 0.02
add a constant 0.02 eV to Im Σ. this broaden spectral functions that have small imaginary parts near their peak, rendering a cleaner display. 
se band@fn=syml2 nw=2 range=−4,4
Instructs the editor to generate spectral functions at the list of k points you specify. band@fn=syml2
generates a list of k points in lines defined through syml2.lscoq. @ is a delimiter separating flags to band (In general the band instruction can have multiple flags). By default band operates in symmetry line mode, which read points symmetry lines from syml2.lscoq. Here we synchronize the k points with the same list for which lmfdmft made the self energy. That way $\Sigma$ will not need to be interpolated in q.  nw=n
Refines the given energy mesh by interpolating Σ to an n multiple of the given energy mesh. n must be an integer.  range=#1,#2
Specifies the energy window [#1,#2] where the spectral function is to be generated.
 band@fn=syml2
Finally, invoke plbnds in spectral function mode to generate an “interacting energy bands” figure.
$ plbnds sp~atop=10~window=4,4 spq.lscoq
$ gnuplot gnu.plt
$ open spf.ps [choose your postscript file viewer]
plbnds will generate a gnuplot script file gnu.plt together with a data file spf.lscoq. gnuplot makes a postscript file, spf.ps.
The open merely displays spf.ps.
Spectrum DOS
The spectrum total DOS is a kintegrated version of the k resolved spectral function. Here we generate the noninteracting DOS from the QSGW hamiltonian, and compare it to the DMFT spectrum DOS.
Making the QSGW noninteracting DOS using lmf in the usual manner, using the dos switch
$ lmf lscoq vnkabc=12 dos@npts=721@ef0@rdm@window=6/13.606,6/13.606 quit=dos
It creates file dos.lscoq. quit=dos tells lmf to stop after the DOS has been written, and the modifies to dos have the following meanings:
 npts=721 specifies number of energy points (the spacing between points is the same as for the spectrum dos)
 window=−6/13.606,6/13.606 energy window for DOS, in Ry. It is 3/2 wider than that generated for the spectrum DOS
 ef0 Energy window to be defined relative to Fermi level
 rdm Write DOS in standard Questaal format for 2D arrays.
To make the spectrum DOS, $\Sigma(\omega)$ is needed on the real axis. The steps below are similar those used earlier in the tutorial for energy bands;
$ lmfdmft lscoq rs=1,0 ldadc=82.2 job=1 pade~nw=121~window=4/13.606,4/13.606~icut=30,75
$ lmfdmft lscoq vnkabc=8 ldadc=82.2 job=1 gprt~rdsigr=sig2~mode=19
$ lmfgws lscoq vnkgw=8 'sfuned~units eV~readsek@useef@ib=30:38~eps 0.02~dos isp=1 range=4,4 getev=12 nw=4~savesea~q'
The main change to the lmfdfmt command is that gprt has no band modifier, since we require the selfenergy on a regular mesh to perform kintegration. lmfdfmt writes $\Sigma(\mathbf{k},\omega)$ to se.lscoq.
The main change to the lmfgws instruction is that dos is used instead of se in the selfenergy editor. Note also the modifier getev=12. This interpolates the QP band structure to a finer mesh (12×12×12 divisions), by recalculating the oneparticle levels for the finer mesh. The selfenergy is linearly interpolated from the given 8×8×8 mesh to the 12×12×12 mesh. It writes file sdos.lscoq.
To compare the two DOS, you can use the fplot utility :
$ fplot vef0=.39 sdos.lscoq ab '(x1ef0)*13.6' ord x2/13.6 lt 2,col=1,0,0 dos.lscoq
The noninteracting DOS are written in Ry, on an absolute energy scale (i.e. not relative to the Fermi level).
ab ‘(x1ef0)*13.6’ shifts by ef0 and scales the abscissa from Ry to eV. ord x2/13.6 scales the ordinate from Ry^{−1} to eV^{−1}.
Note that we should use the Fermi level of the noninteracting DOS (0.264 Ry), but we used the Fermi level of the interacting one (0.392 Ry). The mostly noninteracting DOS near 4 eV line up properly that way.
Edit This Page