VO2 monolayer (2D)¶
Source directory: examples/VO2_monolayer_vasp
This example is a non-charge-self-consistent DFT+DMFT setup for a two-dimensional system: a single VO\(_2\) plane on a square lattice with 20 Å of vacuum, using a V-\(d\) \(+\) O-\(p\) Wannier subspace.
Nothing in DMFTwDFT is specific to three dimensions. A monolayer is treated as a slab in a periodic cell, and the two-dimensional character enters only through the k-meshes, the smearing, the Wannier construction, and the band path. This example exists to record those settings in one place, since each of them differs from the bulk examples for a distinct reason.
The deck is constructed as a two-dimensional counterpart to SrVO3
rather than as a model of any known VO\(_2\) phase. Removing the apical oxygens
and the Sr sublattice from bulk SrVO\(_3\) leaves a VO\(_2\) plane, and the in-plane
lattice constant is held at the bulk SrVO\(_3\) value of 3.84652 Å rather than
relaxed, so that the V-O bond length is identical in the two calculations.
Vanadium is V\(^{4+}\), \(d^1\) in both, the Wannier manifold is V \(d\) plus O \(p\) in
both, and the interaction parameters U and J, the inverse temperature
beta \(= 1/k_BT\), and the double-counting scheme are all unchanged. Both decks
use U = 5.0 eV, J = 1.0 eV, and beta = 30.0 eV\(^{-1}\), the last
corresponding to \(T \approx 387\) K.
Differences between the two runs are therefore attributable to dimensionality
rather than to chemistry, bond length, or parameters.
Included files,
input.tomlpara_com.datINCARKPOINTSKPOINTS.bandPOSCARsubmit.sh
See VASP for general VASP setup requirements, including POTCAR, MPI
launcher, and executable setup. A POTCAR for V and O in that order is required
and is not distributed with DMFTwDFT. The deck assumes V_pv.
System¶
Quantity |
Value |
|---|---|
Cell |
\(a = b = 3.84652\) Å, \(c = 20.0\) Å |
Atoms |
V at \((0,0,\tfrac12)\); O at \((\tfrac12,0,\tfrac12)\) and \((0,\tfrac12,\tfrac12)\) |
Coordination |
Square-planar VO\(_4\), aligned with the Cartesian axes |
Formal valence |
V\(^{4+}\) (\(d^1\)), O\(^{2-}\) |
Wannier manifold |
11 orbitals: V \(d\) (5) + 2 \(\times\) O \(p\) (6) |
Electrons in manifold |
13 = \(d^1\) + 2 \(\times\) \(p^6\) |
The layer is placed at \(z = \tfrac12\) so that it sits at the centre of the cell rather than straddling the periodic boundary.
Removing the apical oxygens of the parent perovskite lowers the site symmetry
from \(O_h\) to \(D_{4h}\), which splits the \(d\) shell into four irreducible
representations. cor_orb reflects this decomposition exactly:
Group |
Orbitals |
Irrep |
|---|---|---|
1 |
|
\(a_{1g}\) |
2 |
|
\(b_{1g}\) |
3 |
|
\(b_{2g}\) |
4 |
|
\(e_g\) |
The bulk \(t_{2g}\) triplet is split into d_xy and the d_xz/d_yz pair, and
the lifting of that degeneracy is the physical signature the calculation should
reproduce. All five \(d\) orbitals appear in cor_orb, so the Hartree-Fock path
is not used. Since the square-planar environment is aligned with the Cartesian
axes, L_rot = [0, 0] is correct and no local axis rotation is required.
Two-dimensional settings¶
These are the parameters that differ from a bulk deck, and the reason for each.
KPOINTSis12 12 1. A single k-point along the non-periodic direction.q = [24, 24, 1]. The DMFT k-sum mesh. In a non-charge-self-consistent run the Hamiltonian is Wannier-interpolated onto this mesh directly, soqis independent of the DFTKPOINTSmesh and need not be a multiple of it. A single point along \(k_z\) is correct for a slab: leaving \(q_z\) at the bulk value of 24 would sample a dispersionless direction 24 times over, which is not wrong, only 24 times more expensive.In charge-self-consistent runs driven through the library interface, the mesh is instead refined from the DFT k-points using
nfine(i) = int(q(i)/mp_grid(i)). Thereqmust be at least as large as the DFT mesh componentwise, since a smaller value givesnfine = 0and collapses that direction silently.num_iter_win = 0. This is the setting most specific to slabs and the one most likely to be missed. See below.ISMEAR = 0,SIGMA = 0.05. Gaussian smearing rather than the tetrahedron method used in the bulk examples. The tetrahedron method is unreliable when only one k-point exists along the third direction.KPOINTS.band. The default k-path is simple cubic and includes \(R = (\tfrac12,\tfrac12,\tfrac12)\), which is meaningless for a slab.
Wannier functions in a slab cell¶
With a single k-point along the vacuum direction, the only finite-difference b-vector available along that direction is \(b_z = 2\pi/c\). The spread and the Wannier centre along \(z\) are then determined by a single overlap phase, and that phase is defined only modulo \(2\pi\), so the \(z\) component of the spread carries no useful information.
Because the layer sits at \(z = c/2\), the phase for a state centred in the layer is \(e^{-i b_z c/2} = e^{-i\pi}\), which lies exactly on the branch cut of \(\mathrm{Im}\ln\). Numerical noise flips it between \(+\pi\) and \(-\pi\) from one in-plane k-point to the next, so the apparent centre alternates between \(+c/2\) and \(-c/2\). Since \(\Omega_D\) along a direction is the k-space variance of the apparent centre, this contributes
which for \(c = 20\) Å is 100 Å\(^2\). This example reports \(\Omega_D = 1098.3\) Å\(^2\) for 11 Wannier functions, or 99.8 Å\(^2\) each.
Two features of wannier90.wout confirm that this is bookkeeping rather than
delocalization. Subtracting 100 Å\(^2\) from each reported spread leaves 0.21 to
0.81 Å\(^2\), the same range as the bulk SrVO\(_3\) run. And the reported \(z\)
centres are small but quantized in steps of \(c/N_k\), because each is the mean of
a \(\pm c/2\) distribution and records only how many k-points chose \(+\pi\) over
\(-\pi\). A genuine position would not be quantized on the k-mesh.
Read \(\Omega_I\) instead. It is the only gauge-invariant part of the spread and the only one that tests whether disentanglement selected a good subspace. This example gives \(\Omega_I = 7.90\) Å\(^2\) over 11 Wannier functions, or 0.718 Å\(^2\) each, against 0.716 Å\(^2\) each for bulk SrVO\(_3\). \(\Omega_D\) does not enter \(H(\mathbf{R})\); it is a diagnostic computed from the overlaps.
What does go wrong if the minimization is allowed to run is that it spends its
iterations moving Wannier centres along a direction in which the objective
function is meaningless. Centres drift away from the atomic plane even though
every atom lies in it, and the drift corrupts the orbital character of the
projection. Setting num_iter_win = 0 skips the minimization and keeps the
projected Wannier functions, which is the appropriate choice for a slab and is a
common choice for DFT+DMFT generally. Disentanglement is a separate step and is
still performed; dis_num_iter_win stays at its default.
Note that num_iter_win = 0 does not remove the 100 Å\(^2\) offset, which is
fixed by the cell geometry rather than by the minimizer. It prevents the drift,
not the artifact. Placing the layer at \(z = 0\) would move the phase off the
branch cut and collapse \(\Omega_D\), at the cost of a slab that straddles the
cell boundary; the choice here favours a readable POSCAR over a readable
wannier90.wout.
Running¶
From a copied and edited example directory,
DMFT.py dmft --dft vasp -v
Inspect DMFT/INFO_ITER for convergence. After the DMFT run completes, run
post-processing from inside DMFT,
postDMFT.py ac --average 10
postDMFT.py dos
postDMFT.py bands --plot-plain --omega-points 1000 --band-k-points 1000 \
--normalize --auto-k-path --kpoints ../KPOINTS.band
KPOINTS.band is kept separate from KPOINTS, which holds the
\(12 \times 12 \times 1\) SCF mesh. The --kpoints flag exists so that neither
has to be copied or renamed. The same file serves as the line-mode KPOINTS for
a VASP band-structure run (ICHARG = 11, reusing the SCF CHGCAR) if a
--compare-dft overlay is wanted. In that case pass the SCF OUTCAR to
--outcar, since the absolute energy zero comes from the run that produced the
charge density.
What to check¶
Wannier quality first. Confirm in
wannier90.woutthat the in-plane components of the 11 final Wannier centres sit on the V and O sites, at \((0,0)\) and at \((a/2, 0)\) and \((0, a/2)\). Check \(\Omega_I\) rather than the total spread, for the reasons given above, and ignore the \(z\) components of the centres entirely. A vacuum-localized Wannier function is the characteristic failure mode of a slab calculation and invalidates everything downstream, but it shows up as an in-plane centre off an atomic site or as an inflated \(\Omega_I\), not in \(\Omega_D\).ewin. This is the parameter most likely to need adjustment in 2D. With 20 Å of vacuum the vacuum level sits a few eV above \(E_F\), and free-electron slab states form a dense ladder above it. The window[-7, 7]is a starting point spanning O \(2p\) through the V \(d\) manifold including the strongly antibondingd_x2y2; it may reach the onset of those vacuum states. Inspect the DFT density of states and adjust. If disentanglement selects vacuum states even so, write awannier90.winby hand withdis_froz_minanddis_froz_max, which DMFTwDFT does not write, and pass--no-wintoDMFT.pyso the file is not regenerated.Occupancy.
INFO_ITERshould showNd_lattandNd_impapproaching each other. A persistent gap between them points at the Wannier window or at an unconverged loop rather than at the impurity solver.The dimensional signature. In
ac/Sig.outthe four self-energy groups are no longer degenerate as the bulk \(t_{2g}\) triplet was. Compare the quasiparticle residues against the bulk SrVO\(_3\) run withZ.py --average 5 --cor-orb-index 3and4.
Wannier orbital ordering¶
For postDMFT.py bands --plot-partial -w, the Wannier order follows the
projection block (V:d then O:p), one-based:
Index |
Orbital |
|---|---|
1 |
V |
2 |
V |
3 |
V |
4 |
V |
5 |
V |
6 |
O1 |
7 |
O1 |
8 |
O1 |
9 |
O2 |
10 |
O2 |
11 |
O2 |
The --cor-orb-index argument of plotDMFT.py and Z.py refers to the four
cor_orb groups tabulated earlier instead, and is unrelated to this ordering.