====== Exercise 3: Electronic band structure of monolayer MoS$_\text{2}$ ======

The electronic band structure gives the relation between crystal momentum $\mathbf{k}$, band index $n$ and the energy $\varepsilon_{n \mathbf{k}}$ of an electron in a crystal, according to Bloch's theorem.

The electronic band structure can be computed using DFT in an approximate way via the Kohn-Sham equations

$$
\left(
-\frac{\nabla^2}{2m} + v_\text{ext}(\mathbf{r}) + v_\text{Hartree}(\mathbf{r}) + v_\text{xc}(\mathbf{r})
\right)
\psi_{n\mathbf{k}}(\mathbf{r})
=
\varepsilon_{n\mathbf{k}} \psi_{n\mathbf{k}}(\mathbf{r})
$$

When using the standard exchange-correlation (xc) functionals like PBE, the band gap between the occupied valence bands and the empty conduction bands is usually underestimated with respect to experiment.

Nevertheless, PBE often gives the correct band ordering, dispersions (i.e., curvature as function of $\mathbf{k}$), and orbital character of the bands.

In this exercise, we compute the band structure of monolayer MoS$_2$, a two-dimensional crystal which has been discovered in 2010 ([[https://doi.org/10.1103/PhysRevLett.105.136805|doi:10.1103/PhysRevLett.105.136805]]).

The following input file can be used as a starting point:
<file none DFT_bandstructure.inp>
&GLOBAL
 PROJECT  MoS2
 RUN_TYPE ENERGY
&END GLOBAL
&FORCE_EVAL
 &DFT
   BASIS_SET_FILE_NAME BASIS_MOLOPT
   POTENTIAL_FILE_NAME GTH_POTENTIALS
   &MGRID
     CUTOFF  500
     REL_CUTOFF  100
   &END MGRID
   &QS
     METHOD GPW
     EPS_DEFAULT 1.0E-12
     EPS_PGF_ORB 1.0E-12
   &END QS
   &SCF
     SCF_GUESS ATOMIC
     EPS_SCF 1.0E-9
     MAX_SCF 500
     &MIXING
         METHOD BROYDEN_MIXING
         ALPHA 0.1
         BETA 1.5
         NBROYDEN 8
     &END
   &END SCF
   &XC
     &XC_FUNCTIONAL PBE
     &END XC_FUNCTIONAL
   &END XC
   &KPOINTS
     SCHEME MONKHORST-PACK 8 8 1
   &END KPOINTS
   &PRINT
     &BAND_STRUCTURE
       ADDED_MOS 10
       FILE_NAME bandstructure.bs
       &KPOINT_SET
         NPOINTS 49
         SPECIAL_POINT GAMMA 0.0 0.0 0.0
         SPECIAL_POINT K     0.333333 0.333333 0.0
         SPECIAL_POINT M     0.0 0.5 0.0
         SPECIAL_POINT GAMMA 0.0 0.0 0.0
       &END KPOINT_SET
     &END BAND_STRUCTURE
   &END PRINT
 &END DFT
 &SUBSYS
   &CELL
     ABC                A B C            ! Comment: FILL HERE "LENGTHS" A B C from C2DB (unit: Angström)
     ALPHA_BETA_GAMMA   ALPHA BETA GAMMA ! Comment: FILL HERE "ANGLES" ALPHA BETA GAMMA from C2DB
     PERIODIC XY
   &END CELL
       &KIND S
     BASIS_SET  DZVP-MOLOPT-GTH
     POTENTIAL  GTH-PBE
   &END KIND

   &KIND Mo
     BASIS_SET  DZVP-MOLOPT-SR-GTH
     POTENTIAL  GTH-PBE
   &END KIND

   &COORD
Mo       X Y Z  ! Comment: FILL HERE POSITION OF Mo FROM XYZ FILE FROM C2DB
S        X Y Z  ! Comment: FILL HERE POSITION OF S  FROM XYZ FILE FROM C2DB
S        X Y Z  ! Comment: FILL HERE POSITION OF S  FROM XYZ FILE FROM C2DB
   &END COORD

   &TOPOLOGY
     &CENTER_COORDINATES
     &END
   &END

 &END SUBSYS
&END FORCE_EVAL
</file>

To complete it, the atomic positions and cell size need to be filled for MoS$_\text{2}$.

<note classic>
**Task:** Navigate to the [[https://c2db.fysik.dtu.dk/|Computational 2D Materials Database (C2DB)]] and
search for [[https://c2db.fysik.dtu.dk/material/1MoS2-1|MoS$_\text{2}$]]. Download the "XYZ" file which contains the atomic positions. Also fill in the lengths of the cell and the angles, as given on the CSDB website. Then run CP2K:
<code>
cp2k.psmp DFT_Bandstructure.inp | tee cp2k.out
</code>
</note>

You can validate the output by comparing to the solution available [[https://github.com/JWilhelm/CP2K_Computational_Methods_in_Crystallography/tree/main/3_MoS2_band_structure_solution|here]].

The PBE band structure is contained in the file ''bandstructure.bs''. You can plot the PBE band structure using a plotting script available via github; obtain the script via
<code>
git clone https://github.com/stefabat/cp2k-scripts
</code>
<note classic>
**Task:** Run the plotting script via
<code>
python3 cp2k-scripts/bin/cp2k_plot_bands.py bandstructure.bs --energy_range -2 3
</code>
and compare your plot to the [[https://github.com/JWilhelm/CP2K_Computational_Methods_in_Crystallography/tree/main/3_MoS2_band_structure_solution/bandstructure.bs.png|png file provided in the solution directory]].
</note>

You can also compare the $k$-path we have chosen in the input to the $k$-path of the hexagonal crystal structure available in the appendix of [[https://www.sciencedirect.com/science/article/pii/S0927025610002697|this paper]].

<note classic>
**Task:** Compare your DFT/PBE band structure to reference data, for example from C2DB. Evaluate the convergence parameters, specifically:

  * Basis set size: You can test larger basis sets like ''TZV2P-MOLOPT-PBE-GTH-q6'' (for S) and ''TZV2P-MOLOPT-PBE-GTH-q14'' (for Mo), contained in the basis file ''BASIS_MOLOPT_UZH'' (also found [[https://github.com/cp2k/cp2k/tree/master/data|here]]).
  * $k$-point mesh: You can test larger $k$-grids, for example ''SCHEME MONKHORST-PACK 16 16 1'' instead of ''SCHEME MONKHORST-PACK 8 8 1''
</note>

<note classic>
**Task:** Compare your DFT/PBE band structure to the $GW$ method, a higher level method deriving from
Green's function theory [[https://doi.org/10.3389/fchem.2019.00377|10.3389/fchem.2019.00377]],
for example from C2DB. What is the main difference of the $GW$ bandstructure compared to DFT/PBE?
</note>

In case you would like to execute the $GW$ band structure calculations of 2D materials with CP2K,
you can check out [[https://doi.org/10.48550/arXiv.2507.18411|10.48550/arXiv.2507.18411]] (Fig. 1,4,10 and Table I) and the corresponding input and output files on [[https://github.com/RemiPasquier/Inputs_Outputs_Periodic_Small_Cell_GW_TMDC|github]].

Can you identify what is the effect of spin-orbit coupling on the band structure described in this paper?