Robust wavefront dislocations of Friedel oscillations in gapped graphene

Since the seminal discovery of grapheneNovoselov et al. (2004), two-dimensional materials have attracted wide interest because of their novel physics and great potential applicationsLiu et al. (2020). Usually, two-dimensional materials have high mobility, in which the antiparticles move ballistically and exhibit unconventional quantum tunneling and interferenceKatsnelson et al. (2006); Young and Kim (2009); Cheianov et al. (2007); Beenakker et al. (2009); Wu et al. (2011); Castro Neto et al. (2009). One can intentionally add one or two impurities to form the impurity-design system, this kind of system is charming because it is easily handled theoretically and experimentally, and then it can be regarded as the model system for the exploration of ballistic physicsSettnes et al. (2014a). Experimentally, scanning tunneling spectroscopy (STM) is a proper tool for impurity-design systemBasov et al. (2014). More interesting, the dimension of two-dimensional materials are very unique. On one hand, the bare surface properties of two-dimensional materials are also the bulk properties in contrast to three-dimensional materials. On the other hand, different from the one-dimensional materials, two-dimension materials with two-dimensional parameter space is enough to evolve the global topological quantityDutreix et al. (2020). Therefore, surface-sensitive STM measurement is promising to explore the topological physics of impurity-design two-dimensional materials.

Friedel oscillations (FO) is the quantum interference of electronic waves scattering by the imperfection in crystalline host materialsFriedel (1952). Recently, STM has been demonstrated experimentally to measure the topological Berry phase $\pi$ of monolayer graphene by counting $2$ wavefront dislocations in Friedel oscillations, which is induced by the intentional hydrogen adatom Dutreix et al. (2019). One subsequent experiment shows that for bilayer graphene FOs can exhibit $4$, $2$, or $0$ wavefront dislocations explained by the $2\pi$ Berry phase, the specific sublattice positions of the single impurity and the position of STM tipZhang et al. (2020). Electronic Berry phase as the intrinsic nature of the wave functions is defined in momentum space, which is responsible for many exotic electronic dynamics such as the index shift of quantum Hall effect in monolayer grapheneNovoselov et al. (2005); Zhang et al. (2005) and bilayer grapheneNovoselov et al. (2006), Klein tunnelingKatsnelson et al. (2006); Young and Kim (2009) and the weak antilocalization Wu et al. (2007). The probe of Berry phase usually requires the magnetic fieldXiao et al. (2010); Novoselov et al. (2005); Zhang et al. (2005); Young and Kim (2009). These two experiments not only do not need external magnetic field, but also realize the measurement of Berry phase in real space. However, two experiments both focus on the gapless cases for graphene and emphasize the relation between Berry phase and wavefront dislocation number. This attracts us to explore the effect of gap opening on the interference pattern of FOs.

Gap opening in graphene occurs in various different ways with the substrate coupling as a typical exampleChaves et al. (2020). Most experiments and devices are performed on the substrate-supported graphene, in which lattice mismatch induced inversion symmetry breaking makes Berry phase unquantized multiple of $\pi$Yao et al. (2008). As a result, gap opening should challenge the established correspondence relation between Berry phase and wavefront dislocation number in previous experimentsDutreix et al. (2019); Zhang et al. (2020). In this study, we study the FOs in gapped graphene. Comparing to the gapless graphene, gapped graphene does not have pseudospin-momentum locking, this may prohibit the occurrence of characteristic interference structure (namely, wavefront dislocations) in FOs following the intuitive picture as suggested by the seminal work Dutreix et al. (2019). But the wavefront dislocations in FOs do emerge. Here, we explain the origin of wavefront dislocations by the invariant pseudospin winding number in gaped and gapless graphene, which can be regarded as an updated correspondence relation compatible with the previous experimentsDutreix et al. (2019); Zhang et al. (2020). The wavefront dislocations in FOs of gapped graphene can be verified in the present experimental conditions, and this study helps deepen the understanding of topological physics reflected in the FOs of the impurity-design system.

