Quantum Multicriticality near the Dirac-Semimetal to Band-Insulator Critical Point in Two Dimensions: A Controlled Ascent from One Dimension

We compute the effects of generic short-range interactions on gapless electrons residing at the quantum critical point separating a two-dimensional Dirac semimetal (DSM) and a symmetry-preserving band insulator (BI). The electronic dispersion at this critical point is anisotropic ($E_{\mathbf k}=\pm \sqrt{v^2 k^2_x + b^2 k^{2n}_y}$ with $n=2$), which results in unconventional scaling of physical observables. Due to the vanishing density of states ($\varrho(E) \sim |E|^{1/n}$), this anisotropic semimetal (ASM) is stable against weak short-range interactions. However, for stronger interactions the direct DSM-BI transition can either $(i)$ become a first-order transition, or $(ii)$ get avoided by an intervening broken-symmetry phase (BSP). We perform a renormalization group analysis by perturbing away from the one-dimensional limit with the small parameter $\epsilon = 1/n$, augmented with a $1/n$ expansion (parametrically suppressing quantum fluctuations in higher dimension). We identify charge density wave (CDW), antiferromagnet (AFM) and singlet s-wave superconductor as the three dominant candidates for the BSP. The onset of any such order at strong coupling $(\sim \epsilon)$ takes place through a continuous quantum phase transition across multicritical point. We also present the phase diagram of an extended Hubbard model for the ASM, obtained via the controlled deformation of its counterpart in one dimension. The latter displays spin-charge separation and instabilities to CDW, spin density wave, and Luther-Emery liquid phases at arbitrarily weak coupling. The spin density wave and Luther-Emery liquid phases deform into pseudospin SU(2)-symmetric quantum critical points separating the ASM from the AFM and superconducting orders, respectively. Our results can be germane for a uniaxially strained honeycomb lattice or organic compound $\alpha$-(BEDT-TTF)$_2\text{I}_3$.

We compute the effects of generic short-range interactions on gapless electrons residing at the quantum critical point separating a two-dimensional Dirac semimetal and a symmetry-preserving band insulator. The electronic dispersion at this critical point is anisotropic (E k = ± v 2 k 2 x + b 2 k 2n y with n = 2), which results in unconventional scaling of thermodynamic and transport quantities. Due to the vanishing density of states ( (E) ∼ |E| 1/n ), this anisotropic semimetal (ASM) is stable against weak short-range interactions. However, for stronger interactions the direct Dirac-semimetal to band-insulator transition can either (i) become a fluctuation-driven first-order transition (although unlikely in a particular microscopic model considered here, the anisotropic honeycomb lattice extended Hubbard model), or (ii) get avoided by an intervening broken-symmetry phase. We perform a controlled renormalization group analysis with the small parameter = 1/n, augmented with a 1/n expansion (parametrically suppressing quantum fluctuations in the higher dimension) by perturbing away from the one-dimensional limit, realized by setting = 0 and n → ∞. We identify charge density wave (CDW), antiferromagnet (AFM) and singlet s-wave superconductivity as the three dominant candidates for broken symmetry. The onset of any such order at strong coupling (∼ ) takes place through a continuous quantum phase transition across an interacting multicritical point, where the ordered phase, band insulator, Dirac and anisotropic semimetals meet. We also present the phase diagram of an extended Hubbard model for the ASM, obtained via the controlled deformation of its counterpart in one dimension. The latter displays spin-charge separation and instabilities to CDW, spin density wave, and Luther-Emery liquid phases at arbitrarily weak coupling. The spin density wave and Luther-Emery liquid phases deform into pseudospin SU(2)-symmetric quantum critical points separating the ASM from the AFM and superconducting orders, respectively. Our phase diagram shows an intriguing interplay among CDW, AFM and s-wave paired states that can be germane for a uniaxially strained optical honeycomb lattice for ultracold fermion atoms, or the organic compound α-(BEDT-TTF)2I3.

I. INTRODUCTION
A Dirac semimetal stands as a paradigmatic representative of a symmetry-protected gapless topological phase of matter that, for example, in two spatial dimensions can be realized in pristine monolayer graphene [1][2][3]. In a planar system, such a phase can be envisioned as a bound state of an equal number of vortices and antivortices (with unit vorticity) in reciprocal space, such that the net vorticity is zero. The difference of the vorticities, coined as the axial vorticity (N a ), is however finite and given by N a = 2n, where n is an integer (n = 1 for graphene), that (modulo 2) in turn also defines the integer topological invariant/charge of a twodimensional Dirac semimetal. Still it is conceivable to tune some suitable band parameter to the drive system through a continuous topological quantum phase transition where vortex and antivortex annihilate at a high symmetry point in the Brillouin zone, beyond which the system becomes a trivial band insulator. In the band insulator phase N a = 0. While a Dirac semimetal features linearly dispersing quasiparticles in all directions down to arbitrarily low energy, the critical fermions residing at the Dirac-semimetal to band-insulator quantum critical point (QCP) possess linear and quadratic dispersions along orthogonal directions in momentum space.
The quintessential properties of such a transition can be captured by a simple single-particle Hamiltonian H(k, ∆) = σ 0 vk x τ 1 + bk 2 y + ∆ τ 2 , where two sets of Pauli matrices σ µ = {σ 0 , σ 1 , σ 2 , σ 3 } and τ µ = {τ 0 , τ 1 , τ 2 , τ 3 } respectively operate on the spin and orbital space (and σ 0 = τ 0 is the identity), and we set = 1. Throughout this paper we assume that spin is a good quantum number (neglecting weak spinorbit coupling). In Eq. (1) v and b bear the dimensions of Fermi velocity and inverse mass, respectively, and ∆ has the dimension of energy. The above Hamiltonian represents (i) a Dirac semimetal for ∆ < 0, and (ii) a band insulator for ∆ > 0, as shown in Fig. 1(a). In the Dirac semimetal phase the Dirac points are located at k x = 0, k y = ± −∆/b. 1 The Dirac-semimetal to bandinsulator QCP is located at ∆ = 0, 2 where the quasiparticle spectra, given by E k = ± v 2 k 2 x + b 2 k 4 y , display both linear (along k x ) and quadratic (along k y ) dependence on the two different components of momenta. The , get gradually washed out as the temperature is increased and lattice effects set in at high temperature (T ∼ T0). The zero temperature quantum critical point in the noninteracting system, located at ∆ = 0, is represented by the red dot. The critical regime can also be exposed by frequency and magnetic field (see text). The scaling of physical observables inside the critical fan, within the Dirac semimetal, or within the band insulator side of the transition is reported in Table I. (b) Hopping pattern in a uniaxially strained honeycomb lattice. Strong (weak) hopping amplitudes are represented by thick (thin) lines. The sites on the two sublattices are shown in red (A) and blue (B). The anisotropic semimetal is realized when t2 = 2t1. density of states of the anisotropic semimetal (ASM) vanishes as (E) ∼ √ E. We compute key thermodynamic and transport properties of the noninteracting ASM [see Table I].
In this paper, we study the effects of generic shortrange electronic interactions on such an ASM. We assume that long-range Coulomb interactions [4,5] are screened, e.g. via a proximate gate. Due to the vanishing density of states, the ASM is stable against sufficiently weak shortrange interactions. By contrast, we show that the ASM can undergo a continuous quantum phase transition at strong interaction coupling through a multicritical point (tuned via interaction and anisotropy strengths) and enter into various broken-symmetry phases. At the multi- FIG. 2: Schematic phase diagrams of a two-dimensional Dirac material, residing in close proximity to the Dirac-semimetal to band-insulator quantum critical point (the red dot), in the presence of generic short-range interactions. Critical fermions possessing anisotropic dispersion are found along the black line in either subfigure. This direct transition can be avoided in two separate scenarios. Subfigure (a) depicts a strongcoupling scenario in which the critical anisotropic semimetal (ASM) is replaced by a line of first-order transitions (red dashed line). Subfigure (b) depicts an alternative scenario in which the ASM becomes unstable to the formation of a broken-symmetry phase. In this case, the blue dot represents a multicritical point separating the ASM from spontaneous ordering. Charge density wave, antiferromagnet (Néel) and s-wave superconductivity are the prominent candidates for the broken-symmetry phase [see also Fig. 4(b) and Fig. 12]. For weak coupling in (a), the only effect of the interactions is to shift the phase boundary separating the band insulator and Dirac semimetal. In (b), we assume that any such shift is compensated by the bare anisotropy (∆). We do not here discuss instabilities of the Dirac semimetal (see text for discussion on this issue [6]).
critical points the Dirac semimetal, band insulator, ASM and a broken symmetry phase meet [see the blue dot in Fig. 2(b)]. We identify charge density wave (CDW), antiferromagnet (AFM) and spin-singlet s-wave superconductivity as leading candidates for the broken symmetry phase. These conclusions obtain via the renormalization group (RG) controlled by an expansion about the one-dimensional limit, described below. We also consider an alternative scenario that could occur at strong coupling, wherein the Dirac semimetal can be separated from the band insulator by a fluctuation-driven firstorder transition. These two possibilities are schematically displayed in Fig. 2. We do not discuss in this paper quantum phase transitions between the Dirac semimetal and broken symmetry phases. In addition to CDW, AFM and s-wave pairing, a Dirac semimetal can accommodate additional fully gapped orders, such as quantum anomalous/spin-Hall insulator or Kekule valence-bond solid [6]. Both of these arise due to the valley degree of freedom in the Dirac system, which is annihilated at the Dirac-semimetal to band-insulator (anisotropic semimetal) transition. Since the density of states in a Dirac semimetal (ASM) vanishes as (E) ∼ |E| ( √ E), any ordering in a Dirac semimetal tuned close to the transition into the band insulator is expected to be preempted by those in the ASM, allowing us to focus solely on interaction effects in the latter. In a half-filled honeycomb lattice the ASM can be realized by applying uniaxial strain [7][8][9][10] [see Fig. 1(b)], which in real graphene would require an extremely large (and possibly unrealistic) distortion of the lattice. Nonetheless, in an optical honeycomb lattice for ultracold fermion atoms, the ASM can be achieved by tuning the depth of the laser trap, as has been reported in recent experiments [11][12][13]. In addition, the pressured organic compound α-(BEDT-TTF) 2 I 3 [14][15][16][17] as well as black phosphorous (a system with a few layers of phosphorene) [18][19][20] may reside very close to the ASM QCP. A recent ARPES measurement is suggestive of anisotropic dispersion in black phosphorus [19]. Furthermore, the quasiparticles residing at the interface of TiO 2 /VO 2 are also believed to possess such an anisotropic dispersion [21,22]. The QCP separating the Dirac semimetal and band insulator can also be accessed in uniaxially strained artificial/molecular graphene, as the hopping pattern in this system can be tuned quite efficiently [23].
A recent experiment has observed a first-order transition in α-(BEDT-TTF) 2 I 3 [24], also found in our theoretical analysis. Additionally, a semimetal-Mott (most likely AFM state) crossover in an optical honeycomb lattice driven by onsite repulsion has been reported in an experiment [12]. Therefore, our results are germane to a plethora of condensed matter systems, but also directly applicable for a strained optical honeycomb lattice (populated by neutral atoms). At the same time, our theoretical approach to the interacting ASM smoothly connects it to strong correlation physics in one dimension, and might serve as a platform to explore how 1D physics such as spin-charge separation and fractionalization get modified in higher dimensions. We now highlight the central results of our study.

