Bond-Nematic Dirac Semimetal

• Bond-nematic Dirac semimetal is defined by spontaneous C4 to C2 symmetry breaking that splits a quadratic band touching into two distinct Dirac cones with a linearly vanishing density of states.
• The phase is established nonperturbatively via CP-QMC and tensor-network methods, which reveal convergence with DMRG in energy, order parameters, and finite-size scaling.
• This gapless state contrasts with gapped site-nematic and quantum anomalous Hall phases, offering unique thermodynamic signatures and transport anisotropy.

Bond-nematic Dirac semimetal (BNDS) denotes an interaction-driven electronic phase in which lattice rotational symmetry is spontaneously reduced, typically $C_4 \to C_2$, through anisotropy of bond expectation values while translation symmetry is preserved and, in the canonical checkerboard-lattice realization, time-reversal symmetry remains intact. In the two-dimensional quadratic-band-touching problem, bond nematicity lifts a symmetry-protected quadratic band touching (QBT) and splits it into two symmetry-related Dirac cones, producing a gapless Dirac semimetal with linearly vanishing density of states at low energy rather than a gapped nematic insulator (Liu et al., 21 Jul 2025, Lu et al., 2023). Related three-dimensional Dirac-semimetal literature studies closely allied nematic bilinears that also break rotational symmetry by mixing symmetry-distinct components at a Dirac point, but there the order generally opens a gap; this contrast suggests using BNDS most strictly for the gapless QBT-derived phase and for material settings where bond-centered nematicity coexists with Dirac quasiparticles (Zhang et al., 2015, Butler et al., 2022).

1. Symmetry content and defining criteria

The defining symmetry pattern of the checkerboard-lattice BNDS is the spontaneous breaking of lattice fourfold rotation to twofold rotation, $C_4 \to C_2$, without breaking translation symmetry. In the nearest-neighbor $t_1$–$t_2$–$V$ model at half filling, the noninteracting QBT is protected jointly by time-reversal symmetry and $C_4$, with the rotation center at the midpoint between two sites of the same sublattice. Preserving these symmetries forbids the linear terms that would otherwise split the QBT. Bond-nematic order removes this protection by making next-nearest-neighbor hopping expectation values anisotropic along the two primitive directions while leaving the pattern translationally uniform (Liu et al., 21 Jul 2025).

The bond character of the order is expressed by unequal intra-sublattice correlations along $a_1=(1,0)$ and $a_2=(0,1)$,

$$\Delta_{\mathrm{bond}}=\frac{1}{2}\sum_\lambda \left[ \left|\langle c^\dagger_{r,\lambda} c_{r+a_1,\lambda}\rangle\right|-\left|\langle c^\dagger_{r,\lambda} c_{r+a_2,\lambda}\rangle\right| \right].$$

In the CP-QMC study of the nearest-neighbor-repulsion model, this bond anisotropy appears while on-site charge densities remain equal between the two sublattices, so the phase is nematic without sublattice charge imbalance. By contrast, the site-nematic insulator (SNI) has a nonzero sublattice density imbalance and is fully gapped, even though it breaks the same lattice symmetry and also preserves time-reversal symmetry (Liu et al., 21 Jul 2025).

A common confusion is to identify BNDS and SNI by different broken symmetries. In the broader $V_1$–$C_4 \to C_2$0 checkerboard study, both BNDS and SNI break exactly the same $C_4 \to C_2$1 symmetry and preserve time-reversal symmetry; their decisive distinction is spectral. BNDS is gapless, with Dirac cones inherited from the split QBT, whereas SNI is gapped because the Dirac nodes merge and annihilate when the site-nematic component becomes sufficiently strong (Lu et al., 2023).

2. Emergence from quadratic band touching

The canonical microscopic realization is a model of spinless fermions on a two-dimensional checkerboard lattice with two sites per unit cell and primitive vectors $C_4 \to C_2$2, $C_4 \to C_2$3. The Hamiltonian used in the CP-QMC study is

$C_4 \to C_2$4

$C_4 \to C_2$5

$C_4 \to C_2$6

with $C_4 \to C_2$7, $C_4 \to C_2$8, and $C_4 \to C_2$9 (Liu et al., 21 Jul 2025).

At half filling, the noninteracting band structure has a single QBT at $t_1$0. Near that point the dispersion is quadratic, $t_1$1, so the density of states at the Fermi level is finite. That finite density of states makes the QBT especially susceptible to interaction-driven instabilities. Bond-nematic order exploits this instability in a specific manner: because $t_1$2 protects the QBT, spontaneous bond anisotropy removes the protection and splits the single QBT into two symmetry-related Dirac cones along the edge of the first Brillouin zone. The density of states thereby changes from finite at the QBT to linearly vanishing in energy, characteristic of Dirac fermions (Liu et al., 21 Jul 2025).

