# Fileset

[nanomaterials-10-01154-v2.pdf](https://mdr.nims.go.jp/filesets/c8b8e778-df47-4d56-b85a-4bfd856fb9c6/download)

## Creator

Eberth A. Quezada-López, Zhehao Ge, [Takashi Taniguchi](https://orcid.org/0000-0002-1467-3105), [Kenji Watanabe](https://orcid.org/0000-0003-3701-8119), Frédéric Joucken, Jairo Velasco

## Rights

[Creative Commons BY Attribution 4.0 International](https://creativecommons.org/licenses/by/4.0/)

## Other metadata

[Comprehensive Electrostatic Modeling of Exposed Quantum Dots in Graphene/Hexagonal Boron Nitride Heterostructures](https://mdr.nims.go.jp/datasets/1c0b23a5-5c0d-4c7a-bbf7-4045ea7bbf86)

## Fulltext

nanomaterialsArticleComprehensive Electrostatic Modeling of ExposedQuantum Dots in Graphene/Hexagonal BoronNitride HeterostructuresEberth A. Quezada-López 1,†, Zhehao Ge 1,†, Takashi Taniguchi 2, Kenji Watanabe 3 ,Frédéric Joucken 1 and Jairo Velasco, Jr. 1,*1 Department of Physics, University of California, Santa Cruz, CA 95064, USA;ebaqueza@ucsc.edu (E.A.Q.-L.); zge2@ucsc.edu (Z.G.); fjoucken@ucsc.edu (F.J.)2 International Center for Materials Nanoarchitectronics National Institute for Materials Science, 1-1 Namiki,Tsukuba 305-0044, Japan; taniguchi.takashi@nims.go.jp3 Research Center for Functional Materials National Institute for Materials Science, 1-1 Namiki,Tsukuba 305-0044, Japan; watanabe.kenji.aml@nims.go.jp* Correspondence: jvelasc5@ucsc.edu† These authors contributed equally to this work.Received: 11 May 2020; Accepted: 6 June 2020; Published: 12 June 2020�����������������Abstract: Recent experimental advancements have enabled the creation of tunable localizedelectrostatic potentials in graphene/hexagonal boron nitride (hBN) heterostructures without concealingthe graphene surface. These potentials corral graphene electrons yielding systems akin toelectrostatically defined quantum dots (QDs). The spectroscopic characterization of these exposedQDs with the scanning tunneling microscope (STM) revealed intriguing resonances that are consistentwith a tunneling probability of 100% across the QD walls. This effect, known as Klein tunneling,is emblematic of relativistic particles, underscoring the uniqueness of these graphene QDs. Despitethe advancements with electrostatically defined graphene QDs, a complete understanding of theirspectroscopic features still remains elusive. In this study, we address this lapse in knowledge bycomprehensively considering the electrostatic environment of exposed graphene QDs. We thenimplement these considerations into tight binding calculations to enable simulations of the grapheneQD local density of states. We find that the inclusion of the STM tip’s electrostatics in conjunctionwith that of the underlying hBN charges reproduces all of the experimentally resolved spectroscopicfeatures. Our work provides an effective approach for modeling the electrostatics of exposed grapheneQDs. The methods discussed here can be applied to other electrostatically defined QD systems thatare also exposed.Keywords: graphene; quantum dots; p-n junctions; nanoelectronics1. IntroductionThe endeavor to corral graphene’s massless Dirac fermions has led to the development ofmultiple techniques and novel procedures for charge carrier confinement. These techniques includethe use of lithographic patterning [1–4], ultra-high vacuum chemical synthesis [5–8], controlleddeposition of adatoms [9], and application of perpendicular magnetic fields [10–12]. However, suchtechniques require either complicated fabrication procedures or rely on rigid material synthesis schemes.Recently, a flexible procedure was developed to corral graphene charges that employs a scanningtunneling microscope (STM) and a graphene/hexagonal boron nitride (hBN) heterostructure [13,14].This procedure works by using an STM tip to locally induce charge in the underlying hBN, thus creatingpersistent and embedded local gates within hBN. These local gates enable the corralling of electrons inNanomaterials 2020, 10, 1154; doi:10.3390/nano10061154 www.mdpi.com/journal/nanomaterialshttp://www.mdpi.com/journal/nanomaterialshttp://www.mdpi.comhttps://orcid.org/0000-0003-3701-8119http://dx.doi.org/10.3390/nano10061154http://www.mdpi.com/journal/nanomaterialshttps://www.mdpi.com/2079-4991/10/6/1154?type=check_update&version=2Nanomaterials 2020, 10, 1154 2 of 12an exposed circular p-n junction, which effectively behaves as a quantum dot (QD) [14]. Such QDshave been used to develop novel electronic devices such as Berry phase switches [15] and have enabledthe unprecedented visualization of correlated relativistic charges under large magnetic fields [16].Despite the rampant progress on these exposed graphene QDs, their spatially resolved spectroscopiccharacterization still remains not well understood.An important commonality in existing works on exposed graphene QDs is that the influence ofthe STM tip has been excluded in the theoretical modeling of QD states. This is in spite of severalexperiments that show the STM tip itself can induce a small QD [17], and even create coulomb-likeconfinement [18]. In this study, we address the influence of the STM tip on exposed graphene QDsby comparing our visualization of QD states with tight-binding (TB) calculations that include theelectrostatics of the STM tip and underlying charged hBN. Our results demonstrate that accountingfor the influence of the STM tip is necessary to reproduce key features seen in the experimental data.Additionally, we use the insight from our study to show how the tip’s influence can be mitigated bychoosing an appropriately sized QD.2. Materials and MethodsThe experiments we present in this work were performed on heterostructures composed of asingle graphene layer on a 45 nm thick hBN resting on a SiO2/p-doped Si substrate as depicted inFigure 1a. Both graphene and hBN were mechanically exfoliated from bulk crystals. The graphene/hBNheterostructure was assembled using a standard polymer-based transfer technique [19]. Following thisassembly, the surface of graphene was cleared of debris and polymer residue using a Cypher S atomicforce microscope from Asylum (High Wombe, UK) while in contact mode [20]. We perform this step toensure the absence of contaminants that could affect the electronic properties of our QDs. The creationand characterization of graphene QDs were performed with a commercial low-temperature STM(Erligheim, Germany) from Createc operating at 4.8 K.Nanomaterials 2020, 10, x FOR PEER REVIEW 2 of 11  quantum dot (QD) [14]. Such QDs have been used to develop novel electronic devices such as Berry phase switches [15] and have enabled the unprecedented visualization of correlated relativistic charges under large magnetic fields [16]. Despite the rampant progress on these exposed graphene QDs, their spatially resolved spectroscopic characterization still remains not well understood. An important commonality in existing works on exposed graphene QDs is that the influence of the STM tip has been excluded in the theoretical modeling of QD states. This is in spite of several experiments that show the STM tip itself can induce a small QD [17], and even create coulomb-like confinement [18]. In this study, we address the influence of the STM tip on exposed graphene QDs by comparing our visualization of QD states with tight-binding (TB) calculations that include the electrostatics of the STM tip and underlying charged hBN. Our results demonstrate that accounting for the influence of the STM tip is necessary to reproduce key features seen in the experimental data. Additionally, we use the insight from our study to show how the tip’s influence can be mitigated by choosing an appropriately sized QD. 2. Materials and Methods The experiments we present in this work were performed on heterostructures composed of a single graphene layer on a 45 nm thick hBN resting on a SiO2/p-doped Si substrate as depicted in Figure 1a. Both graphene and hBN were mechanically exfoliated from bulk crystals. The graphene/hBN heterostructure was assembled using a standard polymer-based transfer technique [19]. Following this assembly, the surface of graphene was cleared of debris and polymer residue using a Cypher S atomic force microscope from Asylum (High Wombe, UK) while in contact mode [20]. We perform this step to ensure the absence of contaminants that could affect the electronic properties of our QDs. The creation and characterization of graphene QDs were performed with a commercial low-temperature STM (Erligheim, Germany) from Createc operating at 4.8 K. Figure 1a shows a schematic of the STM circuit used. A sample bias ( ) is applied between graphene and the grounded STM tip to incite a tunneling current between them and enable probing of electronic states. The back-gate ( ) connected to the p-doped Si layer is used to remotely tune graphene’s Fermi level ( ). We use these control voltages to create a graphene QD by applying a  pulse between graphene and the STM tip while maintaining  at a constant value [13–16,21]. During the application of this pulse, defects in hBN underneath the tip become ionized with opposite polarity to . The end result of this process is depicted in Figure 1b for the case where > 0. At a fixed positive value of , graphene becomes globally n-doped except for the circular region where the pulse was applied. In this region, an excess of negative charges embedded in the hBN acts as a local back-gate that p-dopes graphene. The resulting spatial variation of the charge neutrality point (also known as the Dirac point in graphene) with respect to  (dashed orange line) gives rise to an electrostatic potential ( ). A profile of this electrostatic potential is outlined in Figure 1b (red curve) and specifies the boundary of the exposed graphene QD.  Figure 1. Schematic layout and potential of an exposed graphene Quantum Dot (QD). (a) Schematic showing the scanning tunneling microscope (STM) circuit and graphene/hexagonal boron nitride Figure 1. Schematic layout and potential of an exposed graphene Quantum Dot (QD). (a) Schematicshowing the scanning tunneling microscope (STM) circuit and graphene/hexagonal boron nitride (hBN)heterostructure. The circular p-doped region outlines the QD created by applying a bias voltage (VS)pulse between the tip and the graphene/hBN heterostructure while holding the back-gate voltage (VG)constant. (b) Top: UD(r) (red curve) is the potential of the QD which is outlined by tracing along theDirac point in each cone. Bottom: Side-view schematic of the QD in a graphene/hBN heterostructurefor VG > 0. The application of a high electric field by the STM tip induces a localized net chargeaccumulation after exciting defects in hBN. Applying VG with opposite polarity to charges in hBNinduces the spatial variation in doping on graphene which forms the QD boundaries.Figure 1a shows a schematic of the STM circuit used. A sample bias (VS) is applied betweengraphene and the grounded STM tip to incite a tunneling current between them and enable probingof electronic states. The back-gate (VG) connected to the p-doped Si layer is used to remotely tuneNanomaterials 2020, 10, 1154 3 of 12graphene’s Fermi level (EF). We use these control voltages to create a graphene QD by applying aVS pulse between graphene and the STM tip while maintaining VG at a constant value [13–16,21].During the application of this pulse, defects in hBN underneath the tip become ionized with oppositepolarity to VG. The end result of this process is depicted in Figure 1b for the case where VG > 0. At afixed positive value of VG, graphene becomes globally n-doped except for the circular region wherethe pulse was applied. In this region, an excess of negative charges embedded in the hBN acts as alocal back-gate that p-dopes graphene. The resulting spatial variation of the charge neutrality point(also known as the Dirac point in graphene) with respect to EF (dashed orange line) gives rise to anelectrostatic potential UD(r). A profile of this electrostatic potential is outlined in Figure 1b (red curve)and specifies the boundary of the exposed graphene QD.3. ResultsTo map and visualize the electronic properties of the states in the exposed graphene QD,we use scanning tunneling microscopy (STM) and spectroscopy (STS). Figure 2a shows an STMtopographic map of graphene’s ultra-flat surface after the creation of a QD. For the case whereVG > 0, the approximate regions where graphene is p and n-doped are colored brown and green,respectively. With STS we can obtain the differential conductance (dI/dVS) at a specific point ongraphene underneath the STM tip. This dI/dVS signal is proportional to graphene’s local density ofstates (LDOS) [22]. By performing this measurement at different points within the QD, we reveal thecoarse spatial dependence of the QD states. In Figure 2b we plot dI/dVS as a function of VS taken atthe center (black curve) and boundaries (blue and orange curves) of the QD at points corresponding tothe colored crosses in Figure 2a. These curves clearly display differences between the signals recoveredat the center and edges of the QD. At the edges, the dI/dVS curves have prominent peaks while at thecenter, these peaks are suppressed and broadened.To attain a more comprehensive understanding of the spatial dependence of the exposed grapheneQD’s LDOS, we obtain dI/dVS(VS) curves at multiple points along the dashed cyan line in Figure 2a.Figure 2c–e show the compiled dI/dVS(VS) curves plotted as a function of distance, where the origin isdefined at the center of the graphene QD. Additionally, each of these image plots are taken at differentvalues of VG, which offsets the global graphene doping. Following the schematic in Figure 1b, as VGchanges, EF also changes relative to the QD’s potential (UD(r)). For Figure 2c, EF is near the top of thepotential UD(r), which creates a shallow QD. From the image plots in Figure 2d,e, it is apparent that asVG decreases (hole density increases), the QD gains depth and width as the difference between EF andUD(r) increases.The patterns and features observed in Figure 2b–e can be explained by considering the behavior ofmassless Dirac fermions corralled within a circular and harmonic electrostatic potential. Due to Kleintunneling, a p-n junction on graphene perfectly transmits quasiparticles at normal incidence to thejunction but reflects them at larger incident angles [23–25]. Therefore, in a circular potential, electronswith high angular momenta have oblique incidence with the barrier and become internally reflected.This leads to charge carrier trapping and the formation of resonant states [26–30] with pronouncedintensities near the boundary of the circular potential, in agreement with Figure 2b. Additionally,in Figure 2b we see evidence of the internal reflection due to Klein tunneling. This can be seen inthe manifestation of differing peak widths between the dI/dVS curve taken at the QD’s center andcurves taken at the QD’s edges. As electrons with high angular momenta get trapped near the edges,these states exhibit longer trapping times and thus narrower spectroscopic peaks [14,17]. Moreover,the bright nodal features in Figs 2c–e can be attributed to the eigenstates of the exposed grapheneQD [14]. The profile of the confinement potential in this QD is parabolic, akin to that of a harmonicoscillator. However, unlike Schrödinger fermions in a harmonic potential, these nodal patterns areunevenly spaced in energy. Instead, the nodal patterns formed by these massless Dirac fermionsbecome more closely packed as VS decreases (see. Figure 2d) [14,17,28–30].Nanomaterials 2020, 10, 1154 4 of 12Nanomaterials 2020, 10, x FOR PEER REVIEW 4 of 11   Figure 2. Scanning Tunneling Spectroscopy (STS) of an exposed graphene QD. (a) Topography of graphene acquired with the scanning tunneling microscope. Regions where graphene is n and p-doped are indicated with green and brown color scales, respectively. (b) Differential conductance spectra ( / ) obtained at the corresponding crosses indicated in (a). Each spectrum is offset for clarity and was taken with  = 22 V. (c–e) Spatial dependence of /  spectra with different , which is indicated in each panel. These images map the spatial dependence of graphene QD states obtained along the dashed cyan line in (a). Dashed vertical lines in (d) correspond to spectra displayed in (b). Tunneling parameters:  = −0.1 V,  = 1 nA,  = 2 mV. In addition to the well understood features described above, there are some features that lack explanation. For example, all three plots in Figure 2c–e show a bright skirt-like feature around the edges of the QD (~  = −100 mV). A clear downward bending of the QD states is also visible near the QD boundaries for all values of . This bending effect is particularly pronounced in Figure 2c, where the strong distortion of states creates an envelope-like feature. As we will soon show, these features can be reproduced after comprehensive consideration of the QD electrostatic environment. 4. Discussion To study the effect of the STM tip on exposed graphene QDs, we use simplified electrostatics and a numerical tight binding model. We first discuss our considerations for the electrostatics. In our experiments, spatial variations in doping across graphene originate from localized hBN defect charges (see Figure 1b) and inadvertent gating from the STM tip. To first approximation, the localized hBN defect charges will create a fixed spatially varying doping profile. On the other hand, the STM tip will create a mobile doping profile that changes with the STM tip’s position. An exact solution for the charge density in graphene would require inaccessible experimental parameters such as the Figure 2. Scanning Tunneling Spectroscopy (STS) of an exposed graphene QD. (a) Topography ofgraphene acquired with the scanning tunneling microscope. Regions where graphene is n and p-dopedare indicated with green and brown color scales, respectively. (b) Differential conductance spectra(dI/dVS) obtained at the corresponding crosses indicated in (a). Each spectrum is offset for clarity andwas taken with VG = 22 V. (c–e) Spatial dependence of dI/dVS spectra with different VG, which isindicated in each panel. These images map the spatial dependence of graphene QD states obtainedalong the dashed cyan line in (a). Dashed vertical lines in (d) correspond to spectra displayed in (b).Tunneling parameters: VS = −0.1 V, I = 1 nA, Vac = 2 mV.In addition to the well understood features described above, there are some features that lackexplanation. For example, all three plots in Figure 2c–e show a bright skirt-like feature around theedges of the QD (~VS = −100 mV). A clear downward bending of the QD states is also visible near theQD boundaries for all values of VG. This bending effect is particularly pronounced in Figure 2c, wherethe strong distortion of states creates an envelope-like feature. As we will soon show, these featurescan be reproduced after comprehensive consideration of the QD electrostatic environment.4. DiscussionTo study the effect of the STM tip on exposed graphene QDs, we use simplified electrostaticsand a numerical tight binding model. We first discuss our considerations for the electrostatics. In ourexperiments, spatial variations in doping across graphene originate from localized hBN defect charges(see Figure 1b) and inadvertent gating from the STM tip. To first approximation, the localized hBNdefect charges will create a fixed spatially varying doping profile. On the other hand, the STM tipwill create a mobile doping profile that changes with the STM tip’s position. An exact solution for thecharge density in graphene would require inaccessible experimental parameters such as the spatialNanomaterials 2020, 10, 1154 5 of 12distribution of hBN defect charges as well as the in-situ STM tip’s geometry. To circumvent thesedifficulties we proceed by making a set of simple approximations for the doping profiles due to thehBN defects and the STM tip.We first focus on the doping profile due to charged defects in the underlying hBN. The potentialprofile of a graphene QD can be extracted by tracking the spatial evolution of the region with reduceddI/dVS intensity. This region corresponds to the spatially varying Dirac point and can be seen inmeasurements similar to those shown in Figure 2c–e. Figure 3a shows an example of an extractedpotential profile UD(r) (red curve). This UD(r) is then converted to a doping profile n(r) by therelation n(r) = sgn[UD(r)] ∗UD(r)2}2v2Fπ; with vF = 1× 106 m/s, where UD(r) replaces the energy term [31].The resulting plot after smoothing is shown in Figure 3b. We note that the STM data in Figure 3a, whichis the source of our estimate for the hBN defect potential, necessarily include the effect of the STM tip.To remove this effect so that we may treat it separately, we preemptively reduce the lateral extent of thepotential profile in Figure 3b to 65% while leaving the energy scale unchanged. The resulting adjustedcharge density profile is shown in Figure 3c.Nanomaterials 2020, 10, x FOR PEER REVIEW 5 of 11  spatial distribution of hBN defect charges as well as the in-situ STM tip’s geometry. To circumvent these difficulties we proceed by making a set of simple approximations for the doping profiles due to the hBN defects and the STM tip. We first focus on the doping profile due to charged defects in the underlying hBN. The potential profile of a graphene QD can be extracted by tracking the spatial evolution of the region with reduced /  intensity. This region corresponds to the spatially varying Dirac point and can be seen in measurements similar to those shown in Figure 2c–e. Figure 3a shows an example of an extracted potential profile ( ) (red curve). This ( ) is then converted to a doping profile ( ) by the relation ( ) sgn ( ) ∗ ( ) ; with 1 10  m/s, where ( ) replaces the energy term [31].The resulting plot after smoothing is shown in Figure 3b. We note that the STM data in Figure 3a, which is the source of our estimate for the hBN defect potential, necessarily include the effect of the STM tip. To remove this effect so that we may treat it separately, we preemptively reduce the lateral extent of the potential profile in Figure 3b to 65% while leaving the energy scale unchanged. The resulting adjusted charge density profile is shown in Figure 3c.  Figure 3. Charge density profile induced on graphene by charged hBN defects. (a) Experimentally approximated QD potential profile. The underlying two-dimensional (2D) plot shows the spatial dependence of QD states (similar to Figure 2c–e). The red line (same as ( ) in Figure 1b) represents the extracted potential profile for the graphene QD. (b) QD’s charge density profile converted and smoothed from the extracted potential profile in (a). (c) QD’s charge density profile used in our tight-binding (TB) calculation after adjustment. This profile’s x-axis is 65% of the profile in (b). After obtaining the induced charge density profile on graphene due to the charged hBN defects alone, we proceed to estimate the induced charge density due to the STM tip. In our experiments, we use a tip made of tungsten. Since tungsten has a different work function than graphene, there is a finite work function mismatch (ΔΦ) between the STM tip and graphene. In the tunneling regime, the STM tip remains at a distance ~7.5 Å from the graphene surface. Because of the finite ΔΦ, there is a shift of the graphene bands even when  = 0 and  = 0. For |ΔΦ| ≫ | |, the polarity and intensity of the tip induced doping is dominated by ΔΦ. Therefore, in this regime we can acquire an estimate for the doping profile induced by the STM tip by obtaining an approximation of ΔΦ. To get an estimate of ΔΦ and ultimately the electrostatic effect of the STM tip, we measure and plot the dependence of /  as function of  and  on pristine graphene prior to creating a QD (see Figure 4a). The resulting two-dimensional (2D) image plot displays changes in the tunneling current between the tip and graphene as we vary their relative band alignments. In Figure 4a we observe several confined states appearing as a result of the STM tip’s local top gating effect, which were previously reported by Y. Zhao et al. [17]. As we vary  for different values of  we obtain a /  signal that includes the contribution from two channels during the tunneling process. This is a common occurrence for low dimensional systems with low charge density [32,33]. The first channel corresponds to the differential tunneling current between states of the STM tip and states in graphene at an energy determined by the sample bias. Since the horizontal axis  roughly indicates Figure 3. Charge density profile induced on graphene by charged hBN defects. (a) Experimentallyapproximated QD potential profile. The underlying two-dimensional (2D) plot shows the spatialdependence of QD states (similar to Figure 2c–e). The red line (same as UD(r) in Figure 1b) representsthe extracted potential profile for the graphene QD. (b) QD’s charge density profile converted andsmoothed from the extracted potential profile in (a). (c) QD’s charge density profile used in ourtight-binding (TB) calculation after adjustment. This profile’s x-axis is 65% of the profile in (b).After obtaining the induced charge density profile on graphene due to the charged hBN defectsalone, we proceed to estimate the induced charge density due to the STM tip. In our experiments,we use a tip made of tungsten. Since tungsten has a different work function than graphene, there is afinite work function mismatch (∆Φ) between the STM tip and graphene. In the tunneling regime, theSTM tip remains at a distance ~7.5 Å from the graphene surface. Because of the finite ∆Φ, there is ashift of the graphene bands even when VS = 0 and VG = 0. For |∆Φ| �|VS|, the polarity and intensityof the tip induced doping is dominated by ∆Φ. Therefore, in this regime we can acquire an estimatefor the doping profile induced by the STM tip by obtaining an approximation of ∆Φ.To get an estimate of ∆Φ and ultimately the electrostatic effect of the STM tip, we measure andplot the dependence of dI/dVS as function of VS and VG on pristine graphene prior to creating a QD(see Figure 4a). The resulting two-dimensional (2D) image plot displays changes in the tunnelingcurrent between the tip and graphene as we vary their relative band alignments. In Figure 4a weobserve several confined states appearing as a result of the STM tip’s local top gating effect, whichwere previously reported by Y. Zhao et al. [17]. As we vary VS for different values of VG we obtain adI/dVS signal that includes the contribution from two channels during the tunneling process. This is acommon occurrence for low dimensional systems with low charge density [32,33]. The first channelcorresponds to the differential tunneling current between states of the STM tip and states in grapheneNanomaterials 2020, 10, 1154 6 of 12at an energy determined by the sample bias. Since the horizontal axis VG roughly indicates changes inn and because EF ∝ sgn[n] ×√|n| for graphene, tip-induced states exhibit inverted “S”-shape fans inthis channel. To obtain an estimate for ∆Φ, however, we focus on the second channel. This channelcorresponds to the differential tunneling current between the topmost filled states in graphene andstates in the STM tip. Consequently, signal traces in this channel will manifest as lines along whichgraphene retains a constant charge density n [21,34]. One such line appears in Figure 4a enclosed by anorange box. This feature results from the charging of a tip-induced confined state as it aligns with EF.Nanomaterials 2020, 10, x FOR PEER REVIEW 6 of 11  changes in  and because ∝ sgn | | for graphene, tip-induced states exhibit inverted “S”-shape fans in this channel. To obtain an estimate for ΔΦ, however, we focus on the second channel. This channel corresponds to the differential tunneling current between the topmost filled states in graphene and states in the STM tip. Consequently, signal traces in this channel will manifest as lines along which graphene retains a constant charge density  [21,34]. One such line appears in Figure 4a enclosed by an orange box. This feature results from the charging of a tip-induced confined state as it aligns with .  Figure 4. Determination of the charge density profile induced by the STM tip on graphene. (a) Experimentally determined / ( , ) on pristine graphene. This two-dimensional (2D) plot is used to determine the magnitude of the STM tip’s band shifting effect on graphene. The inverted “S” shaped features correspond to direct tunneling into tip-induced confined states. Diagonal features correspond to the charging of spectroscopic features as they coincide with . The dashed orange box encloses a charging feature due to a tip-induced state and the green dashed line outlines a /  suppression belonging to graphene’s charge neutrality point (CNP). The dashed white line indicates where 0. Tunneling parameters:  = −0.2 V,  = 1 nA,  = 2 mV. (b) Charge density profile induced on graphene by the STM tip. The induced charge density intensity is extracted from the extrapolation in (a). The shape of this density profile is determined by assuming the tip to be a charged sphere with a radius of 80 nm. After identifying the origin of the tunneling features in Figure 4a for different band alignment configurations, we are ready to obtain an estimate for ΔΦ. First, we note that for pristine graphene and for ΔΦ 0, graphene’s charge neutrality (CNP) point will cross the Fermi level at 0 mV and 0  in a /  ( , ) plot. However, because ΔΦ 0, graphene’s CNP will cross the Fermi level at nonzero  and  values. To this end, we find the value of  at which graphene’s charge neutrality point (CNP) crosses  at  = 0. With this  value we can estimate the STM tip’s top gating effect due to ΔΦ on the surface of graphene in the absence of gating from below. We follow the slope of the charging feature inside the orange box and trace a green dashed line with the same slope along the furthest dark fringe on the right (see Figure 4a). The suppression of /  along this dark fringe indicates where graphene’s CNP crosses . Finally, we extend this green dashed line downward and find that it crosses  = 0 at  −290 mV. We note that this value for  is greater in magnitude than the  range in our measurements (−100 m   100 mV). Thus, we can reasonably assume that the effect of ΔΦ dominates within our experimental  range used to map the QD states. With an estimate for the graphene band shift due to ΔΦ between the STM tip and graphene we can approximate the profile for the charge induced on graphene by the STM tip. If we assume the STM tip’s apex and graphene act as a parallel plate capacitor with a 7.5  separation, then the  −290 mV offset corresponds to a maximum tip induced charge density of ~ 2.14 10  cm . We calculated the shape of the tip’s doping profile by using a standard Poisson solver; where the tip is modeled by a charged sphere with an 80 nm radius that is placed 7.5  away from a metal surface. Figure 4. Determination of the charge density profile induced by the STM tip on graphene.(a) Experimentally determined dI/dVS(VG, VS) on pristine graphene. This two-dimensional (2D)plot is used to determine the magnitude of the STM tip’s band shifting effect on graphene. The inverted“S” shaped features correspond to direct tunneling into tip-induced confined states. Diagonal featurescorrespond to the charging of spectroscopic features as they coincide with EF. The dashed orange boxencloses a charging feature due to a tip-induced state and the green dashed line outlines a dI/dVSsuppression belonging to graphene’s charge neutrality point (CNP). The dashed white line indicateswhere VS = 0. Tunneling parameters: VS = −0.2 V, I = 1 nA, Vac = 2 mV. (b) Charge density profileinduced on graphene by the STM tip. The induced charge density intensity is extracted from theextrapolation in (a). The shape of this density profile is determined by assuming the tip to be a chargedsphere with a radius of 80 nm.After identifying the origin of the tunneling features in Figure 4a for different band alignmentconfigurations, we are ready to obtain an estimate for ∆Φ. First, we note that for pristine grapheneand for ∆Φ = 0, graphene’s charge neutrality (CNP) point will cross the Fermi level at VS = 0 mVand VG = 0 V in a dI/dVS (VG, VS) plot. However, because ∆Φ , 0, graphene’s CNP will cross theFermi level at nonzero VG and VS values. To this end, we find the value of VS at which graphene’scharge neutrality point (CNP) crosses EF at VG = 0. With this VS value we can estimate the STM tip’stop gating effect due to ∆Φ on the surface of graphene in the absence of gating from below. We followthe slope of the charging feature inside the orange box and trace a green dashed line with the sameslope along the furthest dark fringe on the right (see Figure 4a). The suppression of dI/dVS along thisdark fringe indicates where graphene’s CNP crosses EF. Finally, we extend this green dashed linedownward and find that it crosses VG = 0 at VS ≈ −290 mV. We note that this value for VS is greater inmagnitude than the VS range in our measurements (−100 m < VS <100 mV). Thus, we can reasonablyassume that the effect of ∆Φ dominates within our experimental VS range used to map the QD states.With an estimate for the graphene band shift due to ∆Φ between the STM tip and graphene wecan approximate the profile for the charge induced on graphene by the STM tip. If we assume the STMtip’s apex and graphene act as a parallel plate capacitor with a 7.5 Å separation, then the VS ≈ −290 mVoffset corresponds to a maximum tip induced charge density of ∼ 2.14×1012 cm−2. We calculated theshape of the tip’s doping profile by using a standard Poisson solver; where the tip is modeled by acharged sphere with an 80 nm radius that is placed 7.5 Å away from a metal surface. The tip radiusNanomaterials 2020, 10, 1154 7 of 12and distance to graphene that we used are both consistent with values found in the literature [18,35,36].With this analysis we acquire the tip induced doping profile shown in Figure 4b.We obtained approximations for the charge density profiles induced on graphene due to chargeddefects in the underlying hBN and the STM tip. We proceed by adding these contributions to obtaina potential profile resulting from the summation of the induced charge densities. After adding thecharge densities induced on graphene from Figures 3d and 4b, we convert the resulting charge densityn(r) into a potential profile UD(r); where UD(r) = sgn[n(r)] ∗ }vF√π∣∣∣n(r)∣∣∣. Figure 5a–c shows 2Dmaps of the potentials resulting from the cumulative charge densities of the hBN defects and STM tip.The color scale corresponds to the potential value; where red and blue indicate high and low values,respectively. Notably, each of the potential maps differ because the position of the tip changes betweenthem. As a comparison, we also show a 2D potential map without the effect of the STM tip (Figure 5d).When the STM tip is at the center (Figure 5a) or 50 nm away from the center (Figure 5b), the potentialmap has the highest value at the STM tip’s location, as indicated by the red dot at the correspondinglocations. Additionally, we plot line cuts of the potential maps in Figure 5e–g. Here the QD’s potentialreveals a distorted profile with a prominent peak. For maps with the STM tip 100 nm away from thecenter (Figure 5c), the potential profile has two separate peaks (Figure 5g); effectively becoming anasymmetric double QD system.Nanomaterials 2020, 10, x FOR PEER REVIEW 7 of 11  The tip radius and distance to graphene that we used are both consistent with values found in the literature [18,35,36]. With this analysis we acquire the tip induced doping profile shown in Figure 4b. We obtained approximations for the charge density profiles induced on graphene due to charged defects in the underlying hBN and the STM tip. We proceed by adding these contributions to obtain a potential profile resulting from the summation of the induced charge densities. After adding the charge densities induced on graphene from Figures 3d and 4b, we convert the resulting charge density ( ) into a potential profile ( ); where ( ) sgn ( ) ∗ | ( )|. Figure 5a–c shows 2D maps of the potentials resulting from the cumulative charge densities of the hBN defects and STM tip. The color scale corresponds to the potential value; where red and blue indicate high and low values, respectively. Notably, each of the potential maps differ because the position of the tip changes between them. As a comparison, we also show a 2D potential map without the effect of the STM tip (Figure 5d). When the STM tip is at the center (Figure 5a) or 50 nm away from the center (Figure 5b), the potential map has the highest value at the STM tip’s location, as indicated by the red dot at the corresponding locations. Additionally, we plot line cuts of the potential maps in Figure 5e–g. Here the QD’s potential reveals a distorted profile with a prominent peak. For maps with the STM tip 100 nm away from the center (Figure 5c), the potential profile has two separate peaks (Figure 5g); effectively becoming an asymmetric double QD system.  Figure 5. Simulation of the local density of states (LDOS) of an exposed graphene QD with an STM tip at a fixed position. (a–c) Spatial 2D map of graphene’s Dirac point energy with respect to  after adding the contributions from hBN defects (Figure 3c) and the STM tip (Figure 4b). The STM tip’s location varies for each of these maps. (d) Spatial 2D map of the graphene Dirac point energy with respect to  determined by only considering the contribution from hBN defects (Figure 3c). (e–h) Respective potential profile line cuts along the yellow dashed lines in (a–d). (i–l) Calculated LDOS distributions for each corresponding potential map. These simulations reveal the effect of the STM tip on the QD states for tip positions placed at different distances from the QD’s center. We now discuss our numerical TB calculations, which use the potential profiles from Figure 5a–d. These calculations allow us to simulate the QD’s LDOS in the presence of a fixed STM tip. Figure 5i shows the calculated LDOS distribution of a graphene QD when the STM tip is fixed at the center. In this image we note several distinct nodes that correspond to graphene QD states [26–30]. In Figure 5j,k we show the calculated LDOS distributions when the STM tip is fixed 50 nm and 100 nm away Figure 5. Simulation of the local density of states (LDOS) of an exposed graphene QD with an STMtip at a fixed position. (a–c) Spatial 2D map of graphene’s Dirac point energy with respect to EF afteradding the contributions from hBN defects (Figure 3c) and the STM tip (Figure 4b). The STM tip’slocation varies for each of these maps. (d) Spatial 2D map of the graphene Dirac point energy withrespect to EF determined by only considering the contribution from hBN defects (Figure 3c). (e–h)Respective potential profile line cuts along the yellow dashed lines in (a–d). (i–l) Calculated LDOSdistributions for each corresponding potential map. These simulations reveal the effect of the STM tipon the QD states for tip positions placed at different distances from the QD’s center.Nanomaterials 2020, 10, 1154 8 of 12We now discuss our numerical TB calculations, which use the potential profiles from Figure 5a–d.These calculations allow us to simulate the QD’s LDOS in the presence of a fixed STM tip. Figure 5ishows the calculated LDOS distribution of a graphene QD when the STM tip is fixed at the center. In thisimage we note several distinct nodes that correspond to graphene QD states [26–30]. In Figure 5j,kwe show the calculated LDOS distributions when the STM tip is fixed 50 nm and 100 nm away fromthe QD’s center. As a comparison we also show the calculated LDOS of a QD that excludes the effectof the STM tip (Figure 5l). When the tip is fixed at the QD’s center (Figure 5i) or 50 nm away fromthe center (Figure 5j), several new states with higher LDOS appear at the tip’s location. We also notethat the LDOS distribution and intensities away from the tip’s location (in Figure 5i,j) are similar tothe calculated LDOS that excludes the tip’s effect (Figure 5l). When the tip is fixed 100 nm awayfrom the QD’s center, we observe a new state with a much higher LDOS intensity at the tip’s position(Figure 5k). Similar to the previous cases, the LDOS distribution and intensities away from the tip’slocation remains unaffected.After demonstrating that graphene QD states are affected by an STM tip at a fixed position withour calculations, we study the case for a movable STM tip, which is more akin to our experiment.In measurements such as those shown in Figure 2c–e, each vertical array of pixels in the imagecorresponds to a dI/dVS curve acquired at the location of the STM tip. Consequently, to compare ourexperimental results with our simulations, we calculate the QD’s LDOS with the STM tip located ateach position along a line that crosses the QD. After obtaining the LDOS distribution from each profile(similar to those in Figure 5e–g), we compile the single dI/dVS curves calculated specifically at theSTM tip’s location for each point within the QD and along the line that crosses the QD. In Figure 6a weshow the result of this compilation process.We now consider how varying VG affects the QD states in our calculation. By changing VG theglobal electron and hole densities in graphene are offset. We simulate this effect by shifting the chargedensity profile due to the hBN defects (see Figure 3c). Specifically, this profile is shifted up to simulatean increase in hole density in graphene. Onto this shifted profile we add the unchanged STM tip’scharge density profile (Figure 4b) and perform the sequence of calculations as described for Figure 6a.Figure 6b,c show the complete results for two additional values of VG.To highlight the importance of the STM tip’s influence on our exposed graphene QDs, we showLDOS calculations at different VG values that omit the tip’s presence (Figure 6d–f). Clearly, bycomparing these simulations with measurements from Figure 2c–e and the simulations in Figure 6a–c,it is evident that a tip-inclusive model achieves better agreement with our experimental results. Similaragreement can be seen with other experimental results as well [14–16]. Moreover, the model that omitsthe STM tip (Figure 6d–f) lacks several key features from our experiment. For example, the deflectionof states near the edges of the QD are missing in Figs 6d–f. On the other hand, the tip-inclusive modelcaptures this feature consistently. In addition, we note that the presence of a continuous bright line thatwraps around the edge of the QD’s profile for higher VG values is shared by the experiment (Figure 2c)and the comprehensive model (Figure 6b), but absent in the model that ignores the tip. Finally, for thegraphene QD from the experiment we note that the QD states are less distorted at higher values ofVG, see for example Figure 2e. This insightful trend is also displayed in the comprehensive modelindicating that the effect of the tip can be mitigated when a sufficiently large QD is achieved.Nanomaterials 2020, 10, 1154 9 of 12Nanomaterials 2020, 10, x FOR PEER REVIEW 8 of 11  from the QD’s center. As a comparison we also show the calculated LDOS of a QD that excludes the effect of the STM tip (Figure 5l). When the tip is fixed at the QD’s center (Figure 5i) or 50 nm away from the center (Figure 5j), several new states with higher LDOS appear at the tip’s location. We also note that the LDOS distribution and intensities away from the tip’s location (in Figure 5i,j) are similar to the calculated LDOS that excludes the tip’s effect (Figure 5l). When the tip is fixed 100 nm away from the QD’s center, we observe a new state with a much higher LDOS intensity at the tip’s position (Figure 5k). Similar to the previous cases, the LDOS distribution and intensities away from the tip’s location remains unaffected. After demonstrating that graphene QD states are affected by an STM tip at a fixed position with our calculations, we study the case for a movable STM tip, which is more akin to our experiment. In measurements such as those shown in Figure 2c–e, each vertical array of pixels in the image corresponds to a /  curve acquired at the location of the STM tip. Consequently, to compare our experimental results with our simulations, we calculate the QD’s LDOS with the STM tip located at each position along a line that crosses the QD. After obtaining the LDOS distribution from each profile (similar to those in Figure 5e–g), we compile the single /  curves calculated specifically at the STM tip’s location for each point within the QD and along the line that crosses the QD. In Figure 6a we show the result of this compilation process.  Figure 6. Simulation of the LDOS in a graphene QD with a moving STM tip for different  configurations. (a–c) Simulated LDOS spectra along a line that crosses the center of a graphene QD. This simulation includes the influence of the STM tip. The LDOS spectra at each position is calculated with the tip fixed at that location. (d–f) Simulated LDOS spectra along a line that crosses the center of the graphene QD without considering the effect of the STM tip. The potentials induced by hBN defects used in (a–c) are the same as those used in (d–f), respectively. By comparing our results that include the effect of the STM tip (a–c) with those that exclude it (d–f), we find that a comprehensive treatment of the QD electrostatic environment is necessary to achieve agreement between theory and the experimental results from Figure 2c–e. Figure 6. Simulation of the LDOS in a graphene QD with a moving STM tip for different VGconfigurations. (a–c) Simulated LDOS spectra along a line that crosses the center of a graphene QD.This simulation includes the influence of the STM tip. The LDOS spectra at each position is calculatedwith the tip fixed at that location. (d–f) Simulated LDOS spectra along a line that crosses the centerof the graphene QD without considering the effect of the STM tip. The potentials induced by hBNdefects used in (a–c) are the same as those used in (d–f), respectively. By comparing our results thatinclude the effect of the STM tip (a–c) with those that exclude it (d–f), we find that a comprehensivetreatment of the QD electrostatic environment is necessary to achieve agreement between theory andthe experimental results from Figure 2c–e.5. ConclusionsIn conclusion, we showed that incorporating the STM tip’s electrostatics in conjunction with that ofthe underlying hBN charges enables a more complete understanding of the experimental spectroscopicfeatures of exposed graphene QDs. We compared experimental STM data obtained on grapheneQDs with simulations that include the tip-induced potential as well as with simulations that neglectthis potential. The agreement between experiments and simulations is greater when the simulationsinclude the influence of the tip. In particular, the experimentally observed bright envelope of thepotential and the deflection of states close to the QD edge are only reproduced when the tip-inducedpotential is included. Our results highlight the importance of considering the effect of the STM tipwhen interpreting spectroscopic characterization of exposed graphene QD states. Our analysis alsoreveals the intriguing possibility of studying the interplay between states confined by the potential dueto hBN defects and the potential due to the STM tip (see Figure 5k). Studies that seek to reduce suchinterplay may use insights from our simulations to mitigate the tip’s effect by tuning VG. Additionally,the interaction between these two QDs could potentially be used to emulate relativistic molecularbehavior or other complex coupled QD systems [37].Nanomaterials 2020, 10, 1154 10 of 12Author Contributions: Conceptualization, J.V.J., E.A.Q.-L., Z.G., and F.J.; Methodology, E.A.Q.-L. and Z.G.;Resources, K.W., T.T., and J.V.J.; Software, Z.G., E.A.Q.-L. and F.J.; Investigation, E.A.Q.-L., F.J. and Z.G.; DataCuration, E.A.Q.-L., F.J. and Z.G.; Formal Analysis, Z.G., E.A.Q.-L. and J.V.J.; Original Draft Preparation, E.A.Q.-L.and Z.G.; Review & Editing, J.V.J., E.A.Q.-L. and Z.G.; Supervision, F.J. and J.V.J.; Project Administration, J.V.J.;Funding Acquisition, J.V.J. All authors have read and agreed to the published version of the manuscript.Funding: This research was funded by the National Science Foundation, DMR-1753367, and the Army ResearchOffice, W911NF-17-1-0473, and the Elemental Strategy Initiative conducted by the MEXT, Japan, JPMXP0112101001,and the JSPS KAKENHI, JP20H00354, and the CREST, JPMJCR15F3.Acknowledgments: We thank the Hummingbird Computational Cluster team at UC Santa Cruz for providingcomputational resources and support for the numerical tight binding calculations performed in this work.All LDOS calculations were performed using the Pybinding package [38].Conflicts of Interest: The authors declare no conflict of interest.References1. Ponomarenko, L.A.; Schedin, F.; Katsnelson, M.I.; Yang, R.; Hill, E.W.; Novoselov, K.S.; Geim, A.K. ChaoticDirac Billiard in Graphene Quantum Dots. Science 2008, 320, 356–358. [CrossRef] [PubMed]2. Schnez, S.; Güttinger, J.; Huefner, M.; Stampfer, C.; Ensslin, K.; Ihn, T. Imaging Localized States in GrapheneNanostructures. Phys. Rev. B 2010, 82, 165445. [CrossRef]3. Todd, K.; Chou, H.T.; Amasha, S.; David, G.G. Quantum Dot Behavior in Graphene Nanoconstrictions.Nano Lett. 2009, 9, 416–421. [CrossRef]4. Han, M.Y.; Özyilmaz, B.; Zhang, Y.; Kim, P. Energy Band-Gap Engineering of Graphene Nanoribbons.Phys. Rev. Lett. 2007, 98, 206805. [CrossRef]5. Subramaniam, D.; Libisch, F.; Li, Y.; Pauly, C.; Geringer, V.; Reiter, R.; Mashoff, T.; Liebmann, M.; Burgdörfer, J.;Busse, C.; et al. Wave-Function Mapping of Graphene Quantum Dots with Soft Confinement. Phys. Rev. Lett.2012, 108, 046801. [CrossRef] [PubMed]6. Hämäläinen, S.K.; Sun, Z.; Boneschanscher, M.P.; Uppstu, A.; Ijäs, M.; Harju, A.; Vanmaekelbergh, D.;Liljeroth, P. Quantum-Confined Electronic States in Atomically Well-Defined Graphene Nanostructures.Phys. Rev. Lett. 2011, 107, 236803. [CrossRef]7. Phark, S.H.; Borme, J.; Vanegas, A.L.; Corbetta, M.; Sander, D.; Kirschner, J. Direct Observation of ElectronConfinement in Epitaxial Graphene Nanoislands. ACS Nano 2011, 5, 8162–8166. [CrossRef]8. Lu, J.; Yeo, P.S.E.; Gan, C.K.; Wu, P.; Loh, K.P. Transforming C 60 Molecules into Graphene Quantum Dots.Nat. Nanotechnol. 2011, 6, 247–252. [CrossRef]9. Wang, Y.; Wong, D.; Shytov, A.V.; Brar, V.W.; Choi, S.; Wu, Q.; Tsai, H.Z.; Regan, W.; Zettl, A.; Kawakami, R.K.;et al. Observing Atomic Collapse Resonances in Artificial Nuclei on Graphene. Science 2013, 340, 734–737.[CrossRef]10. Jung, S.; Rutter, G.M.; Klimov, N.N.; Newell, D.B.; Calizo, I.; Hight-Walker, A.R.; Zhitenev, N.B.; Stroscio, J.A.Evolution of Microscopic Localization in Graphene in a Magnetic Field from Scattering Resonances toQuantum Dots. Nat. Phys. 2011, 7, 245–251. [CrossRef]11. Freitag, N.M.; Chizhova, L.A.; Nemes-Incze, P.; Woods, C.R.; Gorbachev, R.V.; Cao, Y.; Geim, A.K.;Novoselov, K.S.; Burgdörfer, J.; Libisch, F.; et al. Electrostatically Confined Monolayer Graphene QuantumDots with Orbital and Valley Splittings. Nano Lett. 2016, 16, 5798–5805. [CrossRef] [PubMed]12. Freitag, N.M.; Reisch, T.; Chizhova, L.A.; Nemes-Incze, P.; Holl, C.; Woods, C.R.; Gorbachev, R.V.; Cao, Y.;Geim, A.K.; Novoselov, K.S.; et al. Large Tunable Valley Splitting in Edge-Free Graphene Quantum Dots onBoron Nitride. Nat. Nanotechnol. 2018, 13, 392–397. [CrossRef] [PubMed]13. Velasco, J.; Ju, L.; Wong, D.; Kahn, S.; Lee, J.; Tsai, H.Z.; Germany, C.; Wickenburg, S.; Lu, J.; Taniguchi, T.; et al.Nanoscale Control of Rewriteable Doping Patterns in Pristine Graphene/Boron Nitride Heterostructures.Nano Lett. 2016, 16, 1620–1625. [CrossRef] [PubMed]14. Lee, J.; Wong, D.; Velasco, J.; Rodriguez-Nieva, J.F.; Kahn, S.; Tsai, H.Z.; Taniguchi, T.; Watanabe, K.; Zettl, A.;Wang, F.; et al. Imaging Electrostatically Confined Dirac Fermions in Graphene Quantum Dots. Nat. Phys.2016, 12, 1032–1036. [CrossRef]http://dx.doi.org/10.1126/science.1154663http://www.ncbi.nlm.nih.gov/pubmed/18420930http://dx.doi.org/10.1103/PhysRevB.82.165445http://dx.doi.org/10.1021/nl803291bhttp://dx.doi.org/10.1103/PhysRevLett.98.206805http://dx.doi.org/10.1103/PhysRevLett.108.046801http://www.ncbi.nlm.nih.gov/pubmed/22400872http://dx.doi.org/10.1103/PhysRevLett.107.236803http://dx.doi.org/10.1021/nn2028105http://dx.doi.org/10.1038/nnano.2011.30http://dx.doi.org/10.1126/science.1234320http://dx.doi.org/10.1038/nphys1866http://dx.doi.org/10.1021/acs.nanolett.6b02548http://www.ncbi.nlm.nih.gov/pubmed/27466881http://dx.doi.org/10.1038/s41565-018-0080-8http://www.ncbi.nlm.nih.gov/pubmed/29556008http://dx.doi.org/10.1021/acs.nanolett.5b04441http://www.ncbi.nlm.nih.gov/pubmed/26852622http://dx.doi.org/10.1038/nphys3805Nanomaterials 2020, 10, 1154 11 of 1215. Ghahari, F.; Walkup, D.; Gutiérrez, C.; Rodriguez-Nieva, J.F.; Zhao, Y.; Wyrick, J.; Natterer, F.D.; Cullen, W.G.;Watanabe, K.; Taniguchi, T.; et al. An on/off Berry Phase Switch in Circular Graphene Resonators. Science2017, 356, 845–849. [CrossRef]16. Gutiérrez, C.; Walkup, D.; Ghahari, F.; Lewandowski, C.; Rodriguez-Nieva, J.F.; Watanabe, K.; Taniguchi, T.;Levitov, L.S.; Zhitenev, N.B.; Stroscio, J.A. Interaction-Driven Quantum Hall Wedding Cake–like Structuresin Graphene Quantum Dots. Science 2018, 361, 789–794. [CrossRef]17. Zhao, Y.; Wyrick, J.; Natterer, F.D.; Rodriguez-Nieva, J.F.; Lewandowski, C.; Watanabe, K.; Taniguchi, T.;Levitov, L.S.; Zhitenev, N.B.; Stroscio, J.A. Creating and Probing Electron Whispering-Gallery Modes inGraphene. Science 2015, 348, 672–675. [CrossRef]18. Jiang, Y.; Mao, J.; Moldovan, D.; Masir, M.R.; Li, G.; Watanabe, K.; Taniguchi, T.; Peeters, F.M.; Andrei, E.Y.Tuning a Circular P-N Junction in Graphene from Quantum Confinement to Optical Guiding. Nat. Nanotechnol.2017, 12, 1045–1049. [CrossRef]19. Zomer, P.J.; Dash, S.P.; Tombros, N.; Van Wees, B.J. A Transfer Technique for High Mobility Graphene Deviceson Commercially Available Hexagonal Boron Nitride. Appl. Phys. Lett. 2011, 99, 232104. [CrossRef]20. Goossens, A.M.; Calado, V.E.; Barreiro, A.; Watanabe, K.; Taniguchi, T.; Vandersypen, L.M.K. MechanicalCleaning of Graphene. Appl. Phys. Lett. 2012, 100, 073110. [CrossRef]21. Velasco, J.; Lee, J.; Wong, D.; Kahn, S.; Tsai, H.Z.; Costello, J.; Umeda, T.; Taniguchi, T.; Watanabe, K.; Zettl, A.;et al. Visualization and Control of Single-Electron Charging in Bilayer Graphene Quantum Dots. Nano Lett.2018, 18, 5104–5110. [CrossRef] [PubMed]22. Tersoff, J.; Hamann, D.R. Theory of the Scanning Tunneling Microscope. Phys. Rev. B 1985, 31, 805–813.[CrossRef] [PubMed]23. Katsnelson, M.I.; Novoselov, K.S.; Geim, A.K. Chiral Tunnelling and the Klein Paradox in Graphene. Nat. Phys.2006, 2, 620–625. [CrossRef]24. Cheianov, V.V.; Fal’ko, V.I. Selective Transmission of Dirac Electrons and Ballistic Magnetoresistance of N-pJunctions in Graphene. Phys. Rev. B 2006, 74, 041403. [CrossRef]25. Shytov, A.V.; Rudner, M.S.; Levitov, L.S. Klein Backscattering and Fabry-Pérot Interference in GrapheneHeterojunctions. Phys. Rev. Lett. 2008, 101, 156804. [CrossRef] [PubMed]26. Wu, J.S.; Fogler, M.M. Scattering of Two-Dimensional Massless Dirac Electrons by a Circular Potential Barrier.Phys. Rev. B 2014, 90, 235402. [CrossRef]27. Schulz, C.; Heinisch, R.L.; Fehske, H. Scattering of Two-Dimensional Dirac Fermions on Gate-DefinedOscillating Quantum Dots. Phys. Rev. B 2015, 91, 045130. [CrossRef]28. Chen, H.Y.; Apalkov, V.; Chakraborty, T. Fock-Darwin States of Dirac Electrons in Graphene-Based ArtificialAtoms. Phys. Rev. Lett. 2007, 98, 186803. [CrossRef]29. Matulis, A.; Peeters, F.M. Quasibound States of Quantum Dots in Single and Bilayer Graphene. Phys. Rev. B2008, 77, 115423. [CrossRef]30. Bardarson, J.H.; Titov, M.; Brouwer, P.W. Electrostatic Confinement of Electrons in an Integrable GrapheneQuantum Dot. Phys. Rev. Lett. 2009, 102, 226803. [CrossRef]31. Zhang, Y.; Brar, V.W.; Girit, C.; Zettl, A.; Crommie, M.F. Origin of Spatial Charge Inhomogeneity in Graphene.Nat. Phys. 2009, 5, 722–726. [CrossRef]32. Rutter, G.M.; Jung, S.; Klimov, N.N.; Newell, D.B.; Zhitenev, N.B.; Stroscio, J.A. Microscopic Polarization inBilayer Graphene. Nat. Phys. 2011, 7, 649–655. [CrossRef]33. LeRoy, B.J.; Kong, J.; Pahilwani, V.K.; Dekker, C.; Lemay, S.G. Three-Terminal Scanning Tunneling Spectroscopyof Suspended Carbon Nanotubes. Phys. Rev. B 2005, 72, 075413. [CrossRef]34. Chae, J.; Jung, S.; Young, A.F.; Dean, C.R.; Wang, L.; Gao, Y.; Watanabe, K.; Taniguchi, T.; Hone, J.;Shepard, K.L.; et al. Renormalization of the Graphene Dispersion Velocity Determined from ScanningTunneling Spectroscopy. Phys. Rev. Lett. 2012, 109, 116802. [CrossRef] [PubMed]35. Zhang, Y.; Brar, V.W.; Wang, F.; Girit, C.; Yayon, Y.; Panlasigui, M.; Zettl, A.; Crommie, M.F. GiantPhonon-Induced Conductance in Scanning Tunnelling Spectroscopy of Gate-Tunable Graphene. Nat. Phys.2008, 4, 627–630. [CrossRef]36. Holdman, G.R.; Krebs, Z.J.; Behn, W.A.; Smith, K.J.; Watanabe, K.; Taniguchi, T.; Brar, V.W. Dynamic BandStructure and Capacitance Effects in Scanning Tunneling Spectroscopy of Bilayer Graphene. Appl. Phys. Lett.2019, 115, 181601. [CrossRef]http://dx.doi.org/10.1126/science.aal0212http://dx.doi.org/10.1126/science.aar2014http://dx.doi.org/10.1126/science.aaa7469http://dx.doi.org/10.1038/nnano.2017.181http://dx.doi.org/10.1063/1.3665405http://dx.doi.org/10.1063/1.3685504http://dx.doi.org/10.1021/acs.nanolett.8b01972http://www.ncbi.nlm.nih.gov/pubmed/30035544http://dx.doi.org/10.1103/PhysRevB.31.805http://www.ncbi.nlm.nih.gov/pubmed/9935822http://dx.doi.org/10.1038/nphys384http://dx.doi.org/10.1103/PhysRevB.74.041403http://dx.doi.org/10.1103/PhysRevLett.101.156804http://www.ncbi.nlm.nih.gov/pubmed/18999625http://dx.doi.org/10.1103/PhysRevB.90.235402http://dx.doi.org/10.1103/PhysRevB.91.045130http://dx.doi.org/10.1103/PhysRevLett.98.186803http://dx.doi.org/10.1103/PhysRevB.77.115423http://dx.doi.org/10.1103/PhysRevLett.102.226803http://dx.doi.org/10.1038/nphys1365http://dx.doi.org/10.1038/nphys1988http://dx.doi.org/10.1103/PhysRevB.72.075413http://dx.doi.org/10.1103/PhysRevLett.109.116802http://www.ncbi.nlm.nih.gov/pubmed/23005662http://dx.doi.org/10.1038/nphys1022http://dx.doi.org/10.1063/1.5127078Nanomaterials 2020, 10, 1154 12 of 1237. Walkup, D.; Ghahari, F.; Gutiérrez, C.; Watanabe, K.; Taniguchi, T.; Zhitenev, N.B.; Stroscio, J.A. TuningSingle-Electron Charging and Interactions between Compressible Landau Level Islands in Graphene.Phys. Rev. B 2020, 101, 035428. [CrossRef]38. Moldovan, D.; And̄elković, M.; Peeters, F. pybinding v0.9.4: A Python package for tight-binding calculations.Zenodo 2017. [CrossRef]© 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open accessarticle distributed under the terms and conditions of the Creative Commons Attribution(CC BY) license (http://creativecommons.org/licenses/by/4.0/).http://dx.doi.org/10.1103/PhysRevB.101.035428http://dx.doi.org/10.5281/zenodo.826942http://creativecommons.org/http://creativecommons.org/licenses/by/4.0/. Introduction  Materials and Methods  Results  Discussion  Conclusions  References