A. Noninteracting system
We begin by noting some hallmark signatures of the noninteracting ASM, the scaling behavior of thermodynamic and transport quantities that can directly be observed in experiments. Even though the ASM can only be found at ∆ = 0 when T = 0, its imprint can be realized over the wide quantum critical regime shown in Fig. 1(a) at finite temperature (up to a temperature T 0 ∼ v 2 /(k B b) beyond which details of the lattice become important). The power-law scaling of various thermodynamic quantities, shown in Table I (center column), is essentially governed by the scaling of density of states, which in the ASM vanishes as (E) ∼ √ E. Consequently the compressibility scales as κ ∼ √ T and the specific heat vanishes as C v ∼ T 3/2 (at charge neutrality or for T |µ|, where µ denotes the chemical potential). If ∆ < 0, the scaling of thermodynamic quantities displays a smooth crossover from the ones in an ASM to those in a Dirac semimetal (displayed in the right column of Table I) around a crossover temperature T * ∼ |∆|/k B [25]. On the other hand, for ∆ > 0, the thermodynamic quantities display activated scaling for T < T * . The quantum critical regime shown in Fig. 1(a) can also be exposed by frequency (ω) as long as ω ω 0 ∼ v 2 /b and ω > |∆|. For example, the Drude conductivity for noninteracting electrons in the ASM scales as

Physical observable ASM DSM
while the interband component of the optical conductivity goes as where a x ≈ 0.06262, a y ≈ 0.20602, k B T 0 = ω 0 = v 2 /b. In the last two equations j = 1, 2 respectively corresponds to x and y. The scaling of the two universal functions F 1 (x) and F 2 (x) is shown in Fig. 3(a), and the detailed calculation of the optical conductivity is presented in Appendix A. Note that as frequency is lowered the interband Here µ is the chemical potential and T is the temperature. (b) Scaling of the function f (x) versus the dimensionless argument x (the ratio of thermal de Broglie wavelength to magnetic length) that controls the diamagnetic susceptibility of an anisotropic semimetal at finite temperature [see Eq. (5)]. component of the optical conductivity displays a smooth crossover from the scaling highlighted in Eq. (3) to the one for massless Dirac fermions (σ xx,yy ∼ e 2 /h without any leading power-law dependence on ω) for ω < |∆| and ∆ < 0. By contrast, the optical conductivity displays activated behavior at low frequency (ω < ∆) when ∆ > 0.
The diamagnetic susceptibility of the ASM at T = 0 is given by which diverges as B −1/3 , where A ≈ 0.075. A detailed calculation is presented in Appendix B. By contrast, the diamagnetic susceptibility of a two-dimensional Dirac semimetal diverges as χ 0 ∼ B −1/2 [26], while that for three-dimensional Dirac [27] and Weyl [28] semimetals scales as χ 0 ∼ log(B/B 0 ), where B 0 ∼ 1/a 2 and a is the lattice spacing. The power-law dependence of the diamagnetic susceptibility, shown in Eq. (4), can be observed in experiments only when B > ∆ 2 , so that we expose the critical regime associated with the noninteracting QCP [ Fig. 1(a)] by an external magnetic field. 3 Otherwise, with decreasing strength of the magnetic field the diamagnetic susceptibility displays a smooth crossover from B −1/3 to B −1/2 scaling when B < ∆ 2 and ∆ < 0. At finite temperature the diamagnetic susceptibility assumes the following universal scaling form (see Ap-pendix B for details) where 3/4 is the thermal de Broglie wavelength for the critical fermions. 4 The scaling of the universal function f (x) is shown in Fig. 3(b). Our proposed scaling at finite temperature (χ T ) is valid as long as B < λ T h < ξ ∆ ∼ 1/∆. 5

B. Electron-electron interactions
Next we turn to the effects of generic short-range electronic interactions on the ASM. We focus on the particlehole symmetric system at charge neutrality with µ = 0, although our results will also strongly influence the finite temperature physics when T |µ|. Scaling of the density of states yields a negative scaling dimension for any four-fermion interaction coupling g, namely [g] = −1/2. Consequently, the noninteracting ASM QCP separating the Dirac semimetal and the band insulator remains stable against sufficiently weak short-range interactions. However, at stronger interaction coupling we predict that one of two possible scenarios occurs: Either (i) the Dirac-semimetal to band-insulator transition becomes a fluctuation-driven first-order transition (see Sec. V A), or (ii) the direct transition gets masked by an intervening broken-symmetry phase (see Sec. VI). These two scenarios are schematically depicted in Fig. 2(a) and Fig. 2(b), respectively. Demonstrating the first-order transition requires a nonperturbative approach, which we here implement using large-N mean field theory (where N is the number of flavors of critical fermions; N = 2 for spin-1/2 electrons in the strained honeycomb model). The onset of broken symmetry takes place through a continuous quantum phase transition. The nucleation of broken symmetry through a continuous transition at intermediate coupling can be established by an RG analysis.
To control the RG calculation, we deform the Hamiltonian from Eq. (1) to where n is an even integer, so that this transformation does not alter the symmetry of the system (discussed 4 By contrast, the thermal de Broglie wavelength for relativistic and nonrelativistic fermions scales as λ T h ∼ T −1 and T −1/2 , respectively. 5 The scaling of the diamagnetic susceptibility for the anisotropic semimetal at finite temperature and weak magnetic field (no Landau quantization) can be estimated as follows: T h , yielding χ(T ) ∼ T −1/2 . On the other hand, for a two-dimensional Dirac semimetal the diamagnetic susceptibility at finite-T scales as χ(T ) ∼ B −1/2 ∼ B ∼ λ th , yielding χ(T ) ∼ T −1 . Here, U and V respectively correspond to onsite and nearest-neighbor interactions, see Eq. (8). The 1D extended Hubbard model has been extensively studied, see e.g. [32][33][34][35]; here we limit attention to the well-known predictions that obtain from bosonization [30,31], expected to hold for sufficiently weak coupling. In 1D, electronic quasiparticles are ill-defined for arbitrarily weak interactions. The interacting phases are shown in panel (a) and consist of (i) spin density wave (SDW), (ii) charge density wave (CDW), (iii) Luther-Emery liquid (LE), and (iv) coexisting spin and charge Luttinger liquids (LLs). In the continuum limit these phases exhibit complete spin-charge separation [31]; the charge is gapped in the CDW and SDW phases, while the spin is gapped in the CDW and LE phases. The latter is a 1D precursor to superconductivity. We obtain a qualitatively similar phase diagram for the ASM, shown in (b). The key differences are that the quasi-long-range-ordered SDW and LE phases are respectively replaced with true long-range-ordered antiferromagnetism (AFM) and s-wave superconductivity (SC) (at zero temperature), and that the CDW, AFM, and SC orders appear at finite coupling due to the vanishing density of states. In (b) we set n = 10 [Eq. (6)] to parametrically suppress higher dimensional quantum fluctuations (see Sec. V for a detailed discussion); here U and V are dimensionless couplings measured in units of . The phase diagram (b) is obtained by simultaneously solving the RG flow equations for four-fermion coupling constants and order parameter source terms (see Sec. VII). For V = 0 and sufficiently attractive U < 0, simultaneous nucleation of CDW and s-wave pairing occurs (red/cyan dashed line). Along this line pseudospin SU(2) symmetry [36][37][38] gets spontaneously broken. Various quantum critical fixed points FPj (see Table IV) control different segments of the phase boundaries, indicated in panel (b). Fixed points FP2,4 exhibit emergent pseudospin SU(2) symmetry, while this is broken at FP3. In the CDW phase (red region) in panel (b), the green dashed line is a crossover boundary. Above the line of coexistence with the CDW, the s-wave superconductivity exhibits short-range ordering that smoothly vanishes at the green dotted line. This is governed by the bicritical point FP5 (see Table IV), where the ASM-CDW phase boundary displays a kink (green dot).
in Sec. II B). The density of states in the deformed system scales as (E) ∼ E 1/n , and consequently the scaling dimension for any four-fermion interaction coupling g becomes [g] = −1/n. In the limit n → ∞, the density of states becomes finite. With purely local interactions, this limit corresponds to a decoupled collection of one-dimensional systems composed of massless Dirac fermions with spectra ε k = ±vk x . All local interactions have vanishing engineering dimension in this onedimensional limit. Hence, n → ∞ sets the marginality condition for generic short-range interactions in the ASM, and the perturbative RG calculation can be controlled in terms of a small parameter = 1/n, following the spirit of the -expansion [29] with [g] = − . Thus our RG calculation is performed in an effective dimension d * = 1 + , and precisely at = 0 local four-fermion interactions are marginal.
The above construction allows us to further control quantum fluctuations arising from the finite quasiparticle dispersion along the k y -direction by a parameter 1/n, entering through loop corrections in the perturbative RG calculation. Thus, our RG analysis is simultaneously controlled by two small parameters, (which sets the engineering dimension of the interaction coupling g) and 1/n (which controls the quantum fluctuations). In this framework, the RG flow equations take the form where l is the logarithm of the RG length scale, and the summation over µ, ν runs over all symmetry-allowed, linearly independent coupling constants {g µ }. To leading order, A µ,ν and B µ,ν are numerical (n-independent) matrices that we compute. We emphasize that we treat and 1/n as independent parameters in our RG scheme, although we should set these equal at the end. There are physical and technical reasons for this. Physically, enters as the engineering dimension for interaction couplings [g µ ] = − [Eq. (7)]. This is due entirely to dimensional analysis. While > 0 makes all such interactions irrelevant at weak coupling, it preserves key aspects of the physics special to one dimension, such as spin-charge separation. By contrast, the explicit 1/n corrections in Eq. (7) obtain from loop integrations involving the k y -dispersion, and thus encode quantum fluctuations beyond 1D. The Hubbard model phase diagram in Fig. 4(b) (discussed below) obtains only after quantum corrections are included, since these select between competing orders that are degenerate when > 0 but 1/n = 0 [see Table V and Sec. VIII for details]. Technically, our calculation works like an -expansion, 6 since we expand order-by-order in the interaction strengths. The expansion of the one-loop corrections (see Appendix E) in 1/n is not technically required, but is performed for physical clarity. 7 As a benchmark, we will show that if we set = 0 and send n → ∞ in Eq. (7), we recover well-established results for interacting fermions in one dimension. For example: (i) the β-function vanishes for spinless fermions (described by a two-component Dirac fermion and a single Luttinger current-current interaction coupling g), and (ii) the above set of flow equations displays spin-charge separation, with independent charge and spin Kosterlitz-Thouless transitions for spin-1/2 electrons [30,31].
The central outcome of our analysis can be summarized in the form of a phase diagram for an extended Hubbard model appropriate for uniaxially strained graphene, tuned to the Dirac-semimetal to band-insulator QCP. The interactions are defined via Here, U and V respectively denote the onsite and nearestneighbor (NN) interaction strengths, and n σ ( X) is the fermion number operator at position r = X with spin projection σ =↑, ↓. For the sake of simplicity, we assume that the strength of NN repulsion is the same among all three NN sites in uniaxially strained graphene. 8 The extended Hubbard model in one dimension has been exten- 6 An additional advantage of our -scheme is that it is performed in fixed d = 2 spatial dimensions; no analytic continuation of the Clifford algebra is required. 7 A comparison can be drawn between our combined -and 1/nexpansion with a conventional double-expansion in and 1/N , where N is the flavor number of the original problem [39]. Typically captures the deviation from the marginal dimension. However, the parameter 1/N only controls the quantum loop corrections and does not enter in the engineering dimension of the coupling [g] = − . By contrast, we here gain control over loop corrections by tuning the band curvature in the ky-direction (through the parameter 1/n, which does not enter the scaling dimension of g, namely [g] = − , in d * = 1+ ) instead of artificially increasing the flavor number. 8 The Dirac-semimetal to band-insulator QCP in the uniaxially strained honeycomb lattice is achieved at the cost of C 3v sym-sively studied by analytical [30,31] and numerical [32][33][34][35] methods (see also Sec. VIII A). The phase diagrams for one and two dimensional systems are respectively shown in Fig. 4(a) and Fig. 4(b). A detailed analysis of the extended Hubbard model is presented in Sec. VIII. The relevance of Eq. (8) and Fig. 4(b) to experiments is as follows. A strained optical honeycomb lattice for ultracold atoms [11][12][13] constitutes an ideal platform to realize the entire V = 0 axis of the phase diagram, as the strength of repulsive and attractive Hubbard interaction can be tuned quite efficiently, while NN interactions are negligible [40].
In pristine (unstrained) graphene U ∼ 9 eV, while V ∼ 5 eV for κ = 2.5, where κ is the effective dielectric constant of the medium [41]; the hopping strength t ∼ 2.8 eV. When uniaxially strained, the strength of onsite repulsion in graphene is likely unchanged, while the NN interaction gets slightly weaker, thus allowing access to the onsite-repulsion-dominated regime of the phase diagram. A very large strain of order 20% is required to tune to the Dirac-semimetal to band-insulator transition point [9]. A strain around 13% was recently achieved using MEMS (albeit in three-layer graphene [42]), corresponding to t 2 ∼ 1.5t 1 [9] (see Fig. 1(b)). This is still relatively far from the transition at t 2 = 2t 1 , but even larger values of strain may be possible. As emphasized above, it is not necessary to tune precisely to the transition point to potentially realize the Hubbard phase diagram in Fig. 4(b), since transitions out of the ASM are expected to dominate over those out of the Dirac semimetal due to the larger density of states of the former. The parabolic curvature responsible for density of states enhancement turns on at intermediate singleparticle energies even when there is still some splitting between degenerate Dirac semimetal valleys at zero energy (see the wide quantum critical regime for ASM in Fig. 1(a)), and intermediate energies are more relevant for U, V > t.
Another method to achieve the ASM in graphene could exploit a Moiré pattern. Such a pattern can occur in twisted bilayers [43], where it generates satellite Dirac points in the vicinity of the K and K points. A smaller (but still substantial) strain could be used to further merge two satellites along a high-symmetry direction.
The strength of U and V in black phosphorus, at a TiO 2 /VO 2 interface, and in α-(BEDT-TTF) 2 I 3 are not precisely known. Nonetheless, by changing the substrate on which these quasi-two-dimensional systems reside, one can efficiently tune the strength of NN repulsion, with suspended samples endowed with strongest NN repulsion. Furthermore, by changing the distance (ξ) between a gate metry (present in the absence of strain), see Fig. 1(b). Thus, the coupling V is expected to be different for the three NN sites. Such modification can only lead to a non-universal shift of the phase boundaries in Fig. 4(b), without altering its qualitative and universal features. and the sample one can also tune the strength of NN repulsion [44]. The screened Coulomb interaction is where the first (second) result applies to a single (double)-gated sample. In this equation U 0 ∼ e 2 /(κξ); the nearest neighbor repulsion is V = V (a), with a the lattice spacing. Thus different gate configurations can tune the relative strength of U and V , and experimentally access various repulsive regimes and quantum phase transitions of the phase diagram in Fig. 4(b). Black phosphorous (already residing close to the ASM fixed point) in particular can be tuned to the ASM band structure via impurity doping and the giant Stark effect [19], but this also shifts the chemical potential away from charge neutrality. The effects of finite chemical potential on the phase diagram shown in Fig. 4(b) are an important avenue for future work. We hope that the above discussion will motivate first principle calculations in other materials (similar to Ref. [41]) to estimate interaction strengths in various 2D ASMs, and future experiments to search for various broken symmetry phases. Superconductivity in an ASM can also be induced via proximity effect. Proximity-induced superconductivity has been recently achieved in monolayer graphene, when deposited on Rhenium (an s-wave superconductor), and Pr 2−x Ce x CuO 4 (a d-wave superconductor, with coherence length ξ S ∼ 30nm), with respective superconducting transition temperatures T c ∼ 2K [45] and T c ∼ 4.2K [46]. Therefore, it is conceivable to induce superconductivity in an ASM by controlling the distance (ξ) between strained graphene or black phosphorus and bulk superconducting materials (such as Rh, Pr 2−x Ce x CuO 4 ).
We now discuss salient features of the two phase diagrams shown in Figs. 4(a) and 4(b). In 1D, electronic quasiparticles are ill-defined for arbitrarily weak interactions. The 1D extended Hubbard phase diagram reviewed in Fig. 4(a) is that obtained via bosonization, which predicts spin-charge separation [30,31]. The spin and charge sectors can each separately exhibit gapless Luttinger liquid or gapped Mott insulating phases. This leads to the four composite phases shown in Fig. 4(a): (i) spin density wave (SDW), (ii) charge density wave (CDW), (iii) Luther-Emery liquid (LE), and (iv) coexisting spin and charge Luttinger liquids. The charge (spin) is gapped in the CDW and SDW (CDW and LE) phases. The SDW and Luther-Emery phases exhibit quasi-longrange order in the spin and charge sectors, respectively; these are precursors to antiferromagnetism and superconductivity in higher dimensions. The CDW phase exhibits true long-range, sublattice-staggered charge order at zero temperature.
The phase diagram of the extended Hubbard model in the two-dimensional ASM is shown in Fig. 4(b), as obtained via our RG analysis. The 2D ASM (yellow shaded region) remains stable against sufficiently weak local in-teractions. The onset of ordered phases takes place at a finite strength coupling (g µ ∼ ) in the ASM. In terms of the interacting phases, 1D [ Fig. 4 Fig. 4(b)] are qualitatively similar. The key difference is that the quasi-long-range-ordered SDW and LE phases in 1D are respectively replaced with true long-range-ordered Néel antiferromagnetism (AFM) and s-wave superconductivity (SC) in 2D (at zero temperature).