The same mechanism is explicit in the more general mean-field parameterization used in the tensor-network study,

$t_1$3

with $t_1$4 and $t_1$5, which realizes a QBT at $t_1$6. In that framework, a bond-nematic term $t_1$7 splits the QBT into two Dirac cones but does not by itself open a gap, whereas a site-nematic term $t_1$8 can drive the Dirac cones together and gap them out (Lu et al., 2023).

This mechanism is interaction-induced rather than externally imposed. The resulting semimetal is therefore not a weakly perturbed band structure with accidental anisotropy, but a many-body phase selected by repulsion in a QBT system.

3. Nonperturbative establishment and the failure of site-independent Hartree–Fock

The phase diagram of the checkerboard QBT model remained disputed because determinant QMC becomes impractical precisely where the intermediate phase emerges: on torus geometries the average sign drops below $t_1$9 at $t_2$0 for $t_2$1. The constrained-path quantum Monte Carlo (CP-QMC) study addressed this regime by projecting the ground state in Slater-determinant space with auxiliary fields and imposing the overlap constraint $t_2$2, which removes the sign problem at the cost of a trial-state-dependent bias (Liu et al., 21 Jul 2025).

For observables commuting with $t_2$3, the calculation uses the mixed estimator

$t_2$4

while noncommuting operators are evaluated by backpropagation,

$t_2$5

with $t_2$6 sufficient for convergence in the reported calculations. The Hubbard–Stratonovich decomposition is taken in the charge ($t_2$7) channel because it yields a less severe sign problem than the $t_2$8 channel (Liu et al., 21 Jul 2025).

A central technical point is the improved trial-state optimization. The mixed single-particle Green’s function $t_2$9 is iteratively diagonalized and used to rebuild $V$0, but on torus geometries translational degeneracies can trap the self-consistency in local minima. The improved scheme therefore launches multiple self-consistency cycles from different site-independent mean-field solutions and retains the converged state with lowest energy. This refinement is what enabled the CP-QMC calculation to resolve the BNDS on a symmetry-preserving torus and to extrapolate to the thermodynamic limit (Liu et al., 21 Jul 2025).

Cross-validation with DMRG is quantitative. On a small torus $V$1, CP-QMC energies agree with DMRG to within $V$2 over $V$3, and phase boundaries as well as order-parameter jumps coincide. On cylinders such as $V$4, CP-QMC reproduces DMRG energies with $V$5 relative error and the same phase structure. By contrast, site-independent Hartree–Fock stabilizes at intermediate coupling a mixed phase that breaks both time-reversal symmetry and $V$6. The BNDS found by CP-QMC and DMRG instead preserves time-reversal symmetry. The discrepancy identifies BNDS as a fluctuation-driven phase beyond static mean-field treatment (Liu et al., 21 Jul 2025).

4. Phase structure, order parameters, and low-energy diagnostics

At half filling, the CP-QMC phase diagram on the square torus contains three interaction-induced phases: a quantum anomalous Hall (QAH) insulator at weak coupling, a BNDS at intermediate coupling, and an SNI at strong coupling. On the $V$7 torus with periodic boundary conditions, first-order transitions occur at $V$8 for QAH $V$9 BNDS and $C_4$0 for BNDS $C_4$1 SNI, with discontinuous jumps in $C_4$2, $C_4$3, and $C_4$4. DMRG bipartite entanglement entropy shows step-like jumps at the same couplings. For larger tori up to $C_4$5, the BNDS broadens and shifts toward weaker coupling, approximately $C_4$6, and remains finite after extrapolation to the thermodynamic limit (Liu et al., 21 Jul 2025).

The principal order parameters are:

• the QAH current-current correlator $C_4$7, whose square root is used as the QAH order parameter;
• the bond-nematic order $C_4$8, finite in BNDS;
• the site-nematic order

$C_4$9

finite in SNI.

For representative couplings, finite-size extrapolation gives $a_1=(1,0)$0 at $a_1=(1,0)$1 in BNDS, $a_1=(1,0)$2 at $a_1=(1,0)$3 in SNI, and $a_1=(1,0)$4 at $a_1=(1,0)$5 in QAH. On cylinders, CP-QMC again reproduces the DMRG sequence, with QAH $a_1=(1,0)$6 BNDS near $a_1=(1,0)$7 and BNDS $a_1=(1,0)$8 SNI at $a_1=(1,0)$9 (Liu et al., 21 Jul 2025).

Phase	Symmetry and spectrum	Principal diagnostic
QAH	breaks time-reversal; preserves translation and $a_2=(0,1)$0; fully gapped	finite $a_2=(0,1)$1
BNDS	breaks $a_2=(0,1)$2; preserves translation and time-reversal; gapless with two Dirac cones	finite $a_2=(0,1)$3
SNI	breaks the same $a_2=(0,1)$4; preserves time-reversal; fully gapped	finite $a_2=(0,1)$5