Fig. 1 gives the schematic measurement of impurity-induced electronic density oscillation of monolayer graphene by STM. The present experimental technology allows one to intentionally introduce a single-atom vacancy on arbitrary sublattice site of grapheneGonzalez-Herrero et al. (2016), e.g., on the sublattice $A$ as shown by the red dot in Fig. 1(a). Corresponding to the introduction of the single vacancy, FOs occur, and lead to the change of the space-resolved and energy-resolved local density of states (LDOS) $\delta\rho(\mathbf{r,}\varepsilon)$Bena (2008); Zou et al. (2016):

$$\delta\rho(\mathbf{r,}\varepsilon)=-\frac{1}{\pi}\operatorname{Im}[\text{Tr}\delta\mathbf{G}(\mathbf{r},\mathbf{r,}\varepsilon)],$$

(1)

where $\delta\mathbf{G}=\mathbf{G-G}^{0}$ represents the change of the total Green’s function (GF) or propagator $\mathbf{G}$ incorporating the effect of the vacancy relevant to the bare propagator $\mathbf{G}^{0}$ of host system (i.e., graphene in Fig. 1(a)), and has the form

$$\delta\mathbf{G}(\mathbf{r}_{2},\mathbf{r}_{1}\mathbf{,}\varepsilon)=\mathbf{G}^{0}(\mathbf{r}_{2},\mathbf{0,}\varepsilon)\mathbf{TG}^{0}(\mathbf{0},\mathbf{r}_{1},\varepsilon).$$

(2)

Here, the $T$-matrix approach is used to describe the effect of vacancy whose potential is simulated by $V_{0}\delta(\mathbf{r})$, and $T$-matrix isDutreix and Katsnelson (2016)

$$\mathbf{T(}\varepsilon\mathbf{)}=\mathbf{V}[1-\int d^{2}\mathbf{kG}^{0}(\mathbf{k,}\varepsilon)\mathbf{V}]^{-1}.$$

(3)

In the $T$-matrix, $\mathbf{V}$ is usually a matrix and its form depends on the specific position of the vacancy, e.g., in Fig. 1(a) (for the other vacancy configuration in supplementary materials), it is

$$\mathbf{V}=\left[\begin{array}[c]{cc}V_{0}&0\\ 0&0\end{array}\right].$$

(4)

For graphene, there are two Dirac valleys in the Brillouin zoneCastro Neto et al. (2009) (cf. Fig. 1(b)), $\mathbf{K}=(\frac{2\pi}{3},\frac{2\pi}{3\sqrt{3}})$ and $\mathbf{K}^{\prime}=(\frac{2\pi}{3},-\frac{2\pi}{3\sqrt{3}})$, then $\delta\rho(\mathbf{r,}\varepsilon)$ are contributed by the intravalley and intervalley scattering. In graphene, the intravalley scattering contribution to $\delta\rho(\mathbf{r,}\varepsilon)$ has been well understood through many theoreticalCheianov and Falko (2006); Bena (2008); Hwang and Das Sarma (2008); Pereg-Barnea and MacDonald (2008); Bena (2009); Pellegrino et al. (2009); Bacsi and Virosztek (2010); Gomez-Santos and Stauber (2011); Lawlor et al. (2013); Settnes et al. (2014b); Settnes et al. (2015); Dutreix and Katsnelson (2016); Rusin and Zawadzki (2018) and experimentalRutter et al. (2007); Brihuega et al. (2008); Mallet et al. (2012); Clark et al. (2014); Gonzalez-Herrero et al. (2016) efforts, while intervalley scattering contribution attracts people’s attention very recently due to its underlying topological natureDutreix et al. (2019); Zhang et al. (2020); Zhang and He (2020). $\delta\rho(\mathbf{r,}\varepsilon)$ is measured easily by STM. FOs of $\delta\rho(\mathbf{r,}\varepsilon)$ are dominated by the backscattering events along the constant energy contourSprunger et al. (1997); Zhang et al. (2019). Focusing on the intervalley scattering, the corresponding physical process of Eq. 2 is shown in Fig. 1(a), i.e., STM tip emits the electron waves from one valley $\mathbf{K}$, as shown by the propagator $G^{0}(\mathbf{K},\mathbf{0,r,}\varepsilon)$, scattering by the vacancy back to STM tip through the other valley $\mathbf{K}^{\prime}$, as shown by the propagator $G^{0}(\mathbf{K}^{\prime},\mathbf{0,r,}\varepsilon)$. Of course, the conjugate process also exists, i.e., the emission/scattering waves are from $\mathbf{K}^{\prime}/\mathbf{K}$ valley. Here, $G^{0}(\mathbf{K/K}^{\prime},\mathbf{0,r,}\varepsilon)$ is the matrix element of $2\times 2$ matrix $G^{0}$ and is associated with the valley momentum index. When the STM tip is shifted around the vacancy in real space as shown in Fig. 1(a), the contributing momentum states to the propagator change their momentum around $\mathbf{K}$ and $\mathbf{K}^{\prime}$ as shown in Fig. 1(c). As a result, the Berry phase defined in momentum space is measured by STM in real space, and the key is the pseudospin-momentum lockingDutreix et al. (2019). In contrast, we will show that the pseudospin winding number instead of Berry phase is measured by STM in gaped graphene without pseudospin-momentum locking.