(a)] and 2D [
Sufficiently strong repulsive U > 0 drives the ASM through a continuous quantum phase transition and places it into the Néel AFM phase, which breaks spin SU(2) and time-reversal symmetry on the honeycomb lattice. Strong onsite attraction U < 0 with V = 0 induces the simultaneous nucleation of s-wave pairing and CDW orders, stemming from the spontaneous breaking of the exact pseudospin SU(2) symmetry. 9 Nonzero V breaks the pseudospin symmetry, preferring CDW (s-wave superconductivity) for repulsive V > 0 (attractive V < 0) coupling. The CDW and AFM phases respectively exhibit sublattice-staggered charge and collinear spin orders; they preserve translational invariance on the bipartite honeycomb lattice, but break reflection invariance across any line of bonds.
Thus the phase diagram of the extended Hubbard model in the ASM can be considered a controlled deformation (by two small parameters and 1/n) of its counterpart in one spatial dimension. In addition, the results obtained from our leading order RG analysis are consistent with pseudospin SU(2) symmetry [36][37][38] for V = 0. The fixed points FP 2 and FP 4 that govern the transitions along this line exhibit degenerate scaling dimensions for all operators residing in each irreducible representation of this symmetry (see Table V and Fig. 11). The proposed phase diagram in Fig. 4(b) should therefore qualitatively describe the strongly interacting ASM, residing at the phase boundary between the Dirac semimetal and band insulator, for which = 1 2 and n = 2. Approaching from within the s-wave SC and AFM broken symmetry phases, the pairing amplitude and Néel order parameters vanish continuously at the quantum critical points FP 2 and FP 4 , respectively. These transitions can therefore be regarded as the two-dimensional deformations of the 1D quasi-long-range-ordered LE and SDW phases. To the leading order in the -and 1/nexpansion the correlation length exponent ν −1 = is the same across all continuous transitions to broken symmetry phases. 10 This is, however, an artifact of the leading order calculation. Setting = 1/2, the transition or crossover temperature (T c ) between broken symmetry phases and the ASM scales as T c ∼ δ 2 , while that to the Dirac semimetal scales as T c ∼ δ, where δ = (X − X c ) /X c is the reduced control parameter for the zero-temperature quantum phase transition from ASM or Dirac semimetal into the ordered phase, with X = U, V [Eq. (8)] for example, and X c as the critical strength for ordering. Thus, even though the ASM lives adjacent to a Dirac semimetal, the scaling of transition/crossover temperature to these scale distinctly. The quantum phase transition in the former system should precede the one in the Dirac semimetal if the latter is tuned close to the transition into the band insulator, since the ASM density of states vanishes more slowly with energy.
Each of the interacting QCPs (discussed in Sec. VI), controlling various continuous transitions into broken symmetry phases are expected to accommodate strongly interacting non-Fermi liquids that lacks sharp quasiparticle excitations. It would be extremely interesting to look for remnants of spin-charge separation in the twoloop self-energy, which will give an anomalous dimension and lifetime to the fermion field. Our proposed phase diagram and associated quantum critical phenomena can be directly tested by quantum Monte Carlo simulations on the honeycomb lattice with onsite (U ) and nearest-neighbor (V ) interaction, as the extended Hubbard model can now be simulated without encountering the infamous sign problem [47][48][49][50][51][52][53][54]. The phase diagram (or at least some part of it) can also be exposed in future in optical honeycomb lattice experiments with ultracold fermion atoms, wherein the strength of interactions can possibly be efficiently tuned. Our methodology can also be adopted to investigate the effects of electronic interactions in various other itinerant systems, possessing anisotropic fermionic dispersion, such as general Weyl semimetals [55], which we discuss qualitatively in Sec. IX.

C. Outline
The rest of the paper is organized as follows. In Sec. II we establish the field theory description of the noninteracting Dirac-semimetal to band-insulator quantum phase transition and discuss its symmetries. Possible brokensymmetry phases (both excitonic and superconducting) and the associated quasiparticle spectra are discussed in Sec. III. In Sec. IV we describe the minimal models for spinless and spin-1/2 versions of the interacting ASM. In Sec. V we determine the effects of electronic interactions on the ASM for spinless fermions. The RG analysis for spin-1/2 fermions, connection to 1D spincharge separation, existence of various interacting QCPs and associated quantum critical phenomena are discussed in Sec. VI. In Sec. VII we determine the nature of the broken-symmetry phases at strong coupling across various continuous transitions for the interacting spin-1/2 model. Sec. VIII provides a more detailed discussion of the extended Hubbard model [Eq. (8)]. We summarize our findings, discuss applications of our methodology in other correlated systems and highlight prospects for future work in Sec. IX. Additional technical details are relegated to the Appendices.

II. NONINTERACTING SYSTEM: FIELD THEORY AND SYMMETRIES
In this section we construct the field theoretic description of the QCP separating a two-dimensional Dirac semimetal and a band insulator in a noninteracting system, and discuss the symmetries of the ASM separating these two phases.

A. Hamiltonian, Lagrangian, and scaling
The universality class of the Dirac-semimetal to bandinsulator quantum phase transition is captured by the quadratic Hamiltonian from Eq. (1) for ∆ = 0. The four-component spinor basis can be chosen as Ψ where c X,k,σ is the fermion annihilation operator with momentum k, spin projection σ =↑, ↓ and on sublattice X = A/B of the honeycomb lattice. In other compounds, such as black phosphorus or α-(BEDT-TTF) 2 I 3 , X represents orbital degrees of freedom.
The imaginary time (t) action associated with the noninteracting Hamiltonian reads where r ≡ (x, y) and the kinetic operator L 0 is given by Here we define L 0 via the deformation H(k, ∆) → H n (k, ∆) [see Eq. (6)], which leaves all the symmetries of the Hamiltonian unaffected (since n is here an even integer). We now introduce the notion of coarse-graining for this model, under which we send t → e l t, x → e l x, and y → e l/n y; the field transforms as Ψ → e −l(1+1/n)/2 Ψ. Here l denotes the logarithm of the renormalization group length scale. The imaginary time action S 0 remains invariant, so that the Fermi velocity (v) and inverse mass (b) are effectively dimensionless (carry zero scaling dimension [v] = [b] = 0). On the other hand, the scaling dimension of ∆ is [∆] = 1. This is a relevant perturbation at the Dirac-semimetal to band-insulator QCP that controls the transition between these two phases that possess the same symmetry. The scaling dimension of ∆ yields the correlation length exponent for this noninteracting QCP: ν −1 = 1.

B. Symmetries
The Hamiltonian describing the ASM possesses a discrete chiral symmetry [also known as sublattice symmetry (SLS)], generated by a unitary operator σ 0 τ 3 and encoded in the condition Time-reversal symmetry is generated by the antiunitary operator T = σ 2 τ 3 K, where K is complex conjugation, and T 2 = −1. The Hamiltonian H n (k, ∆) also possesses spin SU(2) symmetry generated by S = στ 0 , and x-reflection symmetry defined via This symmetry operation encodes invariance of the system under the exchange of two sublattices or orbitals, and will be denoted by R π [see Table II] [56,57]. We relegate a detailed discussion on the symmetry properties of the ASM to Appendix C. In the next section we discuss various possible ordered phases that break at least one of the above symmetries of the noninteracting system.

III. BROKEN-SYMMETRY PHASES
Since the DOS of the 2D ASM vanishes as (E) ∼ √ E, the noninteracting QCP separating a Dirac semimetal and a band insulator remains stable against sufficiently weak short-range interactions. Nevertheless, beyond a critical strength of the interactions the ASM can become unstable toward the formation of various brokensymmetry phases. Here we consider all possible phases, including both particle-hole or excitonic as well as particle-particle or superconducting orders. However, we restrict ourselves to the momentum-independent (local or intra-unit-cell) orderings. Then the ASM all together supports eight excitonic and four pairing orders, enumerated in Table II. We distinguish these with the labels in the order parameter ("OP") column of this table. The transformation of these fermion bilinears under the symmetry generators of the noninteracting system are also highlighted in Table II. One can add eight particle-hole channel fermion bilinears to the noninteracting Hamiltonian. Of these, only the anisotropy parameter ∆ s 2 that tunes through the Dirac-semimetal to band-insulator transition preserves all discrete and continuous symmetries. The charge density ∆ s 0 , x-current ∆ s 1 , and x-spin-current ∆ t 1 each break at least one discrete or continuous symmetry, but these conserved symmetry currents do not represent true order parameters. The remaining four bilinears are order parameters for the following broken symmetry phases: sublattice-staggered charge density wave (CDW) ∆ s 3 , sublattice-staggered collinear (Néel) antiferromagnet No ∆t TABLE II: All possible local (momentum independent or intra unit-cell) order parameters in the particle-hole or excitonic channel (first 8 rows) and particle-particle or superconducting channel (last 4 rows), and their transformation properties under various discrete and continuous symmetry operations, discussed in Sec. II B. Under spin SU (2), bilinears transform in either the singlet (0) or triplet (1) representations. Note that the anisotropy parameter (AP) preserves all microscopic symmetries, and respectively gives rise to gapless and gapped spectra for ∆ < 0 and ∆ > 0 on the corresponding sides of the Dirac-semimetal to band-insulator transition. The second to the last column shows whether the quasiparticle spectrum inside a given broken-symmetry phase is fully gapped or not. Notice that only three ordered phases, namely the charge density wave (CDW), Néel antiferromagnet (AFM), and s-wave superconductor (SC) yield a fully gapped spectrum. These three orders are energetically most favored, at least within the framework of the weak coupling renormalization group analysis [see, for example, the phase diagram of the extended Hubbard model in Fig. 4(b)]. The last column labels the corresponding order-parameter (OP) source term [see also Sec. VII]. The density, x-current, and x-spin-current are included here for completeness, although these are not true order parameters. We note that due to the lack of valley degrees of freedom a two-dimensional anisotropic semimetal does not allow any topological order such as quantum anomalous/spin Hall insulator, which are odd under the exchange of two valleys, unlike the situation in a Dirac semimetal [6].
Only the CDW and AFM order parameters induce a fully gapped spectrum in the ordered phase. We therefore anticipate CDW and AFM to be the dominant orders for the ASM with strong repulsive interactions [see Figs. 4(b) and 12], since they optimally minimize the free energy (at least at T = 0). Nevertheless, it is worth understanding the effects of other fermionic bilinears in the ordered phase. For example, a nonzero vacuum expectation value for ∆ s 1 renormalizes the Fermi velocity in the x-direction, while a ferromagnetic compensated ASM is realized for ∆ t 0 = 0. In this phase, perfectly nested electron-like and hole-like Fermi surfaces arise for opposite projections of the electron spin. Thus, such a compensated semimetallic phase can undergo a subsequent BCS-like weak coupling instability towards the formation of an AFM phase [58][59][60][61]. However, the order parameter in this AFM phase will be locked in the spin-easy-plane, perpendicular to the ferromagnetic moment in the compensated semimetal. For ∆ t 1 = 0, the Fermi velocity in the x-direction acquires opposite corrections for the two projections of electronic spin, while a nonzero value of ∆ t 2 yields both the Dirac semimetal and band insulator, but for opposite spin projections.
Sufficiently strong attractive interactions could induce various superconducting phases. When we restrict ourselves to local (momentum independent or intra-unitcell) pairings, all together the ASM supports only three spin-singlet and one spin-triplet pairings. Among them only the singlet s-wave pairing ∆ s gives rise to a fully gapped quasiparticle spectrum. Thus we expect rest of the pairings to be energetically inferior to the singlet s-wave superconductor. We will demonstrate that for strong attractive interactions the ASM becomes unstable towards the formation of this s-wave state [see Figs. 4(b) and 12].

IV. INTERACTING THEORY: MINIMAL MODEL
We now consider the minimal model of electronic interactions in the 2D ASM separating the Dirac semimetal and band insulator, compatible with the symmetry of the system, discussed in Sec. II B. We discuss spinless and spin-1/2 fermions separately. 11

A. Spinless fermions
For spinless fermions, generic local four-fermion interactions are described by four quartic terms, In this notation, g µ > 0 represents a repulsive electron-electron interaction that promotes condensation of the corresponding Hermitian fermion bilinear ψ † τ µ ψ . However, not all quartic terms are linearly independent. The Fierz identity mandates that generic interactions for spinless fermions can be captured by only 11 For the clarity of presentation we here denote the spinor for spinless (two-component) and spinful (four-component) fermions as ψ and Ψ, respectively.
one local four-fermion interaction (see Appendix F 1). Nevertheless, the final outcome regarding the fate of the critical fermions in the ASM at strong coupling depends upon the bilinear used to decouple the interactions (in mean field theory), as discussed in Sec. V.
Since we throughout assume spin to be a good quantum number, a collection of interacting spin-1/2 electrons in the ASM in the presence of generic short-range interactions is described in terms of eight quartic terms Again, not all quartic terms are independent. As shown in Appendix F 2, a maximum of four independent interactions are possible for spin-1/2 fermions. With a judicious choice (see Sec. VI) we define the interacting Hamiltonian as In Sec. VI we will perform an RG calculation with H int to understand the effect of electronic interactions in the spinful 2D ASM. It is instructive to express the interacting theory in terms of a different set of four couplings that will allow us to connect with well-established results for spin-1/2 electrons in one dimension, which in our formalism is achieved by setting = 0 and taking n → ∞ in the RG analysis. We can cast the interacting Hamiltonian H int in terms of four-fermion interactions that capture the physics of spin-1/2 electrons in one dimension most efficiently, according to where U N , U A , W , and X are coupling constants respectively for (i) an SU(2) 1 spin current-current interac- The precise definition of these operators is relegated to Appendix D, see Eq. (D2).
The correspondence between the four coupling constants, appearing in Eqs. (17) and (16) is or This mapping will prove useful to connect with the physics of interacting electrons in 1D, as discussed in Sec. VI.

V. INTERACTING SPINLESS FERMIONS
We first discuss the effects of electronic interactions on spinless fermions. As shown in the previous section, generic local interactions for spinless fermions are described by only one quartic term. However, we will demonstrate that for stronger interaction the direct transition between the Dirac semimetal and band insulator can either (i) get replaced by a first-order transition or (ii) get avoided by an intervening broken-symmetry phase, depending on how we address the effects of strong electron-electron interactions.
A. Topological first-order transition We first illustrate the possibility of the first-order transition between the Dirac semimetal and band insulator at strong interaction, for which we consider the interacting Hamiltonian with local interaction g 2 [see Eq. (14)]. This outcome can be demonstrated from the following mean-field or large-N free energy where N is the number of two-component spinors, obtained after performing a Hubbard-Stratonovich transformation of the four-fermion interaction proportional to g 2 in favor of a bosonic field Σ = ψ † τ 2 ψ , and subsequently integrating out the critical fermions. The parameter N (flavor number) should not be confused with n [controlling the degree of anisotropy in ASM, Eq. (6)].
To proceed with the calculation, let us now define a set of variables according to vk x = ρ cos θ, bk 2 y = ρ sin θ, with 0 ≤ θ ≤ π and 0 < ρ < E Λ , where E Λ is the highenergy cutoff up to which the quasiparticle dispersion is anisotropic (given by E k ). In terms of dimensionless variables, defined as (here σ should not be confused with Pauli matrices) and the dimensionless free-energy density f = F (2π Λ ), we can express the above equation as We numerically minimize the free energy and the resulting phase diagram is shown in Fig. 5(a). Nonetheless, the salient features of the phase diagram can be appreciated by expanding the above free energy in powers of σ, after shifting the variable σ +m → σ. It should, however, be noted that such expansion of the free-energy density in powers of σ is not a well-defined procedure, and for j > 4 the f j s are non-analytic functions. But, most importantly, the free energy f contains all odd powers of σ.
When the interaction is sufficiently weak the σ field does not condense, and the profile of the free energy possesses a single global minimum, as shown in Fig. 5(b). Consequently, the nature of the direct, continuous transition between the Dirac semimetal and band insulator remains unchanged at weak coupling. In this regime only the term linear in σ is important and δ = λ 2 f 1 defines the phase boundary between these two phases; this is represented by the black solid line in Fig. 5(a).
On the other hand, beyond a critical strength of the interaction λ 2 > λ * 2 ∼ 1, the σ field acquires a finite vacuum expectation value and all powers (including the odd ones) of σ in f become important. In this regime the profile of the free energy accommodates two inequivalent minima, as shown in Fig. 5(c) and Fig. 5(d), respectively on the band-insulator and Dirac-semimetal sides of the phase diagram. As a result, the direct transition between these two phases at strong coupling becomes a fluctuation-driven first-order transition along the red dashed line in Fig. 5(a) 12 . When the Dirac semimetal is separated from the band insulator by a first-order transition, there exists no ASM at the transition. A similar first-order transition at strong coupling has also been predicted in two-and three-dimensional topological insulators [63][64][65][66] and also in a three-dimensional Weyl material residing in proximity to the semimetal to bandinsulator QCP [55].
It is worth mentioning that a first-order transition has recently been observed in the pressured organic compound α-(BEDT-TTF) 2 I 3 [24]. Thus our theory provides a possible explanation of this experimental observation.

B. Continuous quantum phase transition
Next we will demonstrate the second possibility regarding the fate of the direct transition between the 12 Even though we here do not account for any local or spatial fluctuation of the order parameter Σ or σ field, the terminology "fluctuation driven first-order transition" is justified in this case since it can only be captured upon incorporating quantum corrections to the free-energy density (f ) due to electronic interaction g 2 or λ 2 . In particular, in the large-N (flavor number) limit such corrections arises from the Feynman diagrams (b) and (d), shown in Fig. 6. See also Ref. [62] for discussion on similar issue in one spatial dimension. Dirac semimetal and band insulator, when strong interactions causes nucleation of a broken symmetry phase that masks the direct transition between these two symmetrypreserving phases [see Fig. 2(b)]. For spinless fermions, only the CDW phase gives rise to a fully gapped spectrum; inside this ordered phase ψ † τ 3 ψ = 0. Thus CDW order is energetically advantageous at temperature T = 0 for sufficiently large g 3 > 0 [Eq. (14)]. Next we will establish the onset of such a CDW phase by carrying out a RG calculation on the interacting model with local interaction g 3 . The imaginary-time action for such interacting system reads After accounting for perturbative corrections to the quadratic order in g 3 (computing the one-loop diagrams shown in Fig. 6) and integrating out the fast Fourier modes within the Wilsonian shell E Λ e −l < ω 2 + v 2 k 2 x < E Λ and 0 < |k y | < ∞, we arrive at the following RG flow equations, where∆ = ∆/E Λ and λ 3 = g 3 E Λ /(8π 2 vb ) are dimensionless parameters. Here E Λ is the energy cutoff for ω and v|k x |, and ω is the Matsubara frequency. The functions appearing in the RG flow equations are (see Appendix E for details) Therefore, as n → ∞ only the contribution from f 1 (n) survives, while f 2 (n) and f 3 (n) capture subleading logarithmic divergences. The above coupled flow equations support only two fixed points: 1. ∆ , λ 3 = (0, 0) corresponds to the noninteracting system, the ASM critical point that connects the Dirac semimetal to the symmetry-preserving band insulator.
, which, on the other hand, stands as an interacting multicritical point in the ∆ , λ 3 plane, where the Dirac-semimetal, bandinsulator, ASM and CDW phases meet. Therefore, this point controls a continuous quantum phase transition of critical excitations residing at the phase boundary between the Dirac semimetal and band insulator towards the formation of the CDW phase. The RG flow and the phase diagram are respectively shown in Fig. 7(a) and Fig. 7(b).
We now show how to compute the location of the interacting multicritical point from Eq. 25. Respectively, d∆/dl = 0 and dλ 3 /dl = 0 leads to where quantities with subscript " * " denotes their fixed point values. Since, λ 3, * ∼ , for a systematic extraction of∆ * to the leading order in (interaction mediated shift of the band parameter), we only need to account for the first term in the above expression for∆ * , yielding [67] ∆ leading to the result announced above. The remaining terms in the expression of∆ * from Eq. (27) yields contributions at least O( 2 ). Therefore, while determining the location of the multicritical point, RG flow and phase diagram in the ∆ , λ 3 plane, we only accounted for the terms proportional to∆ and λ 3 , and neglected the contribution proportional to∆λ 3 in d∆/dl. We follow the same procedure for spin-1/2 electrons discussed in the next section.
From the leading order in -and 1/n-expansions we can only estimate the correlation length exponent (ν) for the ASM-CDW continuous transition, yielding where λ * 3 = n 4 represents the strength of four-fermion interactions at the multicritical point. Hence, for the physically relevant situation with = 1/2 and n = 2, we obtain ν = 2. On the other hand, when the Dirac semimetal undergoes a continuous transition into the CDW phase ν = 1 (to one-loop order) [56]. Thus, the universality class of the critical excitation-CDW transition is fundamentally a new one.
However, we note that the multicritical point describing the transition is not strictly controlled by a small parameter as λ * 3 ∼ n ∼ 1. This outcome is specific (and to certain extent pathological) for spinless two-component fermions. For spin-1/2 electrons, we show in the next section that all of the interacting QCPs and the associated critical properties are simultaneously controlled by two small parameters and 1/n.
The spinless case gives an order one coupling strength at the multicritical point because the correction in Eq. (25b) vanishes in the n → ∞ limit. This is true to all orders in λ 3 , since the interaction in Eq. (24) is equivalent to a U(1) current-current (Thirring [29]) perturbation. This is an exactly marginal Luttinger interaction in the 1D limit [30,31], which obtains for = 0 and n → ∞. Thus the ASM-CDW quantum phase transition is solely captured by the subleading divergence arising from f 2 (n), which can directly be tested at least in quantum Monte Carlo simulation with only NN interaction in a uniaxially stressed honeycomb lattice [53,54] [with t 2 = 2t 1 , see Fig. 1(b)]. By contrast, for isotropic two-component massless Dirac fermions, the RG flow equation vanishes to at least order λ 3 3 [68], although it is believed that there exists a possibly continuous semimetal-insulator transition at finite interaction strength [69,70]. Thus far this issue remains unsettled. On the other hand, our leading order RG calculation predicts a continuous anisotropic semimetal-CDW transition.

C. First vs. second-order transitions
We have identified two incompatible scenarios for the strong coupling physics in the spinless case, wherein the ASM is either replaced by a first-order transition, or the CDW phase. These conclusions are reached via calculations with different control parameters: the number of flavors N for the first-order transition, and the curvature of the dispersion along k y (characterized by the n) for the intervening CDW phase. Since we must set N = 1 and n = 2 for the spinless model, both results obtain at strong coupling and we cannot make a rigorous mathematical argument that one scenario is more likely. However, microscopic considerations provide some physical intuition which can favor one or the other.
For a spinless ASM obtained from a microscopic strained honeycomb lattice model that is free of frustration (e.g., a "t-V " model with nearest-neighbor hopping and density-density interactions only), the CDW order is expected to preempt the topological first-order transition between the Dirac semimetal and band insulator, since the CDW order completely gaps out critical excitations, producing a uniform mass gap in the spectrum. This can be tested in a quantum Monte Carlo simulation [53,54]. However, we cannot completely rule out the first-order transition, since bandwidth renormalization and/or the suppression of the quasiparticle weight may become crucial to determine the ultimate fate of the strongly interacting ASM. These effects are not included in the one-loop RG.
For spin-1/2 electrons generic local interactions live in a four-dimensional coupling constant space (see the next section). We can again gain some intuition by considering a particular microscopic model, e.g. the strained honeycomb lattice extended Hubbard model [see Fig. 4(b)]. Again this model is "frustration free," and the fully gapped mass orders CDW, AFM, and s-wave pairing are expected to occur at strong coupling. Nucleation of these mass orders via a continuous phase transition is energetically superior over the first order transition out of the Dirac semimetal, since the ordered phases fully gap the spectrum (maximal gain of condensation energy).
We therefore conclude that if the extended Hub-bard model can be realized in strained graphene or an anisotropic honeycomb optical lattice for ultracold fermions, the first-order transition will be unlikely. For a material with different microscopics [e.g. the organic compound α-(BEDT-TTF) 2 I 3 ], the Dirac-semimetal to band-insulator transition could occur at strong coupling in a corner of the four-dimensional interaction parameter space where the first-order transition is actually preferred. It is a nontrivial task to identify a microscopic model that can exhibit the proposed first order transition in this high-dimensional coupling space. Finally, it must be emphasized that the first-order transition we discuss here between two symmetry-preserving phases (Dirac semimetal and band insulator) is fundamentally different from the one between two distinct broken symmetry phases. Therefore, such a fluctuation-driven first-order transition cannot be captured by an RG calculation and the computation of free-energy density (non-perturbative analysis in the large-N limit) is necessary [62].

VI. SPIN-1/2 ELECTRONS
Next we perform the RG calculation to the leading order in and 1/n starting from the interacting H int [defined in Eq. (16)] for spin-1/2 electrons. We made a judicious choice in selecting the linearly independent coupling constants in H int . Note that we did not choose g s 2 as one of the four independent couplings, since as we know when this coupling constant is strong enough the Dirac-semimetal to band-insulator quantum phase transition becomes a fluctuation-driven first-order transition [following the analysis presented in Sec. V A but now for spinful electrons, which can be accomplished by taking N → 2N ]. Nonetheless, as we will discuss in the next section, the source term ∆ s 2 (corresponding to the anisotropy parameter) never displays the leading divergence. Thus, we can proceed with the following RG analysis tailored to address continuous phase transitions into broken symmetry phases, leaving aside the possibility of a first-order transition as an alternative scenario. Following the same procedure described for spinless fermions in Sec. V B (see Fig. 6 for the relevant Feynman diagrams) we arrive at the following RG equations where λ a µ = g a µ E Λ /(π 2 vb ) is the dimensionless coupling constant and∆ = ∆/E Λ is the dimensionless band parameter, and recall that here n can only take even integer values. The above set of RG flow equations can also be expressed in terms of the coupling constants introduced in Eq. (17) (namely U N , U A , W and X) according to In the last set of equations all coupling constants are dimensionless. Next we systematically analyze the above set of RG flow equations. We addressed the effect of RG flow of∆ in the phase diagram of the interacting ASM for the spinless case. Unless otherwise mentioned we will work in the hyperplane defined by∆(l) = 0, where∆(l) is the renormalized band parameter. Therefore, all interacting fixed points with λ s,t j = 0 [see Tables III and IV, and Fig. 9] have one additional unstable direction, which is∆, and all of them are truly multi-critical in nature. The functions f 1 (n), f 2 (n) and f 3 (n) have already been defined in Eq. (26).
A. Emergent one-dimensional system: spin-charge separation In the limit n → ∞, Eq. (31) reduces to If we set the engineering dimension of all the four-fermion couplings = 0, then we recover the well-known RG flow equations for these coupling strengths in one dimension. In this limit, the anisotropic semimetal with local interactions reduces to a decoupled collection of spin-1/2 onedimensional "wires". Such systems exhibit spin-charge separation [30,31], and this is reflected in Eq. (32). The only spin sector interaction U N couples to the SU(2) 1 current-current perturbation [see Eq. (D2) and the surrounding text for a precise definition]. For = 0, U N resides on the inflowing (outflowing) part of the spin sector Kosterlitz-Thouless separatrix for U N < 0 (U N > 0) [30,31]; spin SU(2) symmetry is preserved everywhere along the separatrix.
The remaining parameters couple to interaction operators that perturb the charge sector, as can be seen from the bosonization of the latter [Eq. (D3)]. U A and X couple to the U (1) current-current and U (1) stress tensor operators, respectively; these modify only the Luttinger parameter and charge velocity of the free boson description, as shown in Eq. (D4). By contrast, the umklapp interaction W couples to a sine-Gordon perturbation. The flow equations for U A and W with = 0 are the charge sector Kosterlitz-Thouless equations. The RG flow in the (U A , W ) plane is shown in Fig. 8.
Since spin and charge are independent, Eq. (32) with = 0 implies that either sector can be a gapless Luttinger liquid or gapped Mott insulator. The spin sector becomes massive when U N → +∞, while the charge sector becomes massive when U A → +∞ and W → ±∞. The product of these gives four different composite phases that must be interpreted through a microscopic model. The one-dimensional extended Hubbard model described by the phase diagram in Fig. 4(a) is an example; it is further reviewed in Sec. VIII, below.
The gapless charge sector Luttinger liquid phase is a stable line of fixed points, shown in Fig. 8. The tuning parameter is the interaction U A or equivalently, the charge sector Luttinger parameter K c [Eq. (D4)]. Up to logarithmic corrections [30] that we ignore, the spin Luttinger liquid is a fixed point with U N = 0 and spin sector Luttinger parameter K s = 1.
B. From 1D to 1 + Moving beyond the strict 1D limit, it is instructive to reinstate the dispersion in the y-direction in two steps. First, we consider the influence of the nonzero scaling dimension (− ) for all local four-fermion interactions in Eq. (32). With > 0, these RG flow equations give birth to six fixed points, tabulated in Table III and schematically shown in Fig. 9(a). All of them are located in a hyperplane defined by X = 0 or λ s 1 = 0. Note that inclusion of the nonzero scaling dimension does not destroy spin-charge separation.
The gapless charge sector Luttinger liquid phase [stable fixed line in Fig. 8] gets replaced by a stable noninteracting fixed point and two QCPs, as shown in Fig. 9(b). An additional quantum critical point appears in the spin sector at U N = .
Now we characterize all six fixed points in the hyper-plane∆(l) = 0, cataloged in Table III. The fixed point FP 1 represents the noninteracting ASM, stable against    (32). The noninteracting semimetal is labeled FP1 and has only one unstable direction; it controls the direct quantum phase transition between the Dirac semimetal and band insulator for sufficiently weak interactions. By contrast, each of FP2,3,4 has two unstable directions and FP5,6 has three unstable directions, and these correspond to interacting multicritical points. The penultimate rightmost column shows the value of the band parameter∆ at each such fixed point. If we follow the RG trajectory along which the renormalized band parameter is kept fixed, i.e.,∆ = 0 then FP2,3,4 represent interacting quantum critical points with only one unstable direction and FP5,6 become bicritical points with two unstable directions. A schematic plot of these fixed points is shown in Fig. 9. The last column shows the number of unstable directions (UDs) at a given fixed point in the∆ = 0 hyperplane as well as in the five dimensional coupling constant space (shown in parentheses).
sufficiently weak generic short-range interactions. Each of FP 2 , FP 3 and FP 4 is characterized by only one unstable direction (in addition to∆). These three fixed points correspond to interacting QCPs and are analogous to the interacting multicritical point shown in Fig. 7. They describe continuous transitions from the ASM to broken symmetry phases.
To determine the actual nature of the broken symmetry across the various QCPs in the ASM requires inclusion of 1/n corrections to the RG flow equations (30) or (31). These corrections are due to quantum fluctuations beyond one dimension and eliminate certain special symmetries specific to 1D, such as spin-charge separation. In Sec. VII, we will combine the RG with a scaling analysis of fermion bilinear susceptibilities in order to pin the pattern of symmetry breaking at strong coupling in an unbiased fashion.
The remaining two fixed points FP 5 and FP 6 each possess two unstable directions, as shown in Fig. 9(a). In thê ∆(l) = 0 hyperplane they represent interacting bicritical points, separating the basins of attraction for the various interacting QCPs. As shown in Fig. 9(a), FP 5 separates the domain of attraction for FP 2 and FP 3 , while FP 6 separates FP 4 and FP 2 . The basins of attraction for FP 3 and FP 4 are separated by FP 1 , which can be seen more transparently in Fig. 9(b).  Table III in the n → ∞ limit. It is quite challenging to find the location of these fixed points analytically. Nevertheless, upon numerically locating the fixed points for various values of n, we can extract the functional dependence on 1/n. At least for large n (namely for n ≥ 8) we find that the position of all the fixed points are well approximated functions of 1/n, and they are tabulated in Table IV. All fixed points are located in the hyperplane λ s 1 = 0. Therefore, this coupling constant does not change any outcome qualitatively or quantitatively.
The location of the six fixed points in the four dimensional coupling constant space spanned by (λ s 1 , λ s 3 , λ t 2 , λ t 3 ) or (U N , U A , W, X) are presented in Table IV. The nature of all the fixed points FP j for j = 1, · · · , 6 remains unchanged from that discussed in the previous section after 1/n corrections are included (see the last column of Table IV). From the leading order in and 1/n expansion we can also allude to some emergent quantum critical phenomena at the interacting QCPs (namely FP 2,3,4 ), describing continuous transitions from an ASM to various broken symmetry phases. For example, the correlation length exponent at the interacting QCPs is given by ν −1 = , and thus for n = 2 (physical situation) ν = 2, similar to the situation discussed previously for spinless fermions [see Eq. (29) for definition of ν]. The fact that correlation length exponent is same at all QCPs is, however, likely an artifact of the leading order calculation.
The residue of the quasiparticle pole of critical excitations (determined by fermionic anomalous dimension η Ψ ) is expected to vanish smoothly as these QCPs are approached from the ASM side of transition. (To oneloop order η Ψ = 0, so that a nontrivial fermionic anomalous dimension requires the computation of self-energy diagrams to two-loop order. We leave this exercise for future investigation.) On the other hand, if we subscribe to an appropriate order-parameter theory (also known as the Yukawa formalism), where the order-parameter or bosonic field is coupled with gapless fermions, a leading order calculation yields non-trivial anomalous dimensions for both fermionic and bosonic fields [71]. The residue of the quasiparticle pole can serve as the order parameter in the semimetallic phase. The interacting QCPs describe strongly coupled non-Fermi liquids where the notion of sharp quasiparticle excitations becomes moot.

VII. SUSCEPTIBILITY OF SOURCE TERMS
To pin the actual nature of the broken symmetry, we track the enhancement or suppression of two-point correlations amongst Hermitian fermion bilinears (order parameters). We do this by computing the flow equations for their corresponding source terms. The Hamiltonian including all of the source terms in the particle-hole or excitonic channel reads The physical interpretation of all fermion bilinears was summarized in Table II. The leading order perturbative corrections to the source terms can be derived after computing the Feynman diagrams shown in Fig. 10. The resulting RG flow equations are   Table III after accounting for 1/n corrections in Eqs. (30) and (31).
The results obtain by utilizing the large-n expansions for f1(n), f2(n) and f3(n) [Eq. (26)]. The inclusion of 1/n corrections does not change the number of fixed points (at least for large enough n ≥ 4), but eliminates the spin-charge separation. Thus the -expansion, augmented by 1/n-expansion allows a controlled route to access some strong coupling phenomena (such as the quantum phase transition between the anisotropic semimetal and a broken symmetry phase) in two spatial dimensions. The coefficients of 1/n appearing in the locations of the fixed points are extracted numerically. The last column shows the number of unstable directions (UDs) at a given fixed point in the∆ = 0 hyperplane as well as in the five dimensional coupling constant space (shown in parentheses).
We also track the RG flow of the source terms corresponding to various local (momentum independent or intra unit-cell) superconducting orders, presented in Table II. The Hamiltonian capturing all local pairings assumes the form The leading RG flow equations for all the source terms associated with local pairings are given by The quantities on the right hand side of Eqs. (34) and (36) are the scaling dimensions of the corresponding Abelian spin current along x η3 (σ1, σ2, σ3) τ1 Spin bond-density-wave η0 (σ1, σ2, σ3) τ2 1 + − 1 4 + 0.092 source terms. The scaling dimensions associated to all fermionic bilinear source terms at the three interacting QCPs (namely FP 2 , FP 3 and FP 4 ) are displayed in Table V. If the source term has scaling dimension y s , the anomalous dimension of the associated fermion bilinear y b is given by using the fact that the operator and coupling constant scaling dimensions add to the scaling dimension of spacetime [72]. The latter is equal to 2 + 1/n according to Sec. II A. Since y b controls the decay of correlations (in the x-direction) of the bilinear two-point function at the quantum critical point [72], a larger source dimension y s means stronger correlations (smaller y b , slower decay).

A. Scaling dimensions of fermion bilinears: spin and pseudospin symmetries
In Table V, the bilinears previously cataloged in Table II are identified with 8 × 8 matrices, obtained by combining physical spin (σ) and sublattice/orbital (τ ) with particle-hole (η) degrees of freedom. The matrix representations of the various fermionic bilinears in Ta In this basis the three generators of physical spin ( S) and pseudospin ( P S) [36][37][38]73] symmetry are given by Note that the pseudospin generators do not involve genuine spin degrees of freedom, and one of its generators, namely P S 3 , is the charge density operator. The remaining two generators of pseudospin symmetry are the real and imaginary components of the chiral singlet pairing 2 . See Appendix C for symmetries in the Nambu basis.
We now discuss the results. First, we consider the 1D (n → ∞) limit, as a check. Then Eq. (26) implies that we should set f 1 (n) = 1 and f 2 (n) = 0 in Eqs. (34) and (36). We see that the charge bond density wave/anisotropy, CDW, spin bond density wave, and antiferromagnet source fields ∆ s,t 2,3 receive an anomalous correction +U A /4, while the singlet and triplet SC terms ∆ s,t receive the anomalous correction −U A /4, when U N = W = X = 0. The remaining terms receive no corrections in this limit. This is consistent with the FIG . 10: (a) The bare vertex associated with a source term Ψ † M Ψ (particle-hole channel) or ΨM Ψ (particle-particle channel). M and N are four-dimensional Hermitian matrices [see Eqs. (33) and (35)]. Renormalization of the bare vertex (a) to the leading order in interaction couplings arises from diagrams (b) and (c). Note there is no mixing of bilinear operators since all transform differently under the symmetries, see Table II and Appendix C.
predictions of bosonization, summarized in Table VII and using the Luttinger parameter in Eq. (D4).
Next, note that at FP 2 , the scaling dimensions for CDW and s-wave pairing are largest and exactly equal, reflecting the underlying pseudospin SU(2) symmetry among these two orders, see Fig. 11(I). Since the bipartite lattice Hubbard Hamiltonian possesses pseudospin SU(2) symmetry at half filling [36][37][38], we expect that FP 2 can be accessed by tuning the strength of onsite attraction (see Sec. VIII). On the other hand, the scaling dimension for CDW order is largest amongst all possible orders at the interacting QCP FP 3 . Hence, it is natural to anticipate that repulsive NN interaction in the strained honeycomb model can induce the transition across this fixed point. By contrast, the AFM order possesses the largest scaling dimension at the interacting QCP FP 4 , indicating that the onset of AFM for strong repulsive Hubbard interaction is controlled by this fixed point. In the following section we will substantiate these observations from the phase diagram of an extended Hubbard model. The fact that a two-dimensional ASM supports only six mass matrices representing CDW (1), AFM (3) and s-wave pairing (2), where the quantities in the parentheses denotes the number of matrices required to define a specific order parameter [see Table V], can also be established from the Clifford algebra of real matrices, as shown in Appendix G.
The pseudospin symmetry between s-wave pairing and CDW order at FP 2 extends to three other bilinear "triads". Each triad contains three bilinears that transform in the pseudospin triplet representation, as explained in Fig. 11. At FP 2 and FP 4 , the bilinears within each triad share a common scaling dimension, indicating that pseudospin SU(2) symmetry is emergent at both of these QCPs, see Figs. 11(II), (III), and (IV). This suggests that FP 2,4 are realized via attractive or repulsive onsite Hubbard interactions in the strained honeycomb lattice model, as we confirm in Sec. VIII. By contrast, pseudospin symmetry is explicitly broken at FP 3 . The perturbative restoration of pseudospin symmetry at inter-  acting QCPs that can be accessed by tuning the Hubbard interaction serves as a good anchoring ground of the analysis [74].

