Critical Behavior of Four-Terminal Conductance of Bilayer Graphene Domain Walls Benjamin J. Wieder,1 Fan Zhang,2 and C. L. Kane1 arXiv:1409.7645v2 [cond-mat.mes-hall] 5 Aug 2015 1 Department of Physics and Astronomy, University of Pennsylvania, Philadelphia, PA 19104 2 Department of Physics, University of Texas at Dallas, Richardson, TX 75080 Bilayer graphene in a perpendicular electric field can host domain walls between regions of reversed field direction or interlayer stacking. The gapless modes propagating along these domain walls, while not strictly topological, nevertheless have interesting physical properties, including valley-momentum locking. A junction where two domain walls intersect forms the analogue of a quantum point contact. We study theoretically the critical behavior of this junction near the pinchoff transition, which is controlled by two separate classes of non-trivial quantum critical points. For strong interactions, the junction can host phases of unique charge and valley conductances. For weaker interactions, the low-temperature charge conductance can undergo one of two possible quantum phase transitions, each characterized by a specific critical exponent and a collapse to a universal scaling function, which we compute. PACS numbers: 74.45.+c, 71.10.Pm, 74.78.Fk, 74.78.Na I. INTRODUCTION A. Background Bilayer graphene1,2 provides an attractive platform for unconventional two-dimensional electronic physics due to the two quadratic band contacts at its Fermi points, and because of the variety of ways by which one can introduce a band gap and produce momentum-space Berry curvature3 . The interlayer nearest-neighbor hopping, γ1 , warps the band structure of the individual graphene layers, repelling two bands away from the Fermi energy and leaving the remaining two dispersing quadratically. This warping is a consequence of the two-step process in which electrons hop between the two low-energy sublattices via the two high-energy ones. The band touching points at inequivalent Brillouin zone corners K or K are protected by the Berry phase ±2π (or winding number ±2) and the required chiral (or sublattice) symmetry between the low-energy sublattices on opposite layers4 . Keeping all of these symmetries, the band touching point can only split, instead of being gapped out, even when trigonal warping and other weak remote-hopping processes are taken into account2 . However, the chiral symmetry between the low-energy sublattices can be intrinsically broken by spin-orbit coupling5 , spontaneously broken by electron-electron interactions3,6,7 , or explicitly broken by an interlayer potential difference2,8–10 . As a consequence, the quadratic band touching is no longer symmetry-protected and gaps open up at valleys K and K . While the first two types of gaps are small in practice4,11–13 , the electric-field-induced gap, which is the focus of this paper, saturates at a large value comparable to the interlayer hopping γ1 ∼ 0.3 eV10,14 . Opening the band gap produces large momentum-space Berry curvature in the quasiparticle bands, with the curvature integral quantized to ±1 over a half Brillouin zone centered at K or K 3,15 . Moreover, for bilayer graphene gapped by an electric field, the sign of this partial Chern number depends on the valley index, the sign of the energy gap (given by the direction of interlayer electric field), and the layer stacking order (i.e., AB or BA)3,16 . Here AB (BA) stacking refers to the case in which γ1 couples the top layer A (B) and bottom layer B (A) highenergy sublattices. In the presence of an interlayer electric field, when the field direction is reversed across a line15–21 or when the field is uniform and the layer stacking switches from AB to BA16,22,23 , the valley-projected Chern number changes by 2 (−2) across the domain wall in valley K (K ). As a result, both types of domain walls host two chiral edge states in each valley with chirality (direction) locked to valley index K or K , as shown in Fig. 1. Similar domain wall states also occur spontaneously due to interactions in the absence of electric fields but at finite temperature24 . Importantly16 , these “Quantum Valley Hall” (QVH) edge states are not strictly topological and can be gapped out by a sufficiently strong, large-momentum scattering which couples the two valleys, even if the underlying symmetries are still preserved. It is therefore crucial that valley index also remains a “good quantum number,” for which we will assume that short-range disorder, interlayer stacking, and electric field direction changes are smooth on the scale of the lattice. Under this assumption, backscattering is prohibited and the system of domain walls provides an attractive platform for Tomonaga-Luttinger liquid physics25 . In this paper, we study the electronic transport properties of a junction where two domain walls intersect (Fig. 2). Such a structure resembles the quantum point contact of the edges of two Quantum Spin Hall (QSH) insulators, which has been studied in Refs. 26–28 and can be probed by four-terminal transport measurements. A domain wall junction can be tuned through a “pinch-off transition” by applying a local field (such as a perpendicular electric field) to the junction region. In Ref. 28, it was found that the corresponding pinch-off transition for QSH systems is controlled by a novel quantum critical point, and that at low temperatures the conductance is 2 described by a universal scaling function across the pinchoff transition. In contrast to the QSH edge states, which have a single time-reversed pair of helical modes, the domain wall states in bilayer graphene host four helical pairs (including electron spin). We find that this leads to several important modifications of the low-energy properties of the junction. Unlike with the QSH edge states, whose forwardscattering interactions are characterized by a single Luttinger parameter, it has been argued that for the domain wall states in bilayer graphene, one should characterize interactions with two independent Luttinger parameters25 . This leads to an expanded phase diagram for the possible stable states of the junction. Moreover, we find that the pinch-off transition is modified. Depending on the interaction strengths, there are two possible regimes for the reduced, two-terminal conductance: one in which it undergoes a single pinch-off transition directly from 0 to 8e2 /h and one in which it undergoes two separate transitions, separated by a stable state with conductance 4e2 /h. We study the critical behavior of these transitions and compute the universal crossover scaling functions for weak interactions. This paper is structured as follows. First we introduce in detail the domain wall states in bilayer graphene and derive low-energy effective field theories for them. Adding interactions, we show how these states are Luttinger liquids described by two independent Luttinger parameters. From there, we characterize the geometry of two intersecting domain wall states in the language of the resulting effective quantum point contact. We then analyze the resulting four-terminal junction in both the context of many-body tunneling in a bosonization framework and with an S-matrix renomalization group using diagrammatic perturbation theory. Combining these analyses, we determine the behavior of the reduced, two-terminal conductance over a range of interaction strengths. B. Measurable Results tance GXX undergoes either a direct transition from 0 to 8e2 /h or an indirect one with an intermediate step to 4e2 /h depending on experimental specifics (Fig. 9). Labeling the direct transition A and the first step of the indirect transition B, we go on to show that at low temperatures, the conductance transitions should collapse onto universal scaling functions GA/B with critical exponents αA/B as functions of the interaction strengths (Figs. 10,11). II. MODEL SYSTEM In this section we introduce our model system of bilayer graphene domain wall modes. First, we begin with the Hamiltonian for a single domain wall and the Luttinger liquid physics which govern it in the presence of interactions. Then, we discuss the four-terminal geometry which arises at the intersection of two domain walls and its equivalence to a quantum point contact. A. Domain Walls in Bilayer Graphene As discussed in the introduction, bilayer graphene domain walls can be created by varying the direction of the interlayer electric field or by varying the interlayer stacking order. In either case, The valley-projected Chern number changes by 2 (−2) across the domain wall in valley K (K ). This necessitates the existence of two domain wall states in each valley, with the states at K having equal and opposite velocities to those at K’. Adopting the notation of Ref. 25, we label these bands 0 and π respectively (Fig. 1). For the purposes of our model, we will assume that the Fermi energy lies exactly in the middle of the bulk gap, which allows the simplification vF,0 = vF,π = vF . This allows us to write down the non-interacting Hamiltonian density † † ψaσ,in ∂x ψaσ,in − ψaσ,out ∂x ψaσ,out H0 = i¯ vF h a=0,π σ=↑,↓ In this paper, we calculate several measurable properties of bilayer graphene domain wall quantum point contacts. In section II A, we find the critical exponent αT which characterizes the low temperature tunneling conductance scaling for a single domain wall, a result previously derived in Ref. 25. In II B, we introduce a diagonal conductance GZZ = 8e2 /h, which is only strictly quantized when valley index is conserved both within individual domain walls and across the junction. This conductance, therefore, stands as a first test of whether experimental samples are in the appropriate disorder regime for the analysis in this paper. We also show in II B that states of exotic charge and valley tunneling character dominate the conductance of the junction under very strong attractive or repulsive interactions (Fig. 4). Finally, in III B 3, we show that the left-to-right conduc- (1) where the indexes “in” and “out” refer to direction and can correspond to electron operators with valley index K or K depending on the orientation of the domain wall. We will see later that in a four-terminal structure, the correspondence between in/out and K/K will alternate with lead index due to the valley-momentum locking of the domain wall modes. In the limit of the Fermi energy resting exactly at the chiral-symmetric point, the band indexes become arbitrary labels for all of the calculations in this paper. A consequence of this is the emergence of a band-indexexchange symmetry in this problem, which we will frequently highlight in calculations throughout this paper. Deviations from this point in the Fermi energy are expected in physical systems, and will lead in general to a 3 (a) V+ = u+ (n0↑ + n0↓ + nπ↑ + nπ↓ )2 V− = u− (n0↑ + n0↓ − nπ↑ − nπ↓ )2 . (b) (c) (2) V+ is the usual two-body forward scattering term which leads to Luttinger liquid physics and V− is a new one which breaks the U(2) symmetry of electron distribution between the 0 and π bands. Both can be effectively tuned by altering the strength of the perpendicular electric field, though for all reasonable numerical estimates, Ref. 25 found u− < u+ and u− harder to tune, which is sensible as only V+ contains contributions from the long-range part of the Coulomb interaction. Other density-density interactions, specifically those which affect electron spin, should be small in practice. In the absence of an external magnetic field, and given the weak spin-orbit interaction in graphene, electron spin should remain SU(2)-symmetric. In this limit, the electron spin sectors of this problem should remain noninteracting and terms which lead to phenomena such as spin-density waves will be marginally irrelevant29 . Returning to our Hamiltonian, we can consider a single domain wall by bosonizing, ψaσ,i = √ 1 eiφaσ,i 2πxc (3) where a = 0, π; σ =↑, ↓; i = in, out; and xc is the short wavelength cutoff. The bosonic fields φaσ,i obey the commutation algebra: FIG. 1. Domain walls in bilayer graphene can be induced by applying a perpendicular electric field and varying either the interlayer stacking (a) or the electric field direction (b). Both kinds of domain walls (the dotted lines) have similar domain wall band structures (c) when the Fermi energy EF is near the chiral symmetric point. Adopting the notation of Ref. 25, the two domain wall states in each direction are labeled 0 and π working from the Brillouin zone edge in. When the Fermi energy is exactly in the middle of the bulk gap, the Fermi velocities are the same for the 0 and π bands and electron direction is set by valley index K/K . z [φaσ,i (x), φbσ ,j (y)] = iπδab δσσ τij sgn(x − y). (4) Under this transformation, the bare Hamiltonian and interactions become: H0 = h ¯ vF 4π (∂x φ0σ,in )2 + (∂x φ0σ,out )2 σ=↑,↓ + (∂x φπσ,in )2 + (∂x φπσ,out )2 h ¯ vF V± = λ± (∂x φ0σ,in − ∂x φ0σ,out ) 8π σ=↑,↓ ± (∂x φπσ,in − ∂x φπσ,out ) relaxation of this symmetry. If the deviations are small, the physics should resemble the predictions of this paper, with corrections of order ∼ 1 − vF,0 /vF,π . These corrections greatly complicate the calculations in this paper and obscure key generalities, and to that end we consider vF,0 = vF,π = vF a desirable simplification of this problem. Calculations by Killi, Wei, Affleck, and Paramekanti indicated that for interacting bilayer graphene systems, the interaction is dominated by two density-density interactions25 : 2 (5) where λ± = u± /π¯ vF . The interacting Hamiltonian h can be simplified by the sum/difference changes of basis: φ±σ,i = φ0σ,i ± φπσ,i φ±c/s,i = φ±↑,i ± φ±↓,i (6) where the c and s sectors are charge and spin respectively. In this basis, all of the interactions are in the charge sector and, as motivated earlier in this section, the spin sector is noninteracting: 4 H0 = H+c + H−c + H+s + H−s h ¯ vF (∂x φ±c/s,in )2 + (∂x φ±c/s,out )2 H±c/s = 8π h ¯ vF 2 V± = λ±c [∂x φ±c,in − ∂x φ±,out ] . 8π (a) 1 (7) The plus and minus charge sectors of H0 are then each renormalized by only V+/− respectively, encouraging us to express the interaction parameter g+/− separately for each charge sector. Therefore we can write down the interacting Hamiltonian for each charge sector, H±c,int = H±c + V± . Diagonalizing this Hamiltonian, the definition of the Luttinger parameters g± arises naturally. The change of basis φ±c,in φ±c,out = 1 2g± 1 + g± 1 − g± 1 − g± 1 + g± ˜ φ±c,in ˜±c,out φ (8) returns our interacting Hamiltonian to the form of one for non-interacting chiral bosons Hint = H+c,int + H−c,int + H+s + H−s h ¯ v± ˜ ˜ (∂x φ±c,in )2 + (∂x φ±c,out )2 H±c,int = 8πg± 1 + 2λ± , g± = 1 . 1 + 2λ± 2 1 4 B 3 A 2 A 4 QPC A B 3 FIG. 2. (a) Two parallel domain walls in bilayer graphene can be created by varying either the interlayer stacking or the perpendicular field direction between regions A and B. (b) Distorting region B such that the walls approach each other results in the equivalent of a Quantum Point Contact (QPC) for the domain wall modes. The numbers 1 − 4 are lead indexes and the two modes displayed for each domain wall are those at 0 and at π. All of the modes shown here are at valley K; a counterpropagating set of modes exists at K and is related by time-reversal symmetry. Including electron spin, there are 4 modes in each valley in each domain wall, for a total of 16 modes to consider for this QPC structure. Four-Terminal Geometry (9) (10) ˜ In the new basis, the charge fields φ±c,i obey the commutation relation: z ˜ ˜ φuc,i (x), φvc,j (y) = iπgu δuv τij sgn(x − y) (b) B B. where v± = vF A (11) where u, v = +, −; i = in, out; and we note that the noninteracting spin sector fields still obey this commutation relation with g± = 1. As in Refs. 26–28, the tunneling density of states for a single edge ρ(E) ∝ E αT is controlled by the interactions. However here, unlike in the QSH case, the critical exponent is a function of two Luttinger parameters, such that in agreement with Ref. 25, A pattern of two domain walls which pass nearly by each other can be formed by varying either electric field direction or interlayer stacking twice (Fig. 2a). If we distort the central region of this picture, we could imagine bringing the two domain walls so close that tunneling between them becomes significant. In this case, the two domain walls have formed the equivalent of a quantum point contact (Fig. 2b). Mapping our work in II A onto a spinful Luttinger liquid for the two-wall system, we deduce in this section, for all interaction regimes, the interaction strengths at which many-body tunneling processes become relevant and alter the junction’s charge and valley conductances. (12) While each lead retains the properties and Luttinger parameters from II individually, we will find it conve˜ nient to limit the usage of the φ± basis to the treatment of isolated domain walls and adopt a new basis with effective charge and valley sectors. The charge sector which arises here is a new degree of freedom and comes from a rotation of the indexes for K and K and propagation direction. We will label these charge and valley sectors ρ and v respectively to differentiate them from the charge and spin sectors which arose in II A, which are labeled c and s respectively. From an experimental perspective, measuring this critical exponent for the tunneling conductance would be a valuable first step in confirming the Luttinger liquid physics of these bilayer graphene domain wall states. The two domain walls in Fig. 2 have opposite helicties, due to being on the top (bottom) of the central region. We can define, for the interacting system, fields labeled by sum/difference, charge/spin, direction, and valley (K, K ): αT = 1 1 (g+ + g− + 1/g+ + 1/g− ) − . 8 2 5 φ±c/s,RK φ±c/s,LK φ±c/s,LK φ±c/s,RK = φ±c/s,in,1 (−x)Θ(−x) + φ±c/s,out,2 (x)Θ(x) = φ±c/s,out,1 (−x)Θ(−x) + φ±c/s,in,2 (x)Θ(x) = φ±c/s,in,3 (x)Θ(x) + φ±c/s,out,4 (−x)Θ(−x) = φ±c/s,out,3 (x)Θ(x) + φ±c/s,in,4 (−x)Θ(−x) (13) where Θ(x) is the Heaviside step function and the indexes 1 − 4 on the noninteracting φ± refer to the individual lead in Figure 2 with which they are associated. The intersection of the two domain walls occurs at x = 0 such that our theory consists of four, isolated domain walls everywhere except at that point. Assuming that the interaction strength is controlled globally such that each lead has the same values of λ± , we can create a ρ/v basis for each +/− sector of the combined Hamiltonian of the two 4 i i domain walls Hint = i=1 α=c,s (H+α,int + H−α,int ): 1 2 1 = 2 1 = 2 1 = 2 φ±c/s,RK = φ±c/s,ρ + φ±c/s,v − θ±c/s,ρ − θ±c/s,v H±,int = h ¯ v± 8π g±αa (∂x φ±αa )2 + α=c,s a=ρ,v φ±c/s,RK φ±c/s,LK φ±c/s,ρ − φ±c/s,v + θ±c/s,ρ − θ±c/s,v φ±c/s,ρ − φ±c/s,v − θ±c/s,ρ + θ±c/s,v (14) where again c, s are the charge and spin sectors which resulted from rotating the indexes for ↑, ↓ and ρ, v are the indexes which have, along with the choice of φ, θ, resulted from rotating indexes for propagation direction (R, L) and valley (K, K ). The new fields are governed by the modified commutation relation: [θuαi (x), φvβj (y)] = 2πiδuv δαβ δij Θ(x − y) (15) where u, v = +, −; α, β = c, s; and i, j = ρ, v. This unitary rotation of the variables effectively changes the sign of the interaction “cross-term” individually for the choice of φ, θ, ρ, and v within the c sector: H±c,int = h ¯ vF 8π (1 + λ± ) (∂x φ±cρ )2 + (∂x φ±cv )2 + (∂x θ±cρ )2 + (∂x θ±cv )2 − λ± (∂x φ±cρ )(∂x φ±cρ ) − (∂x φ±cv )(∂x φ±cv ) − (∂x θ±cρ )(∂x θ±cρ ) + (∂x θ±cv )(∂x θ±cv ) . (16) The previous equation, though diagonal, was left unsimplified and in the form of Eq. 5 such that by the same 1 g±αa g±cρ = g± , g±cv = 1/g± , g±sρ = g±sv = 1 (∂x θ±αa )2 (17) such that Hint now has the form of a spinful Luttinger liquid. In this basis, both the interacting and noninteracting Hamiltonians are diagonal and so the transformation between the interacting and noninteracting θ/φ requires just a simple rescaling by the interaction parameter for each sector. For this geometry, one can probe experimentally by measuring the current Ii at one of the leads in response to an applied voltage on another lead Vj such that a 4×4 conductance matrix characterizes the system, φ±c/s,ρ + φ±c/s,v + θ±c/s,ρ + θ±c/s,v φ±c/s,LK logic as in II A, the form of the simplified diagonalized Hamiltonian, as well as the interactions, can just be read off: Ii = Gij Vj (18) where i = 1 − 4 is a lead index. In the presence of time-reversal and valley symmetries, the number of independent or nonzero parameters in Gij is greatly reduced, as described in detail in the appendix of Ref. 28. For this system, we can then consider a reduced set of voltages and currents: IX IY = GXX GXY GY X GY Y VX VY (19) where IX = I1 + I4 is the left-to-right current and IY = I1 + I2 is the top-to-bottom current. VX and VY are similarly defined such that VX is a bias of leads 1 and 4 relative to leads 2 and 3 and VY is a bias of leads 1 and 2 relative to leads 3 and 4. Therefore GXX and GY Y are the two-terminal conductances measured left-to-right and top-to-bottom respectively. GXY = GY X are skew conductances, equal as a consequence of time-reversal symmetry. In the noninteracting model, this skew conductance is zero as a consequence of artificial spatial symmetries, such as mirror symmetry. Though it may become nonzero under increased interaction strengths, the skew conductance is still negligible along the relevant directions which characterize transitions in this system28 . We can define a final current across the junction IZ = I1 + I3 , which one can probe by applying a voltage VZ which biases leads 1 and 3 relative to leads 2 and 4, with a conductance IZ = GZZ VZ . (20) If valley is conserved, then electrons cannot enter at lead 1 and exit at lead 3, implying that a measurement of an exactly quantized 6 GZZ = 8e2 h (21) would be an experimental confirmation that valleynonconserving disorder is absent and the system is appropriately described by the physics in this paper. The factor of N = 8 in the Landauer prediction G = N e2 /h comes from factors of 2 for band index (0 and π), electron spin degeneracy, and the two incoming leads at K (1 and 3). We can also, in a similar manner, characterize the valley conductance of the system in terms of left-to-right and top-to-bottom parameters GV and GV Y . Before XX Y the introduction of any interactions or tunneling operators, our system consists of two, left-to-right domain wall states and we consider it “fully-open.” For this system, GXX = 8e2 /h and GV = 0 such that it is a left-to-right XX charge conductor, valley conductor, which we will denote as the CC phase. A 90◦ rotation and relabeling (with regards to Fig. 2) or the pinch-off inversion of this phase, for which GXX = GV = 0 and GY Y = 8e2 /h, GV Y = 0, is Y XX considered “fully pinched-off” and is a left-to-right charge insulator, valley insulator, which we denote as the II phase. With this framework established, we can examine perturbatively tunneling processes between the two adjacent domain walls which may lead to differing charge and valley conductances. Using our bosonization work, we can examine the rescaling of the coupling strength for each process, noting the interaction regime in which it dominates the physics of the quantum point contact. Figure 3 illustrates the single-particle valleyconserving processes which can exist in this system within a single spin channel. In general, many-body tunneling processes will also be present and may become relevant. These many-body tunneling processes can be considered products of single-electron-tunneling processes which, in the most general case, may or may not conserve spin or valley indexes. However, restricting ourselves to the set of processes which conserve valley, it becomes apparent that the linear combinations of bosonized operators which can lead to relevant operators can only be achieved through products of single-particle processes which conserve spin. Therefore, in the analysis of many-body processes which may become relevant and drive to phases with different conductance behavior, we can simply consider products of spin- and valley-conserving single-electron tunneling: (a) (b) FIG. 3. Schematic of the valley-preserving single-particle tunneling processes. Many-body processes which conserve spin and valley can be constructed as products of these processes. Among the processes which conserve valley, only spinconserving processes can become relevant and destabilize the fully pinched-off II (t) and fully-open (v) CC phases, due to the nature of the scaling dimension calculation. For each process about the charge and valley conducting phase (a), there is a dual process about the charge and valley insulating phase (b). The diagram here depicts modes for only a single spin direction; the full QPC hosts an additional set of modes related by a spin flip. αβ† β, L → R, α will be expressed as Oσu . Weak tunneling about the CC phase (Fig. 3a) is related to weak backscattering in the left-to-right direction, and is also dual to weak tunneling about the II phase (Fig. 3b). Taking advantage of this duality, we will restrict our discussion to the set of v operators which may destabilize the CC phase, noting that the t tunneling operators about the II phase are related by a duality. This duality is explored in greater detail in Ref. 28. The near-intersection of the two domain walls is a 0+1dimensional object, and therefore the coupling strength va of a given tunneling process flows, to first order, as dva = (1 − ∆(va ))va dl where ∆(va ) is the scaling dimension of the tunneling process Vα . To understand which operators may become relevant, we can first examine the single-electron tunneling processes (Fig. 3a): n † αβ Oσu = ψασRu ψβσLu , Vn−body = vn Oi + H.C. (22) i=1 where α, β = 0, π; σ =↑, ↓; u = K, K ; vn is the coupling strength of the process; and Oi is an arbitrary valley- and spin-conserving single-particle tunneling process. For the sake of condensing notation, tunneling from (23) αβ Oσu + H.C. V 1 = v1 (24) αβ,σ,u where α, β = 0, π; σ =↑, ↓; u = K, K ; and we have restricted ourselves to processes which preserve spin and valley. For single-particle tunneling, 7 B2- * * * ICB * * CIA * II/CC g- B1+ * * B2+ ICA * CIB * * * B1g+ FIG. 4. The regions in interaction space for which tunneling processes become relevant (∆ < 1 in Eqs. (28) and (30)) and the fully open (CC) or pinched-off (II) junction phases are destabilized. The central dot at g+ = g− = 1 is the noninteracting point and the dotted oval is the region of predicted accessible interaction strength in Ref. 25 for a suspended sample. The IC (CI) regions are characterized by four-body tunneling processes which transmit exclusively valley (charge) across the junction. The B regions represent relevant eight-body tunneling processes which are charge insulating (+) or conducting (−) and differ by band-index and valley-transmission character. Regions of overlap between boundaries, denoted with ∗, have multiple relevant operators at different orders and presumably more complicated behavior. In the central region, the fully open (CC) or pinched-off (II) phases remain stable and the conductance is characterized by single-electron tunneling. vant unless the operator in bosonized form only contains c-sector variables. For higher-body tunneling, this factor is greater than or equal to 1, and prevents any process which isn’t pairwise spin conserving and invariant under arbitrary SU(2) spin rotations from becoming relevant. One can view this as pairwise spin conservation allowing αβ decomposition into products of Oσu and SU(2) invariance providing the necessary linear combinations of σ =↑, ↓ to isolate c sector variables. Restricting to processes which conserve valley and obey these spin index constraints, the first operators to become relevant as interactions are increased are therefore a specific set of four- and eight-body tunneling processes. Figure 4 details the region in interaction space for which each class of many-body operators becomes relevant (∆(va ) < 1). When a process becomes relevant, the bosonized operators will become locked into the values which minimize the tunneling operator and become gapped out, altering the conductance of the junction. As an applied voltage only couples to the total electron density, only relevant operators containing θ+cρ or φ+cρ can cause the junction to become charge insulating. In the central region, all operators are marginal or irrelevant, though single-electron tunneling V1 is only the least irrelevant operator close to the non-interacting point g+ = g− = 1. Four groups of four-body processes are present which can become relevant and drive to phases which are completely charge or valley insulating: (1) 00 00 00 00 VICA = vICA O↑K O↑K O↓K O↓K + 0 ↔ π (2) 0π π0 0π π0 + vICA O↑K O↑K O↓K O↓K + 0 ↔ π + H.C. ∆(v1 ) = 1 1 1 g+ + + g− + +4 8 g+ g− (25) such that single-electron tunneling is always marginal or irrelevant (∆(v1 ) ≥ 1) for all possible inter- and intraband scattering processes. In the nearly noninteracting regime, where the least irrelevant operators are V1 , the strength of v1 can be controlled by an external parameter, such as the gate voltage VG for a given set of interaction strengths g± . In this interaction regime, we know that the CC and II phases are stable and that at least one quantum critical point exists to mediate the transition between them. In the subsequent section, Section III, we will use diagrammatic perturbation theory in the interactions about the CC phase to search for the set of possible intermediate phases and quantum critical points which characterize the single-electron-tunneling behavior of the junction. The factor of 4/8 in ∆(v1 ) is due to g±sρ = g±sv = 1 and will remain an obstacle to a process becoming rele- (1) 00 ππ 00 ππ VICB = vICB O↑K O↑K O↓K O↓K + 0 ↔ π (2) 0π 0π 0π 0π + vICB O↑K O↑K O↓K O↓K + 0 ↔ π + H.C. (1) 00† 00† 00 00 VCIA = vCIA O↑K O↑K O↓K O↓K + 0 ↔ π (2) π0† 0π π0† 0π + vCIA O↑K O↑K O↓K O↓K + 0 ↔ π + H.C. (1) ππ† 00 ππ† 00 VCIB = vCIB O↑K O↑K O↓K O↓K + 0 ↔ π (2) 0π† 0π 0π† 0π + vCIB O↑K O↑K O↓K O↓K + 0 ↔ π + H.C. (26) or in bosonized form: 8 (1) (1) (2) VICB = cos(θ+cρ ) (1) vICB cos(θ−cv ) + (2) vICB (1) (2) vCIB (1) ∆(vB ± ) 2 cos(φ−cρ ) VCIA = cos(θ+cv ) vCIA cos(θ−cv ) + vCIA cos(φ−cρ ) VCIB = cos(θ+cv ) 1 1 = (2) ∆(vB ± ) 2 = 8 . g± (30) (2) (1) vCIB (2) ∆(vB ± ) = ∆(vB ± ) = 8g± VICA = cos(θ+cρ ) vICA cos(θ−cρ ) + vICA cos(φ−cv ) cos(θ−cρ ) + cos(φ−cv ) (27) (1/2) where vIC/CI are the coupling constant strengths for the two different choices of interband scattering for each class of tunneling process and have absorbed factors of 2 during bosonization and simplification. The scaling dimensions for these couplings strengths are: (1) (2) (1) (2) (1) (2) (1) (2) ∆(vICA ) = ∆(vICA ) = 2g+ + 2g− ∆(vICB ) = ∆(vICB ) = 2g+ + ∆(vCIA ) = ∆(vCIA ) = ∆(vCIB ) = ∆(vCIB ) = 2 g− 2 2 + g+ g− 2 + 2g− . g+ (28) In the regions where the IC (CI) operators become relevant, the system will be charge (valley) insulating and valley (charge) conducting. Expressed in terms of conductance elements, the IC (CI) phases will have GXX = GY Y = 0, GV = GV Y = 0 (GXX = GY Y = 8e2 /h, XX Y V GV XX = GY Y = 0). These phases are related to the charge and spin insulating phases for the topological insulator QPC28 . However, unlike those phases, the IC and CI phases in the bilayer graphene junction will still be transmitting in the spin sector as long as g±sρ/v are relatively close to noninteracting. The remainder of the II/CC region is bounded by four eight-body processes. Each one can be considered as a product of one of the terms in the previous four-body processes multiplied by a product of selected conjugates of itself as to isolate just a single charge-sector variable. In bosonized variables, these eight-body tunneling operators are: (1) (2) They are only relevant under extremely strong interactions (g± < 1/8 or g± > 8) and are charge insulating (+) or conducting (−). Higher-order many-body tunneling processes are of course also possibly relevant, however due to the nature of the scaling dimension calculation, they will become relevant at much larger values of interaction strength than the boundaries of the shaded regions in Figure 4. III. THE PINCH-OFF TRANSITION As demonstrated in the previous section, under weak interactions the junction is stable in either the open (CC) phase or the closed-off (II) phase, both of which are characterized by single-electron tunneling and are related to each other by both 90◦ rotations and the pinch-off duality. In this section, we expand perturbatively in the interactions about the CC fixed point in search of the quantum critical point(s) which characterize the CC↔II quantum phase transition. In the process, we discover that, in addition to the T0/π = 1/2 critical point, which is expected as a consequence of the pinch-off duality, an additional family of intermediate critical points and phases are also present. For each of the possible paths between the II and CC phases we derive the conductance signatures which characterize the low-temperature transitions as functions of the two interaction strengths and the external gate voltage. First, we show how the general S-matrix characterizing the junction is renormalized by weak interactions, deriving a phase diagram and Renormalization Group (RG) in the case where scattering between the 0 and π bands is disallowed. We then allow for interband scattering, as might be present in the case of relatively smooth disorder, and introduce an S-matrix parameterization incorporating the additional system parameters. Using the results of an extensive renormalization group calculation, detailed in Appendix A, we assert that the most general S-matrix flows back to one with small-momentum conservation. The RG flow on this surface therefore contains all of the characteristic non-trivial quantum critical points for the pinch-off transition of this problem. Finally, we derive the critical exponents and universal scaling functions for the two classes of conductance transitions, up to leading order in the interactions. VB ± = vB ± cos(2θ±cρ ) + vB ± cos(2φ±cv ) 1 VB ± = 2 1 (1) vB ± 2 1 cos(2θ±cv ) + (2) vB ± 2 cos(2φ±cρ ) (29) The strengths of these operators have scaling dimensions: A. Non-Interacting Electrons In the absence of interactions, tunneling through the junction can be characterized by an S-matrix restricted 9 only by time-reversal symmetry and valley-index conservation αβ βσ ασ |ψi,out = Sij δ σσ |ψj,in t r −r† t† (32) where the rows and columns of SK indicate scattering of the incoming modes with valley index K (leads 1 and 3) to the outgoing ones (leads 2 and 4). The matrices r and t live in the 2 × 2 space of band indexes. Identical copies of SK exist for the up and down spins. The elements of SK are otherwise unconstrained if we allow scattering between the 0 and π bands such that SK is an arbitrary U (4) matrix. We can choose to parameterize SK = † U1 0 † 0 U3 t r −r† t† √ t= T0 √0 Tπ 0 , r= U2 0 0 U4 2 1 − T0 0 0 2 1 − Tπ Ui = e e e e . 𝜃1 ϕ21 𝜃2 ϕ41 ϕ23 𝜃4 ϕ43 𝜃3 3 4 FIG. 5. The eight angular variables in the full formulation of the S-matrix, seven of which are linearly independent. θi corresponds to the processes which break small-momentum conservation and scatter electrons between the 0 and π bands at lead i. φij = φi − φj corresponds to the scattering phase acquired tunneling from lead j to lead i at valley K. While there are four independent θi angles, there are only three independent φij as only the relative scattering phase matters, such that φ43 = φ23 −φ21 +φ41 . We have only displayed modes at K; an additional copy of this picture exists for modes at K , related by time-reversal symmetry. This picture is valid for electrons with either spin, as our model system is SU(2) spininvariant for all values of charge-sector interaction strengths. (33) (34) † where the Ui (Uj ) are for now unconstrained U(2) matrices which characterize operations on the outgoing (incoming) electronic wavefunction at lead i (j) and valley index K. In this parameterization, we can choose the tunneling probabilities for each band to be real such that T0/π = |t0/π |2 = 1 − |r0/π |2 . At this point, SK remains characterized by 16 free parameters. Choosing to parameterize the Ui in terms of Euler angles: iφi σ z iθi σ y iαi σ z iξi 2 (31) where i, j are lead indexes 1 − 4; σ, σ are spin indexes ↑, ↓; and α, β are band indexes 0 and π. Given the negligible spin-orbit coupling in graphene, we can consider here that the time-reversal operator τ = K leaves the spin sector unaffected, such that τ 2 = +1. Therefore, βα αβ time-reversal symmetry restricts that Sij = Sji . Additionally, to keep the domain wall states gapless, one must disallow scattering from K to K , which restricts T elements of the S-matrix SK = SK . In this section, we will work the most general allowed S-matrix down to one which is characterized by parameters which have physical meaning. Beginning with the modes in a single valley K: SK = 1 (35) Recognizing that αi and ξi are just U(1)×U(1) transformations at each lead, we can gauge them out. We are then left with nine variables which have physical significance. The four 0 ↔ π mixing angles θi correspond to the breaking of small-momentum-conservation at lead i, the three independent linear combinations of scattering phases φij = φi − φj each correspond to the phase acquired for electrons scattering from lead j to lead i at valley index K, and the two real tunneling probabilities T0 and Tπ characterize the extent to which the junction is pinched off for each band. These seven linearly independent angular variables (illustrated in Fig. 5) and two tunneling probabilities completely span the space of the gauge-independent, time-reversal-symmetric, spinand valley-index-conserving problem with S-matrix, T S = SK ⊕ SK = SK ⊕ SK (36) which has rows and columns characterized by lead indexes i, j. This parameterization of the S-matrix in lead and band-index spaces is restated more explicitly in the beginning of Appendix A 1. In the subsequent sections, we will show how the surface where θi , φij = 0 is not just a welcome simplification of the problem, but also the surface which contains all of the characteristic quantum critical points and their single relevant eigenvectors. 10 B. S-Matrix Renormalization under Weak Interactions As discussed in II A, this system is characterized by two kinds of density-density interactions. Here, we relate the S-matrix to the single-particle Green’s function and then use diagrammatic perturbation theory in those two interactions to find the leading order corrections to the S-matrix and derive RG flow equations for our system parameters28,30 . We calculate the relevant physics on the high-symmetry surface where the interband scattering angles θi = 0. In Appendix A, we further calculate the RG flow for all nine S-matrix parameters, demonstrating that the angles θi , φij either flow back to this surface or are marginal and trivial at the quantum critical points. 1. Gαβ,σσ ij  1  (z, z ) = 2πi δij δ αβ z−z αβ Sij z −z ¯ βα (Sji )∗ z−¯ z δij δ αβ z −¯ ¯ z   δ σσ (38) where z = τ + ix and the rows and columns of Gαβ ij are the in/out indexes such that elements proportional αβ to Sij correspond to in ↔ out. Restricting ourselves to the interactions introduced in II A (Eq. (2)), we can use diagrammatic perturbation theory to calculate the leading order corrections to Gαβ,σσ in the presence of ij weak u+/− , such that the Luttinger parameters g+/− = 1 + +/− and terms are kept up to O( 2 ). +/− The one-loop corrections to the S-matrix are qualitatively similar to those in Ref. 28, with special care taken to properly sum at each vertex over the tensor structure of the interaction Constructing the S-Matrix Renormalization Group αβγδ Scattering processes from one lead, band, and spin to another can be considered in terms of a single-electron thermal Green’s function = +δ αβ γδ δ αβ αβ Sij = Sij + 2 × Λ 1 log 4 E αβ A B (σA Sij σB ) αβ γδ − σz σz (39) where α, β, γ, δ are band indexes. Illustrated in Fig. 6, only two diagrams contribute to the renormalization of the single-particle Green’s function: abαβ,σσ bβσ aασ Gij (x, τ, x , τ ) = −i Tτ ψi (x, τ )ψj † (x , τ ) (37) where Tτ denotes imaginary time ordering and indexes a, b = in, out. In the absence of interactions + G out in,αβ,σσ = αβ 1 Sij δ σσ 2πi z − z ¯ (40) with † Tr[σA Sji σB Sji ] − A,B=0,1 † Sik σA Slk σB Slj αβ † Tr[σA Slk σB Slk ] kl (41) where Λ and E are the ultraviolet and infrared cutoffs respectively. Traces refer to band index space, A/B = +/− , and σA/B = 1, σz (σ0 , σ1 ); tracing over the spin degree of freedom is implicit, which leads to the prefactor αβ dSij 1 = dl 2 αβ A B (σA Sij σB ) of 2. The two terms correspond to the diagrams (a) and (b) respectively in Fig. 6. We can derive flow equations αβ for the elements of Sij by rescaling the cutoff Λ → Λe−l , † Tr[σA Sji σB Sji ] − A,B=0,1 One can immediately observe that this implies that θi = 0 is a fixed surface, since perturbative corrections to 0π 0π Sij = 0 require multiplying by Sij in the above equations. We can thus, for now, simplify our focus to the fixed surface where θi = 0. In Appendix A, we calculate the stability of the quantum critical points on this † Sik σA Slk σB Slj αβ † Tr[σA Slk σB Slk ] . (42) kl surface for a general set of θi and φij , and demonstrate that for physical values of g− relative to g+ , the critical points are stable in all possible out-of-plane, τ - and valley-symmetric directions. αβ Using our parameterization of Sij on the θi = 0 surface, we can exploit the matrix structure of the dia- 11 (a) (a) σ' u v l,in w α i,in σ (b) σ' k,out y γ k,in δ σ α β i,in j,out σ σ (c) (d) σ' w l,out σ v l,out z x u σ' k,in (b) z x y δ β k,out l,in j,out γ σ σ σ FIG. 6. The two non-zero, non-canceling diagrams for O( 2 ) perturbation theory. Note that i − l are spatial lead +/− indexes and α − δ, u − z are band indexes which are summed at each vertex over αβγδ . σ, σ are spin indexes for which the δ σσ within each Green’s function has already been taken into account. Even though the system is spin invariant, the spin index on the loop, here σ , must still be summed over its two values to calculate the physical RG flow. Each diagram contributes a logarithmic correction to the S-matrix which, when renormalized, leads to a term in Eq. (42). grammatic perturbations to obtain flow equations for the transmission probabilities for each band dT0 = −2T0 (1 − T0 ) ( dl +( + − 2 − ) (1 + − 2 − ) (1 + −) 2 (1 − 2T0 ) − 2Tπ ) dTπ = −2Tπ (1 − Tπ ) ( dl +( + + + 2 − ) (1 − 2Tπ ) − 2T0 ) . (43) This system of equations obeys two key symmetries. First and foremost, like the QSH point contact described in Ref. 28, it is invariant under the pinch-off duality T0/π ↔ (1−T0/π ), which we introduced in II B. Eqs. (43) are also invariant under exchange of the band indexes 0 and π. In the context of our renormalization group calculation, 0 and π are arbitrary labels for the band degree of freedom, and so even though this system contains stable fixed points with broken band-index symmetry, the system of flow equations itself must be band-index symmetric. Graphically, these two symmetries manifest themselves as two mirror symmetries in Fig. 7: one about T0 + Tπ = 1 and the other about T0 = Tπ . 2. Fixed Points and Renormalization Group Flow This system of flow equations can have as many as nine fixed points to quadratic order in the interactions (Fig. 7). The two corner fixed points at T0 = Tπ = 0, 1 are stable for all values of − . The central point at T0 = Tπ = 1/2 controls the transition between the CC FIG. 7. RG flow of the variables T0/π , which control the pinching off of the junction, calculated to quadratic order in the interactions +/− (Eq. (43)), panels (a)-(c). Large circles correspond to stable fixed points and small circles indicate nontrivial quantum critical points. The flow is controlled by the ratio of the interaction strengths − / + . When − / + = 1, the 0 and π bands are completely decoupled and each one behaves individually as a copy of the QSH problem in Ref. 28 (a). For − / + > 0, a set of intermediate fixed points exists which allows the 0 and π bands to be pinched off independently (b). When − = 0, the quadratic theory predicts that a fixed line will exist T0 +Tπ = 1 (c). Higher-order corrections about this line, calculated in Appendix B, infer flow along it back to the central quantum critical point T0 = Tπ = 1/2 (d). and II corners and is related to the T = 1/2 critical point in the single-band TI QPC case28 . As in that case, the existence of this central quantum critical point is mandated by the pinch-off duality T0/π ↔ 1−T0/π . The corner fixed points at T0 = 0, Tπ = 1 and T0 = 1, Tπ = 0 represent new, stable, single-electron-tunneling phases where only one of the bands is pinched off. We label these intermediate, mixed-band-character fixed points as “M Phases”. Transitions between the fully open or closed phases and these M phases are controlled by four fixed points which exist for − / + > 0 (Figs. 7a, 7b). For − = + , the 0 and π bands are completely decoupled and each one acts as a (spin-doubled) independent copy of the QSH problem in Ref. 28 (Fig. 7a). When − / + = 0, Eq. (43) predicts that all of the intermediate transitions and stable fixed points will collapse onto a fixed line at T0 + Tπ = 1 (Fig. 7c). This − = 0 case represents a restoration of U(2) band-index symmetry locally at each lead under which the band indexes become trivial. 12 We can examine this fixed line in greater detail by expanding upon our bosonization calculations from II B. In Appendix B, we increase the strength of π ↔ π singleelectron tunneling to drive from the CC phase towards an action about the M-phase corner fixed points for g− = 1. Calculating higher-order correlation functions about this theory and expanding perturbatively again in the interactions, we discover additional fixed points: porates the quadratic-order RG flow and the higher-order corrections near the T0 + Tπ = 1 line. As highlighted by the red and purple arrows in that figure, each of the two classes of nontrivial quantum critical points is characterized by only a single relevant direction in the T0 − Tπ plane. All other out-of-plane θi and φij directions are irrelevant or trivially marginal (Appendix A). 3. 48 1 − T0 = Tπ = 2 π 48 − 3 , T0 = 1 − Tπ = π 2 + − 3 + (44) where the second point is implied by the combination of mirror reflections about the pinch-off and band-indexexchange lines. Taking − → 0, this theory exhibits flow back towards the central quantum critical point T0 = Tπ = 1/2 (Fig. 7d). The simplest assumption would be to postulate that, to lowest order, this flow continues away from the vicinity of the M points without additional fixed points appearing. This implies that the fixed line at 2 − = 0 is simply an artifact of the O( +/− ) perturbation theory and that for extremely small − / + , the M phase is unstable and flow lines in that region point towards the central quantum critical point T0 = Tπ = 1/2. From this information, we obtain Figure 8, a schematic which incor- Conductance Signatures and Universal Scaling Functions For each class of quantum critical point, we can, knowing that its only relevant eigenvector lies on the θi = 0 plane, use the quadratic-order flow equations to calculate the universal conductance scaling and critical exponents to leading order in the interactions. Returning to a discussion of reduced conductance matrixes from II B, we can express the left-to-right, twoterminal conductance GXX in terms of the S-matrix. The elements of the four-terminal conductance G in the lead basis are related to the S-matrix by: Gij = 2e2 † Tr[1 − Sij Sij ] h (45) where i, j are lead indexes 1 − 4 such that the matrix  2 −T+ 0 −(2 − T+ ) 2e  −T+ 2 −(2 − T+ ) 0  G=   0 −(2 − T+ ) 2 −T+ h −(2 − T+ ) 0 −T+ 2  2 where T± = T0 ± Tπ and the factor of 2 on the conductance is due to electron spin degeneracy. The linear combinations of lead currents I1−4 which give IX,Y,Z are a result of, in combination with the requirement i Ii = 0,    I1 IX  IY  = M T  I2  I  3 IZ I4  (47) 1 1  −1 M=  2 −1 1 1 1 −1 −1  1 −1  1  −1 GXX = 4e2 8e2 (T0 + Tπ ) = − GY Y . h h GXX,A (∆VG , T ) = 8  (50) (48) such that GXY Z = M T GM = This confirms the result from II B that, for the valleyconserving problem, GZZ is quantized to be 8e2 /h regardless of the junction state. This reduction also confirms that, in terms of the S-matrix elements, With the structure of the conductance matrix established we can analyze the finite-temperature conductance ∗ transitions near each transition voltage VG . First, we will consider the direct II-CC phase transition, for which GXX scales from 0 to 8e2 /h. We can write the conductance in its scaling form where  (46)  0 0 4e2  T+ 0 2 − T+ 0  . h 0 0 2 (49) e2 ∆VG GA c αA h T (51) ∗ where ∆VG = VG − VG,A in Fig. 9, c is a non-universal constant, and the subscript A on G and α denotes the direct quantum phase transition between the II and CC phases. Observing that infinitesimal movement of T− 13 (a) M 0 𝜋 0 and 𝜋 CC M CC (a) (b) II II M M 0 and 𝜋 0 𝜋 FIG. 8. A schematic phase diagram, in terms of left-to-right conductance, within the T0 −Tπ plane, combining information from Eq. (43) and Appendix B. There are two classes of quantum critical points. The central point controls transitions between the fully open (CC) phase and the fully pinchedoff (II) phase. Four additional critical points on the edges control transitions between the CC/II phases and an intermediate mixed (M) phase in which the two bands have differing conductance contributions. The width of the M phase is O( − / + ). away from 0 is irrelevant (Fig. 8), we can set T− = 0 and characterize the conductance transition from II-CC with a single parameter, dT+ = −2( dl FIG. 9. A reproduction of Figure 8 with dashed lines overlaid to indicate possible voltage curves. As the gate voltage VG winds along a voltage curve, whose exact curvature is dictated by experimental specifics, it passes directly from the II to the CC region along a curve like (a) or indirectly, passing along the way through an intermediate M phase along a curve like (b). At zero temperature, the left-to-right conductance GXX will therefore undergo a direct transition from 0 → 8e2 /h along (a) or one with an intermediate step up to 4e2 /h along (b). This behavior motivates us to search for the finite-temperature scaling of these conductance transitions for ∗ ∗ ∗ VG ∼ VG,A (a), or for VG ∼ VG,B and VG ∼ VG,C (b). where T+ = 2GA such that X ∝ ∆VG /T 2( 2 + + 2 − )T+ (2 − T+ )(1 − T+ ). 0 T+ (2 0 − T+ ) −4( e 0 − T+ )2 2 2 + + − )l = T+ (2 − T+ ) (1 − T+ )2 (53) where we have purposefully left T+ (l) in its implicit form to provide a framework for the more mathematically complicated II-M transition analysis later in this section. 0 As one adjusts the gate voltage, T+ passes through 1 at ∗ 0 VG = VG,A , such that near the transition ∆VG ∝ T+ − 1. To determine the critical behavior, we can therefore ex0 pand T+ around this value. Finite temperature T can be taken into consideration by cutting off the renormalization group flow l at Λe−l ∝ T . Taking ∆VG , T → 0 at an arbitrary ratio, the previous equation can be rewritten 1 GA (1 − GA ) = (2X)2 (1 − 2GA )2 (54) or explicitly inverted GA (X) = 1 X 1+ √ 2 1 + X2 2 2 ++ −) . (52) This equation can be integrated to determine the 0 crossover scaling function. Taking T+ = T+ at l = 0, (1 (b) (55) αA = 2( 2 + + 2 −) (56) is the universal critical exponent for the II-CC quantum phase transition. Figure 10 shows GA at finite temperatures, noting that it collapses onto a step function at T /c1/α = 0 and that it has a crossover point pinned at GA = 1/2 for all values of interaction, making it identical to the T-R scaling function for weak interactions in the related Quantum Spin Hall problem28 . When − = 0 and the only transitions are directly from II-CC, αA = 4αQSH , where αQSH is the T-R critical exponent in Ref. 28. This factor of 4 accounts for the fact that even though here only a single linear combination of indexes is being pinched off, the diagrams which contribute to the flow equations still live in a matrix space which is four times as large as that of the comparable Quantum Spin Hall problem. Obtaining the critical exponent and scaling function for the II-M quantum phase transition is procedurally identical. As an example, we will choose the bottom right transition point in Fig. 9 for our calculation, though that point is restricted to be the same as the one characterizing the finite Tπ II-M transition by band-indexexchange symmetry. The quantum critical point is located at Tπ = 0, T0 = γ/2 where 14 a GA b 1 1 2 GB 1 Γ 2 0 0 4 2 0 2 VG 4 4 2 0 2 VG 4 FIG. 10. The two classes of universal scaling functions as functions of external gate voltage: GA describes the direct II-CC quantum phase transition and GB describes the II-M transition, plotted in panels (a) and (b) respectively. Here, we have plotted using + = 0.212, − = 0.071, such that γ = 1.25. The curves are plotted for increasing temperature, with the red, orange, green, and blue curves representing T /c1/α = 0, 10−5 , 10−2.5 , and 1 V 1/α respectively in equations (54) and (61). Note that the crossover value of GA is fixed to be 1/2, whereas the crossover value for GB is instead at γ/2, where γ varies from 1 to 2 continuously as a function of interaction strength. γ =1+ ( ( − −+ − 2 +) . 2 +) (57) First and foremost, we can note that when γ = 2, = 0 and there is no more available phase space (at quadratic order) for the II-M transition to exist. Near this case, the M phase will exist in a vanishing area of phase space and most transitions will be controlled by the central quantum critical point in Figure 8. However, in Appendix A we have only analytically calculated the stability of the II-M quantum critical point to linear order in − , whereas our analysis of the critical behavior of the II-M transition is up to O( 2 ). We believe it reason− able to assume that this stability extends up to quadratic order in the interactions such that these quantum critical points still describe the relevant physical transitions in this problem. For the II-M transition, the conductance jumps from GXX = 0 → 4e2 /h, with here the T0 axis being the only relevant direction. For this transition then, GXX = 4e2 h T0 . Expressing the conductance in its scaling form − GXX,B (∆VG , T ) = 4 e2 ∆VG GB c αB h T (58) where again ∆VG is the external gate voltage, c is a non-universal constant that may certainly differ from the c in the II-CC transition, and the subscript B denotes that transition between the II-M phases. Taking Tπ = 0, FIG. 11. The universal scaling function GB which characterizes the II-M phase transition, plotted for T /c1/α = 10−3 V 1/α and g+ = 1.212. The blue, green, orange, and red curves are plotted at g− = 1.019, 1.047, 1.071, and 1.212 respectively. The critical value of GB for which the conductance flow changes from the II phase to the M phase occurs at the intersection of each curve with the ∆VG = 0 line, and varies as a function of the ratio of the interaction strengths g− /g+ . When g− = g+ , GB takes on the same functional form as GA in Fig. 10, though with a different critical exponent αB = αA . At T = 0 K, all of these curves collapse onto the same step function; they are increasingly distinguishable as temperature is increased. the conductance transition from II-M is characterized by the flow of a single parameter, dT0 = −2γ( dl + 2 − ) T0 (1 + − T0 ) 1 − 2T0 γ . (59) We can determine the crossover scaling function by 0 integrating this equation. Taking T0 = T0 at l = 0, 1 1 2 2 0 0 T0 (1 − T0 ) γ −1 1− 0 2T0 2 γ e−2γ( 2 ++ −) l = 1 2 −1 γ T0 (1 − T0 ) γ −1 1− γ 2T0 γ 2 γ . 1 2 −1 γ (60) As before, we can cut off the renormalization group flow at finite temperature T ∝ Λe−l and note that the ∗ 0 gate voltage VG = VG,B when T0 passes through γ/2, 0 such that near the transition ∆VG ∝ T0 − γ/2. Again taking ∆VG , T → 0 at an arbitrary ratio, we can finally arrive at an implicit equation for GB (X), (2−γ) 1 2X γ 2 = GB (1 − GB )γ 1− 2GB γ 2 (61) 15 2 where T0 = GB such that X ∝ ∆VG /T ( + + − ) γ(2−γ) . We can note that this equation reduces to Eq. (54) for γ = 1. That case represents − → + , for which T0 and Tπ act as independent copies of the Quantum Spin Hall problem in Ref. 28, but in a higher-dimensional space. Therefore, αB = ( + + − )2 γ (2 − γ) = αA (2 − γ) (62) is the universal critical exponent for the II-M transition. Confirming the relationship to the QSH point contact problem, setting γ = 1 gives αB | − = + = 2αA | − =0 = 8αQSH . Figure 10 shows GB as a function of ∆VG at different temperatures, noting that at zero temperature it is also a step function, indistinguishable from GA , but that at finite temperature it is defined by a crossover value of γ/2 which in general differs from that of GA (Fig. 11). Taking advantage of the duality between the II and CC phases, we can relate the remaining conductance crossover function, one which characterizes the M-CC phase transition, to GB . We can write the conductance in its scaling form GXX,C = 4e2 4e2 + GC h h c ∆VG T αC (63) where c is yet another non-universal constant and C denotes that M-CC transition. By pinch-off symmetry, we know that GY Y for the CC→M transition has to be equivalent to GXX for the II→M transition, therefore for the M→CC transition, ∆VG 4e2 GB −c αB h T GY Y,C = (64) and utilizing Eq. (50), 8e2 4e2 − GB (−X) h h GXX,C = (65) such that we finally deduce GC (X) = 1 − GB (−X) (66) where X ∝ ∆VG /T αC and αC = αB for weak interactions. (67) IV. DISCUSSION In this paper, we have computed the conductance signatures of the four-terminal intersection of two bilayer graphene domain walls. These domain walls can be induced by the presence of a perpendicular electric field and a change in either electric field direction or interlayer stacking. When valley-index is conserved, the domain walls are Luttinger liquids described with two nontrivial interaction parameters g± . The junction is analogous to a point contact and can be analyzed naturally using the language of quantum point contacts. As with a Quantum Spin Hall point contact, the physics of the junction is best understood in terms of reduced, two-terminal conductances. When interactions are strongly attractive (g± < 1/4) or strongly repulsive (g± > 4), the left-toright conductance can be strictly dominated by nonzero charge and valley conductances respectively. For weaker interactions (g± ≈ 1), both left-to-right conductances are nonzero and there exist several stable phases characterized by single-electron tunneling. Transitions between these phases are governed at low temperatures by universal scaling functions and critical exponents, which differ from those in the single-band QSH case and are functions of the two Luttinger parameters. We now briefly discuss the task of experimentally measuring the physics in this paper. First and foremost, the existence of a single domain wall in bilayer graphene requires prohibiting scattering between valleys K and K . Valley-index-breaking perturbations are strongly relevant and will significantly change the physics in both isolated domain walls and for the junction structures at their intersections. To this end, short-range disorder must be kept smooth on the scale of the lattice. Recent experiments have shown promising results in this direction, with fabricated samples showing single domain wall conductances up to nearly the quantized clean limit of 4e2 /h31 . Beyond these results, one should then attempt to verify the Luttinger liquid physics at a single domain wall by measuring the tunneling conductance at several low temperatures and its collapse onto a universal scaling function with critical exponent αT . In junction structures, the conservation of valley index can be confirmed by measuring the quantization of the reduced conductance element GZZ . Creating a multi-terminal junction like the ones we describe poses several fabrication and analysis difficulties which must be overcome to measure the point-contact physics in this paper. Forming a four-terminal junction of electric-field-induced domain walls requires patterning leads on the top and bottom of each of the four bulk regions, as well as a gate on the junction to control the weak-interaction pinch-off transition. Forming a junction from layer-stacking domain walls requires patterning fewer leads, but as clearly demonstrated in the samples in Ref. 22, the three-fold symmetry of the underlying graphene lattice restricts intersections of these domain walls to be six-terminal structures. Conductance transi- 16 tions in these six-terminal structures can be calculated and analyzed using the framework established in this paper, though the task will be algebraically more intensive. Tuning the two Luttinger parameters g± can be accomplished through turning a combination of experimental knobs. The strength of the overall effective Coulomb interaction can be altered for the domain wall states by tuning their widths with the strength of the perpendicular electric field. These changes will be mainly reflected in g+ , as it contains long-range contributions from the overall Coulomb interaction. The other Luttinger parameter, however, can only be adjusted by tuning shortrange interactions. This can be accomplished by testing samples on a variety of substrates. Working in order of increasing dielectric strength, one can work through a suspended sample, a silicon dioxide substrate, or a boron nitride substrate to tune down g− . Several simplifications and assumptions in this paper may not be exactly present under experimental conditions. The assumption that all domain walls in the sample have the same values of the Luttinger parameters requires that the perpendicular electric field strength be globally uniform in magnitude in the bulk regions and change in direction similarly smoothly across electricfield-induced domain walls. It also requires that there are no strong local variations in the dielectric strength and coupling of the underlying substrate. The physics of junctions may be significantly altered if these conditions are not realized experimentally; we have not investigated point contacts at the intersections of two domain walls with differing interaction strengths. We also only calculated universal scaling functions to leading order in the interactions about the noninteracting point g+ = g− = 1. When interactions become stronger and single-electron tunneling is no longer the least irrelevant operator, the critical exponents αA/B in the left-to-right conductance transitions may change significantly in their dependences on the interaction strengths, as they do in Ref. 28. Finally, we assumed that the Fermi energy was exactly at the particle-hole-symmetric point such that vF,0 = vF,π = vF . In practice, it will be quite difficult to exactly tune the Fermi energy to this point, and so for most experimental realizations, 0 ↔ π exchange symmetry in the variables will be broken. For our analysis, the effects of this can be realized by replacing the equal band-index-exchange symmetry which we exploited with one which flips band index and rescales variables by vF,0 /vF,π . This will result in changes to the scaling dimension calculations in Section II B and a relaxation of mirror symmetry about T0 = Tπ in Figure 8. Appendix A: Stability of Quantum Critical Points on the θi = 0 Surface In this appendix, we confirm analytically the stability of the two classes of quantum critical points in Fig. 8 which control the transitions between single-electrontunneling phases. Here, we begin with the general flow equation for elements of the S-matrix (Eq. (42)), leaving free the θi and φij variables in the full S-matrix formulation from III A. As the complexity of this calculation greatly grows with each power of − taken into consideration, we will only here carry out our stability analysis to quadratic order in + and linear order in − . Our methodology can be used to analytically calculate higherorder terms, but we believe it reasonable to truncate the calculation at this order and that the stability should carry over to the O( 2 ) calculation used to produce the − phase diagrams and scaling functions in III B. Additionally, numerical estimates find g− < g+ for a large range of system parameters25 , therefore the limitation to terms of order O( + − , 2 ) is physically motivated. + We can exploit the matrix structure of the diagrammatics by taking several traces of S-matrix products. Derivatives of these traces exploit the matrix structure of the flow equation, effectively closing the external legs of the diagrams in Fig. 6 into additional loops. We will follow through this calculation completely for a single trace † (Tr[Sij σz Sij ]) and then show how taking linear combinations of these traces results in a flow equation for θ1 . All remaining dθi /dl can be obtained by exploiting cyclic reindexing and pinch-off symmetries. 1. Flow Equations for θi Band Mixing Angles and φij Scattering Phases To begin this process, let us rewrite S-matrix elements in the 2 × 2 band space using a simplified version of Eqs. (33), (35): αβ Sij = (Ui† Dij Uj )αβ θi U i = ei 2 σy Dij = aij 1 + bij σz 1 0 aij = d + dπ ij 2 ij 1 0 bij = d − dπ ij 2 ij a  a 0 ta eiφA 0 ra eiφB a  ta eiφa A 0  0 −ra eiφC a  da =  a iφa  0 −r e C 0 ta eiφD  a a −ra eiφC 0 ta eiφD 0 (A1)  ACKNOWLEDGMENTS This work was supported by a Simons Investigator award from the Simons Foundation to Charles Kane. Fan Zhang was supported by UT Dallas research enhancement funds. (A2) where da is a collection of scalars in the band index ij space with a = 0, π, φ0 = −φπ = φij and all matrix ij ij 17 operations only involve the α, β = 0, π indexes. Transmission probabilities are normalized for each band such that (ta )2 + (ra )2 = 1 and the labeling A − C on the phases corresponds to: φA φB φC φD = φ21 = φ2 − φ1 = φ41 = φ23 = φ43 = φB − φA + φC . † † Tr[Sl2 Slk S1k S12 ] = Tr[(a∗ 1 + b∗ σz )(alk 1 + blk σz ) l2 l2 (a∗ 1 + b∗ σz )(a12 1 + b12 σz )]. (A7) 1k 1k (A3) We’ll begin by calculating three traces: † Tr[Sij Sij ] = |d0 |2 + |dπ |2 ij ij † Tr[Sij σz Sij ] = |d0 |2 − |dπ |2 cos θi ij ij † Tr[Sij σx Sij ] = |d0 |2 − |dπ |2 sin θi . ij ij commutivity issues. Consider first the simplest example, † † Tr[Sl2 Slk S1k S12 ]. Without any additional Pauli matrices, the Ui rotations all cancel out pairwise, resulting in: (A4) We can then differentiate and take weighted linear combinations of these traces to obtain explicit flow equations for θi . Specializing to i, j = 1, 2, One might be concerned that converting this to a useful 0/π form, one with dij where S-matrix elements can just be read off, would be a daunting and terrible task. However, 0/π converting to dij basis is equivalent to summing over all of the ways to choose + and − signs for the cross terms, and so for every single one of these four S-matrix traces, all of the cross terms cancel. Our trace is reduced to the quite simple form: † † Tr[Sl2 Slk S1k S12 ] = (da )∗ da (da )∗ da . l2 lk 1k 12 (A8) a=0,π (T0 − Tπ ) d dθ1 † = (cos θ1 ) Tr[S12 σx S12 ] dl dl d † − (sin θ1 ) Tr[S12 σz S12 ] dl (A5) where T0/π = |t0/π |2 . Now, we can examine, in detail, the process of using the diagrammatics (Eq. (42)) to cal† d culate one of these derivatives, namely dl Tr[S12 σz S12 ]. All other trace derivatives, while varying in signs and specifics, follow procedurally from this example. First, by using the cyclic index definition of the trace, we can see that the derivative of the trace is equal to the trace of the chain rule derivative: Worth noting is that this picture is significantly complicated by the addition of Pauli matrices between Smatrix factors, due to the fact that [Ui , σz/x ] = 0. While the complexity doesn’t significantly increase for the addition of a single Pauli matrix, as it can be absorbed into the definition of Dij , it does for two or more Pauli matrices. This can be seen at the level of the two-S-matrix trace, where commutivity issues lead to the addition of a second term: † Tr[σz Sij σz Sij ] = cos θi cos θj |d0 |2 + |dπ |2 ij ij d dS † dS12 † † Tr[S12 σz S12 ] = Tr[ 12 σz S12 ] + Tr[S12 σz ] dl dl dl dS12 † = Tr[S12 σz ] + C.C. (A6) dl where for the second equality we exploited the cyclic nature of the trace to reduce this step to the calculation of one, albeit large, trace. In this case, the trace of the derivative is real, but in the next section where we calculate dφij /dl, it will not be and the addition of the complex conjugate cannot be overlooked. From here, the calculation amounts to taking the traces of terms which contain products of two or four S-matrices. While the products of two, in the form of Eq. (A4), can be calculated by rote algebra without much difficulty, the terms with four S-matrices require a bit more manipulation. We calculate the traces of products of four S-matrices by both recognizing a pattern in the assignment of signs to the products of aij , bij and with a careful treatment of + sin θi sin θj d0 (dπ )∗ + (d0 )∗ dπ . ij ij ij ij (A9) For the four-S-matrix traces, the weighting of the cross terms is altered by the anticommutivity of the Pauli matrices and while only two terms remain for each trigonometric function of θi , they are different than the simple form of Eq. (A8) and contain possible terms which mix 0 and π. As the number of Pauli matrices inserted into these traces increases linearly with the power of − to which we expand, we have chosen for the sake of simplicity and clarity to expand only to linear order in − . Though our analysis in section III B 2 continues to O( 2 ), we be− lieve that the stability deduced here should carry over to higher order terms. With our calculation machinery established, we can produce a flow equation for the 0 ↔ π mixing at lead 1: 18 dθ1 = dl − + − sin 2θ1 [T0 (1 − Tπ ) + Tπ (1 − T0 )] + sin 2θ2 cos (2φ21 ) T0 Tπ [2 − T0 − Tπ ] T0 (1 − T0 )Tπ (1 − Tπ ) + 2 sin 2θ3 cos (2φ31 ) + sin 2θ4 cos (2φ41 ) (1 − T0 )(1 − Tπ )[T0 + Tπ ] . (A10) Flow equations for the remaining three mixing angles can be generated by exploiting underlying symmetries of the problem. The set of nine independent variables which characterizes the S-matrix obeys three symmetries. Two “pinch-off” symmetries exist; the duality between the fully closed (II) and fully open (CC) single-electron phases implies that the system of flow equations is invariant with respect to the exchange T0/π ↔ (1 − T0/π ) and either the exchange of lead indexes 2 ↔ 4 or 1 ↔ 3. The third symmetry is a cyclic relabeling of the lead indexes as well as an exchange of the definitions of pinched off and open, due to the system’s invariance under properlytreated (with respect to valley) 90◦ rotations. Therefore the system is also invariant under the exchange T0/π ↔ (1 − T0/π ) and 1 → 2, 2 → 3,3 → 4, and 4 → 1. Independent calculations of dθi /dl confirm these properties. We can see immediately that for − = 0, all θi are marginal. In this case, the system has U(2) symmetry at each lead and all band index rotations can be gauged out. The dependence of dθi /dl on θj=i can also be suppressed dφA 1 = dl 2 + − tan θ1 − θ2 2 by tuning the scattering phase φij closer to π/4, which amounts to having a π/2 scattering phase difference between the 0 and π bands. Calculating the flow equations for scattering phases φij requires, conversely, tracing over open paths in diagrams, which allows phase to accumulate throughout the summation instead of being canceled out pairwise as frequently occurred in the calculation of dθi /dl. Specializing for the moment towards obtaining dφA /dl, we can note the following: Tr[S12 ] = (t0 eiφA + tπ e−iφA ) cos | Tr[S12 ]|2 = cos2 θ1 − θ2 2 θ1 − θ2 2 T0 + Tπ + 2 T0 Tπ cos 2φA (A11) where the correspondence between the subscripts A−D and the lead indexes i, j comes from Eq. (A3). This can be differentiated, and, with a considerable amount of algebra, used to obtain first order flow equations for the scattering phases. While the specifics of this calculation differ from those of obtaining the dθi /dl, the key point about the summation over + and − possibilities when 0/π converting to the dij basis remains for both the oneand three-S-matrix products here, again greatly simplifying the algebra for the required trace calculations. Utilizing this fact, we obtained an explicit flow equation for scattering phase φA : (sin 2θ1 − sin 2θ2 ) sin (2φA ) T0 Tπ (2 − T0 − Tπ ) + sin 2θ3 (1 − T0 )(1 − Tπ ) (T0 + Tπ ) sin (2φC ) − 2 T0 Tπ sin (2(φA − φC )) − sin 2θ4 (1 − T0 )(1 − T π) (T0 + Tπ ) sin (2φB ) − 2 The flow equations for φB−D can be obtained by exploiting pinch-off and cyclic symmetries as well as the redundancy of φD (Eq. (A3)). As with the θi , one can observe that for the U (2) symmetric case of − = 0, all φij are also marginal and can be gauged away. One can dφij also observe that dl = 0 when θi = θj , as in that case there is no relative interband scattering between leads i and j and φij can be gauged out. Only two unique stability calculations are required, as the four critical points on the boundary of the square are related by pinch-off and band-index-exchange symmetries. 2. T0 Tπ sin (2(φA − φB )) . (A12) Quantum Critical Point Stability We can now utilize the flow equations for θi and, to a lesser extent, φij to determine the stability of the θi = φij = 0 surface quantum critical points in Fig. 8. We will begin by considering the T0 = Tπ = 1/2 central quantum critical point, which controls the direct transition between the fully pinched-off (II) and fully-open (CC) phases. Expanding the θi to linear order:  θ1 d  θ2   =− dl θ3 θ4   θ1  θ2  +M  θ3  θ4  − (A13) 19 Mij = cos 2φij (A14) which gives two Lyapunov exponents of zero and two of λ± = −2 + − 1 ± 1 i,j cos 4φij . Both λ± eigen4 values are either zero or negative for all choices of the independent φA−C . Given the high symmetry of this central quantum critical point, it is unlikely that the λ = 0 marginal directions correspond to instabilities at higher orders. The four external critical points which control the transitions to and from the mixed (M) phase can be handled similarly, but with an attention to the expansions used to arrive at this analysis. As our model is invariant under global exchange 0 ↔ π and the pinch-off symmetries, we perform stability analysis on just one of the four critical points and relate the rest by symmetry. Choosing the critical point at T0 = 0, Tπ = γ/2, we can again construct a stability matrix:  θ1 d  θ2   =− dl θ3 θ4   θ1  θ2  − γN  θ3  θ4  + is a genuine physical phenomenon or an artifact of our O( 2 ) perturbation theory. The first section of this ap± pendix details the derivation of a Euclidean action about an M fixed point, working from the action for the fullyopen CC phase. The second section utilizes that action to calculate flow equations beyond quadratic order in the vicinity of the M phase, resolving the fixed line to additional fixed points for weak-to-moderate interactions. 1. Deriving the Euclidean Action for the T0 = 1, Tπ = 0 M-Phase Fixed Point To examine the higher-order behavior of tunneling processes about the T0 = 1, Tπ = 0 corner fixed point in Eq. (43), we have to first obtain a theory for that fixed point in terms of our existing theory for the T0,π = 1, fully-open fixed point. As we will be focusing on the fate of the quadratic-order fixed line at g− = 1, we will for the purposes of this calculation specialize to g+ = g, g− = 1. Integrating out the φ+/−,c/s,ρ/v sectors in Eq. (17), we can propose as a starting point the Euclidean action about the CC fixed point: (A15) SCC =  0 0 0 cos 2φB γ 0 0 cos 2φC 0  N = 1+ 1 −  . 0 cos 2φC 0 2 cos 2φB 0 0 (A16) Acknowledging that, as we have expanded to linear order in − , γ ≈ 2 − 4 − / +, all of the off-diagonal terms become infinitesimal and we are left with:   θ1 d  θ2    = −2 dl θ3 θ4  θ1  θ2  θ  3 θ4  + − 1− 2 − + (A17) which clearly has strictly negative Lyapunov exponents up to this order of expansion. Appendix B: Resolving the Fixed Line in the U(2) Symmetric Case The flow diagram in Fig. 7, obtained by perturbation theory to quadratic order in the interactions, suitably describes the transport physics of single-electron-tunneling phases for most ranges of infinitesimal interactions + and − . However, when − = 0, Eqs. (43) become symmetric under arbitrary U(2) transformations in band index space and describe a fixed line T0 + Tπ = 1. In this appendix, we use bosonization at the T0 = 1, Tπ = 0 “Mphase” fixed point to determine whether this fixed line vT v |ωn |θα gα θα (B1) α=c,s ωn where  θ+αρ θ  v θα =  +αv  , θ−αρ θ−αv    g 0 0 0 1 1 0  g 0 0 , g =  0  gc =  0 0 0 1 0 s 0 0 0 0 1   1 4πβ (B2) 0 1 0 0 0 0 1 0  0 0 . 0 1 (B3) We can introduce weak tunneling processes, v0/π † ψ0/π,↑RK ψ0/π,↑LK + K ↔ K + ↑↔↓ + H.C. 2 = v0/π cos θ0/π,σu (B4) V0/π = σ,u where σ =↑, ↓ and u = K, K . These processes will generally be present within the vicinity of the T0 = Tπ = 1 fixed point. Both of these processes are just single-electron tunneling, so as described in II B, they are marginal or irrelevant (∆(v0 ) = ∆(vπ ) = 1 8 [g + 1/g + 6]) for all values of g, such that we may arbitrarily increase the coupling strength v0/π without new, additional tunneling processes becoming relevant. Respecting time-reversal, valley-index, and band-index symmetries, V0/π are restricted to live in the plane of 20 Fig. 8 and therefore turning up v0/π represents motion away from the CC fixed point along the T0/π axes respectively. For this reason, we need to return SCC to the 0, π basis so that we can generate an M phase action by taking the tunneling strength vπ → ∞ and “pinch off” just the π band. Additionally, we will find it easier to work in a basis for which the tunneling operators V0/π are just sums of cosines. From our bosonization work in II and utilized earlier in this appendix, we recall that θ±σu = θ0σu ±θπσu . We will also need to transform valley indexes utilizing the property that θaσ,K/K = θaσρ ± θaσv where a = 0, π. Combining these definitions, we create the unitary change-of-basis    θ0αK θ+αρ θ  θ  v θα =  +αv  = Y  0αK  = Yθα θπαK θ−αρ θπαK θ−αv  1 11 Y=  2 1 1 I I I −I 1 −1 1 −1 1 1 −1 −1 1 4πβ  1 −1  −1  1 (B6) T |ωn |θα Ygα Yθα (B7) where Ygc Y separates into blocks: 1 4 A B B A (B8) where 1 1 g+ g +2 g− g 1 1 g− g g+ g +2 = (g + 1/g + 2) I + (g − 1/g) σ x 1 1 g+ g −2 g− g 1 1 g− g g+ g −2 = (g + 1/g − 2) I + (g − 1/g) σ x = A − 4I. (B9) As the spin sector is noninteracting, Ygs Y = I. To complete the transformation into the basis of VCC = V0 + Vπ , . (B11) 1 4πβ Ygc Y + I Ygc Y − I Ygc Y − I Ygc Y + I θ↑ θ↓ ωn (B12) The most general action, to first order in irrelevant electron tunneling processes, is then: SCC = |ωn | θ↑ θ↓ β α=c,s ωn Ygc Y = θ↑ θ↓ Applying this to our action, we can complete this change of basis: S = SCC + and we can note that YT = Y−1 = Y. Under this transformation, the action becomes B= = 0  A= θc θs (B5) where α = c, s and SCC = we now need to address the spin degree of freedom. Utilizing the definitions of the spin and charge sectors (Eqs. (6)), (B10) dτ VCC . τc (B13) To drive towards the T0 = 1, Tπ = 0 fixed point, we need only take vπ → ∞, as we have already established v0 and vπ as the quantities which push back along the T0 and Tπ axes respectively. To massage something useful out of this strongtunneling limit, we can utilize an extreme limit of the Villain approximation, following procedurally a longer formulation of the Kane-Fisher problem32 . Note that the partition function is the product of the contribution from SCC and a term of the form e−vπ dτ cos θ . As vπ → ∞, the entire partition function is zeroed out except when cos θ is a minimum. Acknowledging this, one can make the substitution of e−vπ dτ cos θ → m eimθ where m is now a discrete step in time. The partition function can then be integrated over θ once we complete the square to give a new Gaussian effective action in terms of the integer m. Now, we can define m = ∂t ϕ/2π, such that in the frequency domain, our bare action has the same form and an inverse Luttinger parameter. Hopping events between the minima of cos θ are instantons whose term in the action takes the form t dτ cos ϕ, the same as if we had defined our electron operators as exponentials of ϕ fields and examined tunneling about the bare action. This equivalence between strong electron tunneling in θ and weak quasiparticle backscattering in ϕ is a key feature of the Kane-Fisher problem for single-impurity scattering in Luttinger liquids. Its physical implications, as well as a much more detailed derivation of it, can be found in Ref. 32. For four-terminal quantum point contacts, however, this equivalence manifests itself as a relationship between strong left-to-right electron tunneling and weak top-tobottom electron backscattering28 . Therefore, our ultimate goal is to arrive at ϕπ operators which correspond to the tπ→π single-particle tunneling processes in Fig. 3 about an M point in Fig. 8, for which only the π band is pinched off. . 21 In our problem, the result of taking vπ → ∞ and utilizing this trick is a fair bit more complicated due to the larger set of operators in the action and the matrix nature of the interactions. First, let us start by reducing this problem into manageable blocks. We can simplify using the following definitions: Ygc Y ± I = U± = 1 4 A ± 4I B B A ± 4I = 1 4 SCC = ωn β 0 dτ 2i VCC = τc β  T T θσ U+ θσ + 2θ↑ U− θ↓  (B15) σ=↑,↓ β mT θπσ + πσ 0 σ=↑,↓ θaσ = θaσK θaσK mπσ = (B16) As addressed in Eq. (B13), the full action contains contributions from SCC as well as perturbative tunneling terms. To push towards the corner M phase, we can increase vπ greatly until it becomes valid to replace that S= ωn |ωn | 16πβ → ϕT θπσ πσ ϕπσK ϕπσK 1 2π = ∂t ϕπσK ∂t ϕπσK iω ϕπσ . 2π (B19) a=0,π + 0 σ=↑,↓ iω 2π = T T T θa↑ A− θa↓ + 2θ0↑ Bθπ↓ + 2θ0↓ Bθπ↑ +2 a=0,π β + 16i sgn(ωn ) mπσK mπσK Combining all of these definitions, we arrive at a full expression for the action which has both θ and ϕ fields in it: T T θaσ A+ θaσ + 2θπσ Bθ0σ σ=↑,↓ (B18) such that θσ = θ0σ ⊕ θπσ . Tπ is an instanton tunneling term which enforces the condition that maσK/K is an integer. The factor of 2 in front of it is a normalization requirement from the c, s →↑, ↓ change of basis. As described earlier, we can recognize the maσK/K as discrete steps in time of a new field: where U± = Ygc Y ± I. dτ [V0 + Tπ ] τc (B17) A± B B A± (B14) such that now  |ωn |  4πβ part of the partition function with discrete delta functions in θπσ : dτ [V0 + Tπ ] τc (B20) where we have frequently exploited that A± and B are symmetric. After a good bit of algebra, we can complete the square twice here and integrate out the θπσ , leaving us with an intermediate action SM about the T0 = 1, Tπ = 0 Mphase stable fixed point:  SM = 1 4πβ |ωn | θ0↑ θ0↓ ϕπ↑ ϕπ↓ ωn VM =  θ0↑   gM  θ0↓   ϕπ↑  ϕπ↑ (B22) v0 cos (θ0σu ) + tπ cos (ϕπσu ) (B23) σ,u β S = SM + 0 dτ VM τc (B21) where σ =↑, ↓; u = K, K ; and  2I + ασ x ασ x − sgn(ωn )ασ x − sgn(ωn )ασ x x x x x g−1 ασ 2I + ασ − sgn(ωn )ασ − sgn(ωn )ασ   = , α = sgn(ωn )ασ x sgn(ωn )ασ x 2I − ασ x −ασ x g+1 sgn(ωn )ασ x sgn(ωn )ασ x −ασ x 2I − ασ x  gM (B24) 22 in which the sign on the bottom right block of gM is owed to one of the sets of bosonic fields being evaluated at −ωn and the other at +ωn . The tunneling term tπ for the ϕ fields represents weak, left-to-right tunneling of electrons with band index π; making it very large would return the system back to the CC corner fixed point. Correlation functions about this theory can be calculated using elements of the inverse of this matrix:  2I − ασ x −ασ x sgn(ωn )ασ x sgn(ωn )ασ x 1 −ασ x 2I − ασ x sgn(ωn )ασ x sgn(ωn )ασ x  =   − sgn(ωn )ασ x − sgn(ωn )ασ x 2I + ασ x ασ x 2 x x x x − sgn(ωn )ασ − sgn(ωn )ασ ασ 2I + ασ  2g−1 M where the factor of 2 again comes from the spin-index change of basis and ensures that diagonal correlators like eiθ0↑K (τ ) eiθ0↑K (0) ∼ 1/τ 2 are appropriately marginal. The calculations in the following section are also only appropriate provided that the single-particle, off-diagonal correlation functions are zero. This requirement in enforced for the θ fields for 2 − α > 1 and by 2 + α > 1 for the ϕ variables. However, the range −1 < α < 1 is satisfied by all physical values of g. Therefore, the only bounds on the validity of the SM theory are the regions in Fig. 4 for which many-body processes become relevant. 2. (B25) Cubic-Order Fixed Points about the M Phase for g− ≈ 1 In this section, we will, using the SM theory established in the previous section, develop quartic-order flow equations for v0 and tπ for g− = 1, g+ = g. After establishing the flow in those variables, we will relate v0 and tπ to T0/π in Eq. (43). Expanding g = 1 + + , we will develop the cubic-order corrections to Fig. 7 close to the M-phase fixed points and obtain the schematic phase diagram Fig. 8. There are two sets of correlation functions for each coupling coefficient that we must calculate and rescale to obtain flow equations, as indicated by the two non-zero, off-diagonal elements for each operator and spin in gM . For now, we will work just with the v0 renormalization and then use band-index-exchange and pinch-off symmetries to relate them to the flow of tπ . Starting with the θ0↑K operator, we can calculate terms in the cumulant expansion to cubic order: ∞ 3 δv0 = v0 dτ1 dτ2 Tτ eiθ0↑K (τ ) eiθ0σK σ=↑,↓ (τ1 ) −iθ0σK (τ2 ) e > − eiθ0↑K (τ ) > Tτ eiθ0σK (τ1 ) −iθ0σK (τ2 ) e > −∞ ∞ + v0 t 2 π dτ1 dτ2 Tτ eiθ0↑K (τ ) eiϕπσK σ=↑,↓ (τ1 ) −iϕπσK (τ2 ) e > − eiθ0↑K (τ ) > Tτ eiϕπσK (τ1 ) −iϕπσK (τ2 ) e >. −∞ (B26) As before, Tτ is the time-ordered product and for each of the integrals, the second term is the disconnected ∞ 3 δv0 = 2v0 + dτ1 dτ2 −∞ ∞ 2v0 t2 π −∞ 1 2 τ12 dτ1 dτ2 τ1 τ2 piece. Taking into account all possible time orderings, we can rewrite this: α −1 1 πiα Θ(τ1 < τ < τ2 ) + e−πiα Θ(τ2 < τ < τ1 ) − 1 2 Θ(τ ∈ (τ1 , τ2 )) + e τ12 where τ12 = τ1 − τ2 and Θ(τ, τ1 , τ2 ) is the Heaviside (B27) step function expressed in conditional notation. These 23 integrals can be performed analytically and, after a bit of work, we arrive at an explicit equation for δv0 = v0 − v0 : 3 δv0 = log b 4v0 πα tan vicinity of the possible fixed line in Fig. 7b. Expanding g = 1 + + , α2 = 2 /4 such that Eq. (B29) reduces to + the − = 0 limit of Eq. (43) under the substitution: πα πα − 8v0 t2 sin2 π 2 2 v0 = (B28) where b is the time-integral cutoff. Rescaling b → be−l and taking into account the flow equations’ invariance under combined exchange v0 ↔ tπ , α ↔ −α, we arrive at our higher-order flow equations: πα dv0 πα 3 − 8v0 t2 sin2 = 4v0 πα tan π dl 2 2 πα dtπ πα 2 − 8v0 tπ sin2 = 4t3 πα tan π dl 2 2 + − (1 dTπ = −8 dl + − Tπ 2 + + π2 4 + 48 Equation (B 2) contains a few points of interest. First and foremost, it reflects the pinch-off symmetry in its invariance under the exchange (1 − T0 ) ↔ Tπ . It also − contains fixed points at T0 = 1− 2 + , Tπ = 0 and at T0 = − 0, Tπ = 2 + , in agreement with the small − / + limit of Fig. 7. We specifically normalized our correlations about SM to fix this agreement, such that we could locate at quartic order in + any local quantum critical points which control the transition between the CC/II and M phases. On the line 1−T0 = Tπ , Eq. (B 2) admits an additional fixed point at 1 A. H. Castro Neto, F. Guinea, N. M. R. Peres, K. S. Novoselov, and A. K. Geim, Rev. Mod. Phys. 81, 2 π Tπ . (B30) when 1−T0 1, Tπ 1. We also know that the linear terms in dT0/π /dl must agree with any local expansion in the T0 − Tπ plane, such that we can restore flow to linear order in − near the M-phase by extracting the leading term in Eq. (43): O( + −) dTπ dl (B29) − T0 ) + 4(1 − T0 )2 2 + 4Tπ 1 − T0 , tπ = dT0 dl which are valid only in a perturbative vicinity of this M-phase fixed point and for weak-to-moderate interactions. Now, to answer our initial question, we will examine the consequences of our new flow equations near the d (1 − T0 ) = −8 dl 2 π O( + −) =8 + − (1 = −8 − T0 ) + − Tπ (B31) where we have used the same approximation 1 − T0 1, Tπ 1. Combining this with an expansion of Eq. (B29) to quartic order in + , we obtain, finally, higher-order flow equations about an M-phase fixed point: π2 4 π2 4 + + − 4(1 − T0 )Tπ 2 − + 48 48 π2 4 + − 4(1 − T0 )Tπ 2 − . + 48 2 + + 1 − T0 = Tπ = 48 π2 − 3 . + (B32) (B33) . Relevant flow lines from this point head off towards the central T0 = Tπ = 1/2 critical point, even in the − = 0 case when this point converges with the other critical points at the corner. Therefore, this point stands as a demonstration that the fixed line in Fig. 7b is an artifact of O( 2 ) perturbation theory, at least in the vicinity of + the T0 = 1, Tπ = 0 M point. It is most likely that this flow away from the M point continues, to lowest order, all the way to the central quantum critical point. This infers that the phase diagram of this system for weak-tomoderate interactions is best described by the schematic Fig. 8. 109 (2009). 24 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 E. McCann and V. I. Fal’ko, Phys. Rev. Lett. 96, 086805 (2006). F. Zhang, J. Jung, G. A. Fiete, Q. Niu, and A. H. MacDonald, Phys. Rev. Lett. 106, 156801 (2011). J. Velasco Jr., L. Jing, W. Bao, Y. Lee, P. Kratz, V. Aji, M. Bockrath, C. N. Lau, C. Varma, R. Stillwell, D. Smirnov, F. Zhang, J. Jung, and A. H. MacDonald, Nature Nanotechnology 7, 156 (2012). The related discussions were provided in the supplementary materials. C. L. Kane and E. J. Mele, Phys. Rev. Lett. 95, 226801 (2005). F. Zhang, H. Min, M. Polini, and A. H. MacDonald, Phys. Rev. B 81, 041402(R) (2010); F. Zhang, H. Min, and A. H. MacDonald, Phys. Rev. B 86, 155128 (2012). R. Nandkishore and L. Levitov, Phys. Rev. Lett. 104, 156803 (2010). T. Ohta, A. Bostwick, T. Seyller, K. Horn, and E. Rotenberg, Science 313, 951 (2006). J. B. Oostinga, H. B. Heersche, X. Liu, A. F. Morpurgo, and L. M. K. Vandersypen, Nature Materials 7, 151 (2007). H. Min, B. Sahu, S. K. Banerjee, and A. H. MacDonald, Phys. Rev. B 75, 155115 (2007). W. Bao, J. Velasco Jr, F. Zhang, L. Jing, B. Standley, D. Smirnov, M. Bockrath, A. H. MacDonald, and C. N. Lau, Proc. Natl. Acad. Sci. USA 109, 10802 (2012). H. Min, J. E. Hill, N. A. Sinitsyn, B. R. Sahu, L. Kleinman, and A. H. MacDonald, Phys. Rev. B 74, 165310 (2006). Y. Yao, F. Ye, X. Qi, S. Zhang, and Z. Fang, Phys. Rev. B 75, 041401 (2007). Y. Zhang, T. Tang, C. Girit, Z. Hao, M. C. Martin, A. Zettl, M. F. Crommie, Y. Shen, and F. Wang, Nature 459, 820 (2009). I. Martin, Y. M. Blanter, and A. F. Morpurgo, Phys. Rev. Lett. 100, 036804 (2008). F. Zhang, A. H. MacDonald, and E. J. Mele, Proc. Natl. Acad. Sci. USA 110, 10546 (2013). 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 W. Yao, S. A. Yang, and Q. Niu, Phys. Rev. Lett. 102, 096801 (2009). J. Jung, F. Zhang, Z. Qiao, and A. H. MacDonald, Phys. Rev. B 84, 075418 (2011). Z. Qiao, J. Jung, Q. Niu and A. H. MacDonald, Nano Lett. 11, 3453 (2011). M. Zarenia, J. M. Pereira, Jr., G. A. Farias, and F. M. Peeters, Phys. Rev. B, 84, 125451 (2011). J. Li, I. Martin, M. Buttiker, and A. F. Morpurgo, Nature Phys. 7, 38 (2011). J. S. Alden, A. W. Tsen, P. Y. Huang, R. Hovden, L. Brown, J. Park, D. A. Muller, and P. L. McEuen, Proc. Natl. Acad. Sci. USA 110, 11256 (2013). A. Vaezi, Y. Liang, D. H. Ngai, L. Yang, and E. Kim, Phys. Rev. X 3, 021018 (2013). X. Li, F. Zhang, Q. Niu, and A. H. MacDonald, Phys. Rev. Lett. 113, 116803 (2014). M. Killi, T-C Wei, I. Affleck, and A. Paramekanti, Phys. Rev. Lett. 104, 216406 (2010). C. Y. Hou, E. A. Kim, and C. Chamon, Phys. Rev. Lett. 102, 076602 (2009). A. Str¨m and H. Johannesson, Phys. Rev. Lett. 102, o 096806 (2009). J. Teo and C. L. Kane. Phys. Rev. B 79, 235321 (2009). T. Giamarchi, Quantum Physics in One Dimension (Clarendon Press, Oxford, 2003). K. A. Matveev, D. Yue, and L. I. Glazman, Phys. Rev. Lett. 71, 3351 (1993); D. Yue, L. I. Glazman, and K. A. Matveev, Phys. Rev. B 49, 1966 (1994). L. Ju, Z. Shi, N. Nair, Y. Lv, C. Jin, J. Velasco Jr., C. Ojeda-Aristizabal, H. A. Bechtel, M. C. Martin, A. Zettl, J. Analytis, and F. Wang. Nature 520, 650 (2015). C. L. Kane and M. P. A. Fisher, Phys. Rev. B, 46, 15233 (1992).