The tensor-network study extends the characterization to finite temperature and to a broader interaction space with nearest- and next-nearest-neighbor repulsions $a_2=(0,1)$6 and $a_2=(0,1)$7. There BNDS generically intervenes between SNI and either QAH or stripe order. For $a_2=(0,1)$8, representative ground-state points appear near $a_2=(0,1)$9 for SNI, $$\Delta_{\mathrm{bond}}=\frac{1}{2}\sum_\lambda \left[ \left|\langle c^\dagger_{r,\lambda} c_{r+a_1,\lambda}\rangle\right|-\left|\langle c^\dagger_{r,\lambda} c_{r+a_2,\lambda}\rangle\right| \right].$$0 for BNDS, $$\Delta_{\mathrm{bond}}=\frac{1}{2}\sum_\lambda \left[ \left|\langle c^\dagger_{r,\lambda} c_{r+a_1,\lambda}\rangle\right|-\left|\langle c^\dagger_{r,\lambda} c_{r+a_2,\lambda}\rangle\right| \right].$$1 for QAH, and $$\Delta_{\mathrm{bond}}=\frac{1}{2}\sum_\lambda \left[ \left|\langle c^\dagger_{r,\lambda} c_{r+a_1,\lambda}\rangle\right|-\left|\langle c^\dagger_{r,\lambda} c_{r+a_2,\lambda}\rangle\right| \right].$$2 for stripe order; for $$\Delta_{\mathrm{bond}}=\frac{1}{2}\sum_\lambda \left[ \left|\langle c^\dagger_{r,\lambda} c_{r+a_1,\lambda}\rangle\right|-\left|\langle c^\dagger_{r,\lambda} c_{r+a_2,\lambda}\rangle\right| \right].$$3, QAH disappears but BNDS remains as an intermediate phase between stripe and SNI. The transitions between QAH, BNDS, and SNI are first-order. A distinctive thermodynamic signature is $$\Delta_{\mathrm{bond}}=\frac{1}{2}\sum_\lambda \left[ \left|\langle c^\dagger_{r,\lambda} c_{r+a_1,\lambda}\rangle\right|-\left|\langle c^\dagger_{r,\lambda} c_{r+a_2,\lambda}\rangle\right| \right].$$4 at low temperature in BNDS, consistent with two-dimensional Dirac fermions, whereas SNI, QAH, and stripe phases show activated $$\Delta_{\mathrm{bond}}=\frac{1}{2}\sum_\lambda \left[ \left|\langle c^\dagger_{r,\lambda} c_{r+a_1,\lambda}\rangle\right|-\left|\langle c^\dagger_{r,\lambda} c_{r+a_2,\lambda}\rangle\right| \right].$$5. DMRG with twisted boundary conditions places the Dirac nodes in BNDS near $$\Delta_{\mathrm{bond}}=\frac{1}{2}\sum_\lambda \left[ \left|\langle c^\dagger_{r,\lambda} c_{r+a_1,\lambda}\rangle\right|-\left|\langle c^\dagger_{r,\lambda} c_{r+a_2,\lambda}\rangle\right| \right].$$6, with a minute flux-dependent shift $$\Delta_{\mathrm{bond}}=\frac{1}{2}\sum_\lambda \left[ \left|\langle c^\dagger_{r,\lambda} c_{r+a_1,\lambda}\rangle\right|-\left|\langle c^\dagger_{r,\lambda} c_{r+a_2,\lambda}\rangle\right| \right].$$7 where the single-particle gap is minimized and the fitted central charge peaks; the central charge fit on YC6 gives $$\Delta_{\mathrm{bond}}=\frac{1}{2}\sum_\lambda \left[ \left|\langle c^\dagger_{r,\lambda} c_{r+a_1,\lambda}\rangle\right|-\left|\langle c^\dagger_{r,\lambda} c_{r+a_2,\lambda}\rangle\right| \right].$$8, though not yet fully converged (Lu et al., 2023).

5. Relation to three-dimensional Dirac semimetals and nematic mass generation

In three-dimensional Dirac semimetals, rotationally protected Dirac nodes admit a different but closely related form of nematic instability. A minimal continuum description near the two Dirac points $$\Delta_{\mathrm{bond}}=\frac{1}{2}\sum_\lambda \left[ \left|\langle c^\dagger_{r,\lambda} c_{r+a_1,\lambda}\rangle\right|-\left|\langle c^\dagger_{r,\lambda} c_{r+a_2,\lambda}\rangle\right| \right].$$9 in Na$V_1$0Bi is

$V_1$1

with $V_1$2. The Dirac points are protected by crystalline rotation because the two components at the node transform differently under the protecting $V_1$3. A nematic bilinear breaks that rotation and mixes the symmetry-distinct states at fixed momentum:

$V_1$4

Because $V_1$5 anticommutes with $V_1$6, it opens a gap rather than producing a gapless semimetal. Mean-field energetics favor $V_1$7, and the resulting nematic state breaks rotation and typically time-reversal as well (Zhang et al., 2015).

This is the essential contrast with the checkerboard BNDS: in the two-dimensional QBT problem, bond nematicity splits a quadratic node into two Dirac cones and leaves the system gapless; in the three-dimensional Na$V_1$8Bi-type continuum theory, the analogous rotation-breaking order is a mass term for an already linear Dirac node. The same paper also emphasizes why this instability is specific to Dirac, rather than Weyl, semimetals: a single Weyl node cannot be gapped by any symmetry-preserving bilinear at the same momentum, whereas a Dirac node consists of two symmetry-distinct components at the same momentum that can hybridize once the protecting rotation is broken (Zhang et al., 2015).

A related strong-field problem appears in the quantum limit of Dirac semimetals. Under $V_1$9, the zeroth Landau bands of a Na$C_4 \to C_2$00Bi-type model become highly degenerate, and a translationally invariant excitonic condensate

$C_4 \to C_2$01

with $C_4 \to C_2$02, yields a polarized nematic phase. Its mean-field zeroth-Landau-band spectrum is

$C_4 \to C_2$03

so the instability again gaps the crossing. The field-induced nematic preserves translation, breaks rotational symmetry, and is interpreted as bond-nematic in the effective lattice sense because the excitonic self-energy renormalizes the in-plane kinetic terms anisotropically, producing unequal effective bond or hopping expectation values along orthogonal directions (Song et al., 2017).

6. Materials, spectroscopic signatures, and scope of the concept

Quadratic band touchings occur in Bernal-stacked bilayer graphene, kagome metals, HgTe, pyrochlore iridates, and proposed CrCl$C_4 \to C_2$04(pyrazine)$C_4 \to C_2$05 monolayers. In such systems, if interactions or strain break $C_4 \to C_2$06 while preserving time-reversal symmetry, a BNDS-like state may emerge. The momentum-space signature is the splitting of a single QBT at the zone corner into two Dirac cones along the Brillouin-zone edge, accompanied by anisotropic Dirac velocities and the absence of a gap. Transport should show anisotropic conductivity with preserved time-reversal symmetry, and low-temperature thermodynamics should follow Dirac power laws rather than the finite-density-of-states behavior of a QBT or the activated behavior of an insulator. These diagnostics distinguish BNDS from QAH, which has a full gap, chiral edge modes, and a quantized Hall response, and from SNI, which has a full bulk gap and sublattice charge imbalance (Liu et al., 21 Jul 2025, Lu et al., 2023).

BaNiS$C_4 \to C_2$07 provides a materially different but conceptually connected realization of bond-centered nematicity coexisting with Dirac electrons. STM detects electronic nematicity as two $C_4 \to C_2$08-symmetry striped patterns in the local density of states at approximately $C_4 \to C_2$09 and $C_4 \to C_2$10 meV, and quasiparticle interference identifies a splitting of the Brillouin-zone-boundary electron pockets by $C_4 \to C_2$11 meV. The observed anisotropy is modeled by a uniform $C_4 \to C_2$12 bond order with $C_4 \to C_2$13-form factor

$C_4 \to C_2$14

with $C_4 \to C_2$15 meV. The Dirac cones lie near the nodes of $C_4 \to C_2$16 along $C_4 \to C_2$17–$C_4 \to C_2$18, so they remain almost unaffected; STM and QPI therefore show the coexistence of bond nematicity and Dirac quasiparticles rather than the QBT-splitting mechanism characteristic of checkerboard BNDS. Quantitatively, the surface-projected Dirac node lies at $C_4 \to C_2$19 meV, with velocities $C_4 \to C_2$20 m s$C_4 \to C_2$21 and $C_4 \to C_2$22 m s$C_4 \to C_2$23, while the observed QPI dispersions are renormalized by a factor $C_4 \to C_2$24 relative to correlation-corrected DFT (Butler et al., 2022).

Two broader lessons follow. First, bond-centered nematicity need not imply an insulating state; in the checkerboard QBT problem it stabilizes a symmetry-broken but gapless Dirac semimetal, and in BaNiS$C_4 \to C_2$25 it coexists with nearly unaffected Dirac cones. Second, the phrase “nematic Dirac phase” covers distinct regimes: a gapless Dirac semimetal generated by splitting a QBT, a gapped topological nematic phase generated by adding a nematic mass to a 3D Dirac node, and a field-induced excitonic nematic in the zeroth-Landau-level sector. The precise physical content is therefore fixed by whether rotational symmetry breaking splits, preserves, or gaps the relevant Dirac degeneracy (Zhang et al., 2015, Song et al., 2017).