B. Identification of broken-symmetry phases
The strategy for identifying the actual nature of the broken symmetry at strong coupling is the following. We simultaneously run the flow of the interaction coupling constants and of the source terms. When interactions are weak all of them flow back to zero under course graining, and none of the source terms diverges, indicating the stability of the ASM for weak enough interactions.
As we keep increasing the strength of interaction, beyond a critical strength at least one of the interaction coupling constants diverges. At the same time at least one of the source term diverges. The channel (say ∆ µ ) that diverges fastest determines the pattern of spontaneous symmetry breaking and the resulting phase is characterized with the order parameter ∆ µ = 0. Following this strategy we present various cuts of the phase diagram of the interacting ASM in Fig. 12. We follow the same strategy to arrive at the phase diagram of an extended Hubbard model, discussed in the next section.

C. 1D physics in 2D?
The 1/n corrections to the interaction RG flow in Eq. (31) proportional to f 2 (n) mix the spin (U N ) and charge (U A , W , X) sector interaction coupling constants. This suggests that spin-charge separation is completely destroyed for any finite n; the one-dimensional physics does not survive once quantum fluctuations in 2D are incorporated. On the other hand, the pseudospin-SU(2)symmetric fixed points FP 2 and FP 4 are reminiscent of the 1D spin-charge separated Luther-Emery liquid and SDW phases, respectively. Both the 1D phases and 2D fixed points sit on the precipice of long-range CDW/swave pairing or antiferromagnetic order, and both are expected to lack well-defined quasiparticles. It would therefore be extremely interesting to look for remnants of 1D physics at these fixed points upon incorporating higher order corrections (two loops and beyond).