The physics is essentially the same for gapped monolayer and bilayer graphene, described below using gapped monolayer as an example and the relevant results for gapped bilayer graphene in supplementary materials. The Hamiltonian of gaped monolayer graphene is $H_{0}=v_{F}(\eta\sigma_{x}k_{x}-\sigma_{y}k_{y})+\Delta\sigma_{z}$. $H_{0}$ is expressed in the sublattice basis of $\{A,B\}$, then $\sigma_{x,y,z}$ is the Pauli matrix acting on the pseudospin space. $\pm\Delta$ is the staggered potential on the sublattice $A$ and $B$, which originates from the inversion symmetry breaking, e.g., by the proximity substrateZhou et al. (2007); Chaves et al. (2020). And $\eta=\pm 1$ is valley index for two inequivalent valleys in graphene, $v_{F}=3/2a_{0}t_{0}$ with $a_{0}$ being the carbon-carbon bond length and $t_{0}$ being the nearest-neighbor hopping energy, and $a_{0}$/$t_{0}$ is used as the length/enegy unit in our conventionCastro Neto et al. (2009). The energy spectrum and the spinor wavefunction are, respectively,

$$E_{\xi}(\mathbf{k})=\xi v_{F}\sqrt{\delta^{2}+k^{2}},$$

(5)

and

$$\psi_{\xi,\mathbf{k}}=\frac{1}{\sqrt{1+k^{2}/(\epsilon_{\xi}+\delta)^{2}}}\left[\begin{array}[c]{c}1\\ \frac{\eta k_{x}-ik_{y}}{\epsilon_{\xi}+\delta}\end{array}\right],$$

(6)

where the reduced quantities are defined as $\delta=\Delta/\nu_{F}$, $\epsilon_{\xi}=E_{\xi}/v_{F}$ with $\xi=\pm 1$ for the conduction and the valence band. The GF in momentum space is defined as $G^{0}\left(\mathbf{k,}\varepsilon\right)\equiv(z-H_{0})^{-1}$, and it is

$$\mathbf{G}^{0}\left(\mathbf{k,}\varepsilon\right)=\frac{1}{z^{2}-\Delta^{2}-k^{2}\nu_{F}^{2}}\left[\begin{array}[c]{cc}z+\Delta&\eta v_{F}ke^{i\eta\theta_{k}}\\ \eta v_{F}ke^{-i\eta\theta_{k}}&z-\Delta\end{array}\right].$$

(7)

Here, $z=\varepsilon+i0^{+}$ with $0^{+}$ for the retarded properties of GF and the Fermi level $\varepsilon$ is assumed in the conduction band for brevity. Performing the Fourier transformation to the momentum space GF, we express the real space GFZhu et al. (2011) as $\mathbf{G}^{(0)}(\mathbf{K},\mathbf{r,}\varepsilon)=-e^{i\mathbf{K}\cdot\mathbf{r}}/(2v_{F})^{2}\mathbf{G}^{(0)}(\mathbf{r,}\varepsilon)$ with