VIII. EXTENDED HUBBARD MODEL
We now focus on the phase diagram of the extended Hubbard model in 1D and for the ASM in two spatial dimensions. We discuss these two cases separately.

A. One dimension
The extended Hubbard interaction in 1D is Here j is the 1D lattice site index. Since the single-band 1D system supports two Dirac points at momentum k x = ±k F , we can expand the lattice fermion operator in terms of the low energy Fourier modes according to where R σ and L σ are right-and left-mover operators for modes at ±k F , with spin projection σ =↑, ↓. At halffilling the extended Hubbard interaction from Eq. (39) leads to the continuum form where the interaction operators are defined in Eq. (D2). The initial or bare conditions (values of coupling constants at RG time l = 0) for the 1D extended Hubbard model are The standard 1D RG flow equations in Eq. (32) (with = 0) then determine the phase diagram of the one-dimensional extended Hubbard model, displayed in Fig. 4(a) [30]. The diagonal phase boundaries in Fig. 4(a) are determined by U = ±2V , based on the lowest order Kosterlitz-Thouless equations [Eq. (32) with = 0].

B. Two-dimensional ASM
Since the Dirac-semimetal to band-insulator quantum phase transition takes place when two Dirac points with opposite vorticity merge at a high symmetry point in the Brillouin zone, the ASM in two dimensions lacks any valley (chiral) degeneracy. Thus, we can simply identify the sublattice density operators with Taking the continuum limit of the extended Hubbard Hamiltonian from Eq. (8), we obtain the interaction Hamiltonian H int in the form of Eq. (16). We thereby determine the initial conditions for the extended Hubbard model in the two-dimensional ASM, In this notation X > 0 corresponds to repulsive interaction, while attractive interaction is realized for X < 0, where X = U, V . We follow the strategy described in Sec. VII B for various values of U and V (both attractive and repulsive interactions) to arrive at the phase diagram of the interacting ASM in the U -V plane. Our result is shown in Fig. 4(b). For strong enough onsite repulsion U > 0, the susceptibility in the AFM channel displays the leading divergence, as shown in Fig. 13(a), indicating the onset of collinear antiferromagnetic (Néel) order on the strained honeycomb lattice. On the other hand, for sufficiently strong attractive U < 0, the s-wave and CDW channels diverge in an exactly degenerate fashion, as shown in Fig. 13(b), suggesting simultaneous nucleation of these two orders. Thus the onset of long-range order for strong repulsive and attractive Hubbard interaction takes place through spontaneous breaking of SU(2) spin and pseudospin symmetries, respectively. However, the exact (microscopic) pseudospin SU(2) symmetry gets destroyed for  . The resulting phase diagram in the U -V plane is shown in Fig. 4(b). For the clarity of presentation we here do not present the RG flow of the remaining channels. We have verified that in the entire phase diagram only one of these three channels displays the leading divergence at strong coupling.
Finally, we discuss the role of various interacting QCPs, discussed in Sec. VI, on the phase diagram of the extended Hubbard model. The fixed point controlling the transition into the broken symmetry phase can be identified by scrutinizing the RG flow trajectories of the interaction coupling constants. The continuous transition between the ASM and AFM for repulsive U > 0 is controlled by the interacting QCP FP 4 , while the transition into a coexisting phase of s-wave pairing and CDW order (U < 0, V = 0) takes place through the interacting QCP FP 2 . The ASM-CDW continuous transition driven by NN repulsion is governed by the interacting QCP FP 3 . Notice that for V > 0 as we keep increasing the strength of attractive onsite interaction, such that the broken symmetry phase at strong coupling is always the CDW, the fixed point that controls the ASM-CDW continuous transition switches from FP 3 to FP 2 for an intermediate strength of negative Hubbard-U . At this point ASM-CDW phase boundary displays a cusp, as shown with a dot in Fig. 4(b), and the ASM-CDW transition at this point is controlled by the bicritical point FP 5 . Since such a bicritical point is accessed by holding one of its two unstable directions fixed, the ASM-CDW transition is always continuous.
In the phase diagram of the extended Hubbard model for the two-dimensional ASM [see Fig. 4(b)] we show a green dashed line for repulsive NN but attractive onsite interaction. We now discuss the physical significance of this line. The scaling dimensions for s-wave pairing and CDW order source terms at the interacting fixed points FP 2 , FP 5 , FP 3 respectively read as While the scaling dimensions for these two orders are exactly equal at FP 2 , at the two other fixed points the scaling dimension for the CDW is always bigger than that for s-wave SC. Consequently, we determine that CDW order is favored for strong attractive U < 0 and repulsive V > 0. Also note that, starting from V = 0 and attractive U < 0, if we increase V along the ASM-CDW phase boundary [see the red line in Fig. 4(b)], the transition is first controlled by FP 2 (for weak V > 0), then FP 5 (for a specific strength of V indicated by the green dot), and finally by FP 3 (stronger NN repulsion V > 0). The scaling dimension of the s-wave pairing source term decreases monotonically from FP 2 to FP 5 to FP 3 for any n ≥ 2. Therefore, in the entire regime of the CDW phase [the red shaded region in Fig. 4(b)], s-wave pairing can only be found as a short-range order, without any global long range order. Furthermore, the coherence length (ξ) of such local s-wave pairing decreases monotonically as we keep increasing V > 0 (keeping U < 0), and across the green dashed line ξ vanishes smoothly. Therefore, above the green dashed line in Fig. 4(b), local short-range s-wave pairing loses its support. It must be noted that there is no genuine phase transition across the green line, it only represents a crossover boundary.
Our RG analysis leads to a fairly complete phase diagram of the extended Hubbard model in the uniaxially strained honeycomb lattice, residing at or very close to the Dirac-semimetal to band-insulator QCP. We identified the role of each of the six fixed points (Table IV) in the phase diagram of a simple microscopic interacting model. With recent progress in controlling the strength of electronic interactions in optical honeycomb lattices for ultracold fermionic atoms [11,12] and the realization of an ultracold-atom Fermi-Hubbard antiferromagnet [75], as well as quantum Monte Carlo simulation with onsite or NN interaction [47][48][49][50][51][52][53][54], we believe that the phase diagram of the extended Hubbard model can be established (at least partially) in numerical simulation and/or experiments.
The existence of the three interacting critical points FP 2,3,4 can be argued to survive beyond perturbation theory. For strong repulsive Hubbard interactions on a bipartite lattice in 2D, SU(2) spin rotational symmetry is broken and the system enters into the AFM phase, while for strong attractive Hubbard interactions pseudospin SU(2) symmetry is lifted leading to simultaneous nucleation of s-wave pairing and CDW. By contrast, for strong repulsive nearest-neighbor interactions on a bipartite lattice, the system enters the CDW phase where discrete Z 2 sublattice symmetry is broken. The existence of these phases at sufficiently strong interaction can be anticipated for the strained honeycomb lattice model, which is free of frustration for all of the above orderings. However, due to the vanishing DoS in the ASM ( (E) ∼ √ E) all the broken-symmetry phases are realized only beyond a critical threshold of interaction, and thus through a quantum phase transition. We expect that the transitions are continuous in the strained honeycomb lattice extended Hubbard model, and thus controlled by QCPs. We note that quantum Monte Carlo simulation in the Dirac system (isotropic honeycomb lattice model) also indicate continuous transitions through various QCPs [47][48][49][50][51][52][53][54]. Thus, future numerical works can shed light on the non-perturbative nature of these transitions, which are accessed here from leading order and 1/n expansions.