$$\mathbf{G}^{(0)}(\mathbf{r,}\varepsilon)=\left[\begin{array}[c]{cc}i\varepsilon_{+}H_{0}\left(u\right)&\eta\sqrt{\varepsilon_{+}\varepsilon_{-}}H_{1}\left(u\right)e^{i\eta\theta_{r}}\\ \eta\sqrt{\varepsilon_{+}\varepsilon_{-}}H_{1}\left(u\right)e^{-i\eta\theta_{r}}&i\varepsilon_{-}H_{0}\left(u\right)\end{array}\right],$$

(8)

where $H_{j}$ is $j$-th order Handel function of the first kind, $\varepsilon_{\pm}=\varepsilon\pm\Delta$, $u=r\sqrt{\varepsilon_{+}\varepsilon_{-}}/v_{F}$, and $r$ ($\theta_{r})$ is the module (azimuthal angle) of $\mathbf{r}$. Concentrating on the intervalley contribution, the change of LDOS is

$$\delta\rho(\Delta\mathbf{K,r},\varepsilon)=\delta\rho_{A}(\mathbf{r},\varepsilon)-\delta\rho_{B}(\mathbf{r},\varepsilon),$$

(9)

where the sublattice-resolved LDOS are

	$\displaystyle\delta\rho_{A}(\mathbf{r},\varepsilon)$	$\displaystyle=C\operatorname{Im}\left[t(\varepsilon)H_{0}^{2}(u)\varepsilon_{+}^{2}\right]\cos(\Delta\mathbf{K}\cdot\mathbf{r}),$		(10a)
	$\displaystyle\delta\rho_{B}(\mathbf{r},\varepsilon)$	$\displaystyle=C\operatorname{Im}\left[t(\varepsilon)H_{1}^{2}(u)\varepsilon_{+}\varepsilon_{-}\right]\cos(\Delta\mathbf{K}\cdot\mathbf{r}-\Delta\eta\theta_{r}).$		(10b)
Here, $\Delta\mathbf{K=K-K}^{\prime}$, $\Delta\eta=\eta-\eta^{\prime}$,$\ C=1/(8v_{F}^{4})$ and $t(\varepsilon)=V_{0}/[1-V_{0}G_{AA}^{(0)}(0,\varepsilon)]$ is the matrix element of $T$-matrix induced by the vacancy on the sublattice $A$.

Nontrivially, Eq. 2 for the LDOS contributed by the intervalley scattering has the identical form as that in gapless graphene Dutreix et al. (2019) and it reproduces the result of gapless case when $\Delta\rightarrow 0$. For the vacancy on the sublattice $A$, $\delta\rho_{A}(\mathbf{r},\varepsilon)$ is trivial. $\left|\Delta\eta\right|=2$ for the intervalley scattering, then the phase of $\delta\rho_{B}(\mathbf{r},\varepsilon)$ is singular at $\mathbf{r}=0$. By shifting STM tip around the vacancy, i.e., $\theta_{r}$ is rotated by $2\pi$, there should be $2$ additional wavefronts in the Friedel oscillation pattern of $\delta\rho_{B}(\mathbf{r},\varepsilon)$Dutreix et al. (2019). In Fig. 2, focusing on the intervalley scattering contribution, we show electronic density oscillations around a single-atom vacancy $\delta\rho_{A}(\mathbf{r},\varepsilon)$ and $\delta\rho_{B}(\mathbf{r},\varepsilon)$ on sublattice $A$ (left panel) and on sublattice $B$ (middle panel), respectively. And the sum $\delta\rho_{A}(\mathbf{r},\varepsilon)+\delta\rho_{B}(\mathbf{r},\varepsilon)$ as total electronic density modulation (right panel). As expected, $\delta\rho_{A}(\mathbf{r},\varepsilon)$ exhibit the normal oscillating wavefronts perpendicular to $\Delta\mathbf{K}$ with a wavelength $\lambda_{\Delta\mathbf{K}}=2\pi/\left|\Delta\mathbf{K}\right|=2.60a_{0}\approx 3.69$ $\mathring{a}$, and does not display any topological feature. $\delta\rho_{B}(\mathbf{r},\varepsilon)$ gives two wavefront dislocations at $\mathbf{r}=0$, which accommodates for the phase accumulated along the contour enclosing the singular point of the phase $\theta_{r}$Nye et al. (1974). In the total electronic density modulation (cf. right panel of Fig. 2), $\delta\rho_{A}(\mathbf{r},\varepsilon)$ only shifts the position of dislocations from $\mathbf{r}=0$ along the direction parallel to $\Delta\mathbf{K}$, but does not change the shape and the number of dislocationsDutreix et al. (2019).