IX. DISCUSSION AND CONCLUSION
In this paper we discuss the effect of generic shortrange interactions at the quantum critical point separating a two-dimensional topological Dirac semimetal and a symmetry-preserving band insulator. The critical excitations separating these two phases constitute an anisotropic semimetal. The quasiparticle spectra in the anisotropic semimetal display both linear and quadratic dependence on momenta along two mutually orthogonal directions. Consequently, the quantum critical regime in the noninteracting system displays peculiar power-law scaling of thermodynamic and transport quantities that are distinct from their counterparts in a Dirac semimetal or a band insulator, as shown in Table I. Due to the vanishing density of states ( (E) ∼ √ E) the critical point separating the Dirac semimetal from the band insulator remains stable against sufficiently weak short range electron-electron interactions. However, at strong interaction this direct transition either (i) becomes a fluctuation-driven first-order transition (see Sec. V A), or (ii) gets avoided by an intervening broken-symmetry phase (see Sec. VI). Using renormalization group (RG) analysis we identify charge density wave, antiferromagnet, and spin singlet s-wave pairing as the dominant channels of symmetry breaking for the interacting model containing generic short-range interactions that for example encompass the extended Hubbard model, see the phase diagrams in Fig. 4(b) and Fig. 12.
We demonstrate that (see Sec. VI) the RG can be controlled by expanding about one dimension, for which short-range interactions are marginal at tree level. We recover standard results for 1D systems such as the absence of corrections for a U(1) current-current (Luttinger) perturbation, as well as spin-charge separation and Kosterlitz-Thouless transitions for spin-1/2 electrons. From the leading order calculation we find that the correlation length exponent at all the interacting multicritical points, describing continuous quantum phase transition from anisotropic semimetal to various broken symmetry phases is ν = 2, which is strikingly different from the ones across a Dirac-semimetal to brokensymmetry phase transition for which ν ≈ 1, as predicted in field theoretic works [29, 56, 68-70, 74, 76, 77] and has also been established numerically [48-51, 53, 54].
Throughout we have neglected the long-range tail of the Coulomb interaction and focused only on its shortrange pieces. Recently it has been proposed that infinitesimally weak long-range Coulomb interaction can cause an instability of the critical excitations and drive the system toward the formation of an infrared-stable non-Fermi liquid phase [4,5]. In the future it will be interesting to study the interplay of non-Fermi liquid and broken-symmetry phases. However, in an optical honeycomb lattice (composed of neutral atoms) our study should provide a fairly complete picture of the phase diagram for the interacting uniaxially strained honeycomb lattice. Our results also apply to solid state compounds with screened Coulomb interactions (due e.g. to an external gate).
The RG framework developed here (combined and 1/n expansion) can be applied to many different systems where quasiparticle excitations possess anisotropy, such as a three-dimensional general Weyl semimetal constituted by monopoles and antimonopoles of strength N [55], with N = 1, 2, 3 in a crystalline environment [78][79][80]. In addition, our formalism (with minor modifications) can also be subscribed to address the effects of quenched disorder at the Weyl semimetal-band insulator quantum critical point and possible onset of a metallic phase at strong disorder through a continuous quantum phase transition across a multi-critical point [81]. Fur-thermore, our formalism can also shed light on the effects of electronic interaction in a two-dimensional Dirac semimetal. At the cost of the in-plane rotational symmetry (v x = v y , where v j is the Fermi velocity of massless Dirac fermions in the j direction), which can be achieved by applying a weak uniaxial strain in monolayer graphene (so that the system is sufficiently far from the Dirac-semimetal to band-insulator quantum critical point), we can generalize the proposed and 1/n expansion to address competing orders, emergent quantum critical phenomena, and the global phase diagram of a 2D Dirac system (with the caveat that n can only be odd integer). The interacting theory for a slightly anisotropic Dirac semimetal (with two valleys) is described by 4 (8) independent coupling constants for spinless (spinful) fermions. This exercise will allow us to investigate an intriguing competition among various ordered phases, such as anti-ferromagnet, valence-bondsolid, charge-density-wave, singlet and triplet superconductors as well as topological quantum anomalous/spin Hall insulators, in a controlled renormalization group framework [6]. We leave this problem for future work [82].
Finally, we point out that the dispersion of the ASM studied here closely resembles that of quasiparticles near hot-spots of the Fermi surface in a two-dimensional square lattice system [83]. Recently there has been a surge of theoretical works geared toward understanding the effects of strong forward scattering interactions and quantum critical phenomena within the hot-spot or patch model [84][85][86][87][88]. Prior works are often based on appropriate order parameter field (spin-density-wave or nematic, for example), coupled with gapless fermions, residing close to the hot-spots, through Yukawa coupling. Our approach can provide an alternative route to investigate the effects of strong electronic interactions on the Fermi surface, the role of competition among various incipient orderings, emergent symmetry near quantum critical points, etc. Very recently it has also been argued that a topologically protected, multi-flavored 2D anisotropic semimetal can also be realized in Sr 2 IrO 4 [89]. When electron-doped, Sr 2 IrO 4 is believed to support high-T c dwave superconductivity [90,91]. Therefore, one can (at least in principle) generalize the formalism outlined here to address the intriguing confluence of electronic correlations, exotic broken-symmetry phases, competing orders, quantum critical phenomena, the role of topology (such as in Weyl materials) in a wide variety of systems within a unified perturbatively controlled RG scheme.
Appendix A: Optical conductivity in the anisotropic semimetal Due to the distinct power-law dependence of the quasiparticle dispersion along the x-and y-directions, the critical excitations show distinct scaling of the dynamic (frequency-dependent ac) conductivity along these two directions (denoted by σ xx and σ yy , respectively). We use the Kubo formula to compute the dynamic conductivity. The expression for the polarization bubble is where e is the electronic charge, β is the inverse temperature, ω n , p m are fermionic Matsubara frequencies, J x = vτ 1 ,Ĵ y = 2bk y τ 2 , and N counts the flavor degeneracy. Thus for spinful fermions N = 2. The noninteracting Green's function reads as where µ is the chemical potential and we have used the Lehmann representation in the final expression. The spectral function A(k, ω) is given by A(k, ω) = π 1 + 2 a=1 τ ada δ (ω + µ − E k ) + π 1 − 2 a=1 τ ada δ (ω + µ + E k ) , (A3) where E k = v 2 k 2 x + b 2 k 4 y ,d a = d a / d 2 1 + d 2 2 for a = 1, 2 and d 1 = vk x , d 2 = bk 2 y . The conductivity is defined as σ jj = −Im[Π jj (ω)]/ω. The Drude conductivity (the contribution that comes only from the intra-band piece) evaluates to Eq.
The scaling of F 1 (x) and F 2 (x) is shown in Fig. 3(a). The result for the inter-band component of the optical conductivity is given by Eq. (3), where ω 0 = v 2 /b. Thus as the frequency ω → 0, the dynamic conductivity along the y-(x)-direction vanishes (diverges) as √ ω (1/ √ ω). Such distinct power-law dependence of the conductivity along different directions can be directly measured in experiment to pin the quasiparticle dispersion. Appendix B: Diamagnetic susceptibility Next we will discuss the effects of external magnetic fields on critical excitations, and in particular we focus on the diamagnetic susceptibility. The Landau level (LL) spectrum of the ASM, described by the Hamiltonian H(k, 0), is ±E n (B), where n is the LL index 13 and for n = 0, 1, 2, · · · E n (B) = A 2v 2 b 1/3 ω c n + 1 2 with A ≈ 1.173, ω c = 2eBb is the cyclotron frequency, and B is the strength of the external magnetic field [10]. The free energy at T = 0 is defined as Here C = e/h and CB counts the LL degeneracy. Note that Ω 0 is a divergent quantity, which, however, can be regularized by using the zeta regulator as follows [26,92,93] π s ζ (1 − s) = 2 1−s Γ(s)ζ(s) cos sπ 2 .
(B3) Upon using the above regularization scheme and the definition of the diamagnetic susceptibility χ = ∂ 2 Ω 0 /∂B 2 , we arrive at the expression for the two-dimensional ASM at T = 0 where Notice that the diamagnetic susceptibility diverges as B −1/3 as B → 0. The free energy at finite temperature reads as Ω(T ) = Ω 0 + Ω T , where 13 The LL index n in this Appendix should not be confused with the anisotropy parameter introduced in Eq. (6) Bilinear description PH and B = /(eB) is the magnetic length. After following the standard steps, highlighted above, we arrive at the following expression for the diamagnetic susceptibility at finite temperature χ(T ) = χ 0 + χ T , where χ T = χ 0 f (λ T h /l B ) and ≈ −(11.063) g(x) + 3 10 Here λ T h is the thermal de Broglie wavelength for the critical excitations. Our proposed scaling of the diamagnetic susceptibility at finite temperature is valid as long as l B < λ T h < ξ ∆ , with ξ ∆ ∼ 1/∆. Note f (x) is a universal function of dimensionless argument, and in Fig. 3 Again since the Nambu PH is automatic, Eq. (C8) is an equivalent definition of time-reversal symmetry. The full inventory of 28 Nambu-PH Hermitian fermion bilinears appears in Table V. The properties of these bilinears under SLS, T , x-reflection R π , and spin SU(2) symmetry transformations were summarized in Table II. To this we add the properties under microscopic PH symmetry [Eqs. (C2) and (C7c)], shown in Table VI. We list only particle-hole channel bilinears.
Appendix D: 1D limit and review of spin-charge separation We translate the anisotropic semimetal model into the standard notations for 1D physics. After a τ 2 -rotation that sends τ 1 → τ 3 , the field in Eq. (10) is decomposed into right-R and left-L mover components via Ψ T = [R, L], leading to For the case of spin-1/2 electrons, the chiral components R and L are each two-component spinors.
The four four-fermion interaction operators appearing in Eq. (17) are chosen to exploit spin-charge separation [30,31] in the 1D (n → ∞) limit. These are defined as The operator O N is an SU(2) 1 current-current perturbation; the remaining operators couple only to the U(1) electric charge sector [31].
In the 1D limit, we can bosonize the U(1) charge and spin SU(2) 1 sectors; up to irrelevant operators, these completely decouple [31]. The imaginary time bosonic action for the charge sector can be written as where φ (θ) denotes the polar (axial) field, and the charge velocity v c and Luttinger parameter K c are given by v c = v (1 + X) , K c = 1 − U A /2.
The U(1) current-current operator O A and U(1) stress tensor operator O X are completely absorbed into these  Scaling dimensions of all fermion bilinears in the 1D limit ( = 0, n → ∞), in the stable gapless phase consisting independent spin and charge Luttinger liquids, c.f. Fig. 4(a). The common dimensions for the AP, CDW, spin-BDW, and AFM bilinears are the sum of the spin (1/2) and charge (Kc/2) sector contributions; the same applies to the s-wave SC and triplet SC. The chiral SC1,2 bilinears reside entirely in the charge sector. We ignore logarithmic corrections to the spin sector due to nonzero, but marginally irrelevant UN < 0 [30].
parameters of the free boson, while the umklapp operator 1 2 (O W +Ō W ) becomes a sine-Gordon perturbation. As usual [30,31], the expressions for v c and K c in terms of the microscopic coupling strengths U A and X obtained via perturbative bosonization are valid only to first order. When both charge and spin sectors are gapless, the theory asymptotes to a conformal fixed point in the infrared. This is a product of spin and charge Luttinger liquids [c.f. Fig. 4(a)]. The SU(2) invariant fixed point has U N = W = 0, with U A < 0 (K c > 1). The exact scaling dimensions of all fermion bilinears in Tables II  and V at this fixed point are summarized in Table VII.