The number $2$ of additional wavefronts is regarded as the signature of the Berry phase $\pi$ of graphene. However, the correspondence between Berry phase and the wavefront dislocation number fails since the Berry phase is unquantized multiple of $\pi$ in gapped graphene. Back to Fig. 1(c), in which we do not follow the referencesDutreix et al. (2019); Zhang et al. (2020) to plot the momentum-resolved pseudospin direction, but it still shows synchronous motion of the dominant scattering state contributing to FOs in momentum space and STM tip in real space, e.g., the scattering state rotates clockwise/counterclockwise on the constant energy contour when STM tip is shifted clockwise/counterclockwise around the vacancy. The scattering state can be used to define the pseudospin vector $\mathbf{s=(}s_{x},s_{y}\mathbf{)}$ where

	$\displaystyle s_{x}$	$\displaystyle=\left\langle\psi_{\xi,\mathbf{k}}\left|\sigma_{x}\right|\psi_{\xi,\mathbf{k}}\right\rangle=\frac{2\eta k_{x}(\epsilon_{\xi}+\delta)}{(\epsilon_{\xi}+\delta)^{2}+k^{2}},$		(11a)
	$\displaystyle s_{y}$	$\displaystyle=\left\langle\psi_{\xi,\mathbf{k}}\left|\sigma_{y}\right|\psi_{\xi,\mathbf{k}}\right\rangle=-\frac{2k_{y}(\epsilon_{\xi}+\delta)}{(\epsilon_{\xi}+\delta)^{2}+k^{2}}.$		(11b)
Due to $s_{x}\sim k_{x}$ and $s_{y}\sim k_{y}$, Eq. 11 gives $4$ fixing points with the momentum azimuthal angle $\theta_{k}=0$, $\pi/2$, $\pi$, and $2\pi$, at which the pseudospin directions are fixed. This feature has no dependence on $\delta$, then is robust to the gap opening. To show more visually, we arbitrarily choose a set of parameters to plot the pseudospin texture in Fig. 3. With the scattering state moving on the constant energy contour, the pseudospin direction twists continuously in Fig. 3(a). Clearly, there are $4$ fixing points shown by the red dots (cf. Fig. 3(a)), which implies the invariant pseudospin winding number. Referring to Fig. 3(b), when $\theta_{k}$ evolves from $0$ to $2\pi$, the pseudospin rotates by $2\pi$ corresponding to the winding number $1$. Alternatively, the invariant pseudospin winding number can also be shown mathematicallyFuchs et al. (2010). To rewrite the gapped Hamiltonian of monolayer graphene into the form:

$$H_{0}(\mathbf{k})=\left|E_{\xi}(\mathbf{k})\right|\left[\begin{array}[c]{cc}\cos\alpha&\sin\alpha e^{-i\phi}\\ \sin\alpha e^{i\phi}&-\cos\alpha\end{array}\right],$$

(12)

where the azimuthal $\alpha(\mathbf{k})$ and polar $\phi(\mathbf{k})$ angles on the Bloch sphere are defined as $\cos\alpha=\Delta/\left|E_{\xi}(\mathbf{k})\right|$, $\sin\alpha=f(\mathbf{k})/\left|E_{\xi}(\mathbf{k})\right|$ with $f(\mathbf{k})=v_{F}(\eta k_{x}+ik_{y})$ and $\phi\equiv-$Arg$f$. As a result, we obtain the Berry phase along a closed Fermi surface:

$$\gamma(C)={\displaystyle\oint\limits_{C}}d\mathbf{k\cdot A_{\xi}}=\pi W_{C}[1-\frac{\Delta}{\left|E_{\xi}(\mathbf{k})\right|}].$$

(13)

Here, we have used the Berry connection $\mathbf{A_{\xi}}\mathbf{=}i\left\langle u_{\mathbf{k},\xi}|\nabla_{\mathbf{k}}u_{\mathbf{k},\xi}\right\rangle=-\xi\sin^{2}\frac{\alpha}{2}\nabla_{\mathbf{k}}\phi$. Most importantly, the winding number is introduced

$$W_{C}\equiv-\xi{\displaystyle\oint\limits_{C}}\frac{d\phi}{2\pi}=\eta\xi.$$

(14)

While the Berry phase of Eq. 13 has the $\Delta$ dependence and becomes unquantized multiple of $\pi$, $W_{C}$ of Eq. 14 with the absolute value $1$ is topologically invariant in gapless and gapped graphene. As a result, the $2$ wavefront dislocations also exist in gapped graphene, which should be ascribed to the invariant winding number.

Here, we firstly discuss the experimental implication of our theoretical discovery. When $\Delta\rightarrow 0$ for gapless graphene, one can say that the Berry phase and winding number are equivalent, then Berry phase is measured by the wavefront dislocationsDutreix et al. (2019); Zhang et al. (2020). But a supplementary experiment is necessary to confirm the gapless nature of graphene sample. In previous two experimentsDutreix et al. (2019); Zhang et al. (2020), there is no such supplementary experiment and they do not perform the quantitative comparison between experiments and simulations for the FOs, which both prohibit the discovery of the possible gap opening. The quantitative calculations of FOs, especially properly incorporating the strength and range of impurity potential, is rather important and is worth of further simulations, e.g., by using density functional calculationsNoori et al. (2020). Also, this experimental discussions are also applicable to one recent theoretical proposal, which uses the presence or absence of wavefront dislocations in FOs to distinguish the incompatible low-energy models for twisted bilayer graphenePhong and Mele (2020). In light of present experimental technology, we expect the immediate verification of our theoretical prediction in gapped graphene. In addition, our results may stimulate the theoretical and experimental research interest on the FOs in intrinsic gapped two-dimensional materials, e.g., transition metal dichalcogenidesWang et al. (2012); Schaibley et al. (2016).

The robustness of the wavefront dislocations of FOs to the gap opening dictates the potential applications in graphene-based pseudospintronicsPesin and MacDonald (2012). Each imperfection can be regarded as one vortex, and the sublattice dependence of wavefront dislocation help define the vortex and antivortex since the orientations of the tripod shape of the imperfection on two sublattices are differentDutreix et al. (2019). Then, we can use the vortex/antivortex as 0/1 to construct the memory device. And the the writing of vortex-based memory device can be performed by using STM to shift the imperfection from one sublattice to the other oneGonzalez-Herrero et al. (2016), while its reading is realized by scanning the imperfection-induced FOs to characterize the orientations of the tripod shape of the imperfection. One recent experimentZhang and He (2020) studies the effect of spacing on the wavefront dislocation of two vortices, and shows the irrelevant vortices on the several nanometer scale which favors the high-density of vortex-based memory device. This implies the great potential of impurity-design graphene in pseudospintronics.

In summary, we have studied the Friedel oscillation in gapped graphene. Although the pseudospin-momentum locking is broken, the wavefront dislocations still emerge in intervalley scattering contributed Friedel oscillations. We establish the correspondence between the invariant winding number instead of the Berry phase and the number of wavefront dislocations in gapped graphene, which can be verified by the present experimental technology. This study is helpful to the understanding of correspondence between topological quantity and wavefront dislocations in Friedel oscillations, and broadens the range of materials to observe the wavefront dislocations of FOs, e,g, in transition metal dichalcogenidesWang et al. (2012); Schaibley et al. (2016).