Appendix E: Table of integrals
To proceed with the perturbative RG analysis to the leading order in and 1/n, we need to evaluate the following integrals where ρ = ω 2 + v 2 k 2 x , and E Λ is the ultraviolet energy cutoff. Within the framework of the Wilsonian RG, we integrate out the fast Fourier modes within the shell E Λ e −l < ρ < E Λ , and subsequently integrate over 0 ≤ θ ≤ 2π and 0 ≤ x ≤ ∞. The above integrals are then given by (1 + x 2n ) 2 = a l π csc π 2n 2n 2 = a 1 n l + O 1 n 2 , (E3b) where a ≡ E Λ /(8π 2 vb ). Thus as n → ∞ only the contributions from I 1 and I 2 survive. We note that all integrals over x are convergent for any n ≥ 2.
operators, K j for j = 1, 2, 3, besides the standard identity matrix. The Casimir operators satisfy the quaternionic algebra where δ ij is the Kronecker delta function. From the Casimir operators we can define imaginary matrices E j = iK j for j = 1, 2, 3 and E j s satisfy the SU(2) algebra.
Notice that E j s are purely imaginary Hermitian matrices that commute with H ASM (k). From these three imaginary matrices we can define yet another set of three mutually anticommuting matrices M j = iE j RI 1 for j = 1, 2, 3 that together close a Cl(3) algebra and also anticommute with H ASM (k). Hence, there are all together six purely imaginary matrices (thus mass matrices) that anticommute with H ASM (k) and they can be arranged into two sets according to Together these two sets close a Cl(3) × Cl(3) algebra of mass matrices, as we announced in the main part of the paper. One set we can identify with the three components of Néel AFM order, while the other set is constituted by the CDW and two components (real and imaginary) of s-wave pairing. We note that Clifford algebra of mass matrices in a 2D ASM is identical to that for a three-dimensional line-node semimetal [95].