---
abstract: 'We study interaction-induced broken symmetry phases that can arise in metallic or semimetallic band structures with two nested Weyl or Dirac loops. The odered phases can be of the charge or (pseudo)spin density wave type, or superconductivity from interloop pairing. A general analysis for two types of Weyl loops is given, according to whether a local reflection symmetry in momentum space exists or not, for Hamiltonians having a global PT symmetry. The resulting density-wave phases always have lower total energy, and can be metallic, insulating, or semimetallic (with nodal loops), depending on both the reflection symmetry of the loops and the symmetry transformation that maps one loop onto the other. We extend this study to nested $\mathbb{Z}_2$ nodal lines, for which the ordered phases include also nodal point and nodal chain semimetals, and to spinful Dirac nodal lines. Superconductivity from interloop pairing can be fully gapped only if the initial double loop system is semimetallic.'
author:
- 'Miguel A. N. Araújo$^{1,2}$, and Linhu Li$^3$'
title: 'Broken-symmetry phases of interacting nested Weyl and Dirac loops '
---

Introduction
============

Since the discovery of topological insulators, band structures of fermionic systems with non-trivial momentum-space topology have received much attention in modern condensed matter physics. Their low energy description involves Dirac-like band dispersions, which in some cases imply gapless band structures characterized by the presence of nodal points or lines. Among these, are the nodal line semimetals (NLSMs) [@mullen; @hklattice]. A NLSM has valence and conduction bands touching along one-dimensional (1D) lines in the three-dimensional (3D) momentum space, and feature two-dimensional (2D) “drumhead" surface states surrounded by the nodal lines [@Burkov; @Yang2014; @Chen2015; @Weng2015; @Zhang2016]. Contrary to the well-studied topological insulating phases and nodal point semimetals, the 1D nodal lines of NLSMs provide rich topological structures such as links and knots [@Links_knots1; @Links_knots2; @Links_knots3; @Links_knots4; @Links_knots5], which cannot be described unambiguously by a single sign (e.g. the $Z_2$ number) or a integer (e.g. the Chern number) [@multi_links_Li]. On the other hand, a variety of gapped and gapless topological phases have been predicted in NLSMs (while possibly breaking certain symmetries). For example, a spin-orbit interaction can induce 3D Dirac semimetals from a NLSM [@SOC1; @SOC2], and periodical driving such as linear or circular polarized light may induce different types of nodal points [@light_driven1; @light_driven2; @light_driven3; @light_driven4; @light_driven5; @light_driven6]. By introducing various types of extra gapped terms, a NLSM can also be driven into several different types of topological insulators, including the recently discovered high-order topological insulators [@TI_NL1; @TI_NL2].

Spontaneous symmetry breaking from interactions in three dimensional systems with Weyl/Dirac nodal points or lines have also been addressed. For single nodal loop (NL) systems, superconducting and charge (or spin) density wave instabilities have been investigated using renormalization of fermionic interactions[@Nandkishore2], including also the mean-field description of the ordered phases[@Nandkishore1; @roy; @ryu].

Certain symmetries, such as spatial inversion or time-reversal, imply that Weyl nodes must occur at an even number of Brillouin Zone (BZ) points. Charge and spin density waves, as well as superconducting phases, which arise from nested spherical Fermi surfaces (FSs) in doped (or uncompensated) Weyl/Dirac points have been discussed[@wangye]. Weyl or Dirac NLs, on the other hand, do not necessarily have to exist in pairs. Although two-loop semimetals have not yet been found in nature, pairs of linked NLs (or Hopf-link structures) have been theoretically proposed[@Links_knots3; @Links_knots4; @Links_knots6]. Furthermore, a class of NLs protected by a combination of inversion and time-reversal has recently been discussed[@Z2_loop], which carry $\mathbb{Z}_2$ monopole charges, and must therefore be created or annihilated in pairs.

This has motivated us to address the spontaneous symmetry breaking from short range interactions in two-loop band structures, when the NLs are related through a nesting vector ${{\bf Q}}$ in the BZ. We describe the density wave and superconducting phases, which can be metallic, semimetallic (with double NLs) or fully gapped. A systematic analysis for two-band NL models is given, where the NLs can either satisfy a local reflection symmetry in the loop plane or not. If a global PT symmetry exists, then a symmetry operation can relate the two NLs. These properties combined determine the nature of the ordered phases. We also study specific four-band models that have appeared in the recent literature, such as the $\mathbb{Z}_2$ NLs, and NLs arising from perturbed Dirac points. Superconducting phases arising from pairing of fermions in different loops are also considered for all cases of singlet and triplet gap functions in loop space as well as (pseuso)spin space. But we have restricted our search to order parameters with time-reversal symmetry (TRS) and fully gapped phases, because the latter are expected to be more stable. The possibility of gap functions with a winding number, which break TRS, is not addressed here.

The structure of the paper is as follows. In Section \[modelsec\] we introduce the local $k\cdot p$ Hamiltonians for two-band Weyl NLs, which can either satisfy a local reflection symmetry in the NL plane, or not. The Hubbard interaction and the density wave order parameters associated with the NL nesting vector are also introduced. The density wave phases are described in Section \[densec\], and Section \[examples\] is devoted to example models and to a four-band system that was not included in the general analysis of the previous Sections: the nested $\mathbb{Z}_2$ NLs. In Section \[Diracsec\], we study spin degenerate Dirac NLs and also NLs arising from perturbed Dirac points. The superconducting pairing between nested NLs is studied in Section \[supersec\]. The analysis is focused on interloop pairing and time-reversal symmetric order parameters. In Section \[concsec\] we present our conclusions.

Model {#modelsec}
=====

We consider spinless fermions and let $\boldsymbol\tau$ denote the Pauli matrices acting on the pseudo-spin (orbital) degree of freedom. We assume that the band structure has two degenerate loops. If the Hamiltonian has PT symmetry, both loops involve the same Pauli matrices, $\tau_a\,,\tau_b$, so that each one can be locally described by $k\cdot p$ Hamiltonians:

$$\begin{aligned}
H_0({{\bf k}})&=& v_1\left( p_{\parallel} - p_o \right) \tau_a + v_2 p_\perp \tau_b \,, \\
H_0({{\bf k}}+\ {{\bf Q}})&=& g_1v_1\left( p_{\parallel} - p_o \right) \tau_a +g_2 v_2 p_\perp \tau_b\,,\end{aligned}$$

\[loop1\]

Here, and throughout the paper, $\bf p=\hbar{{\bf k}}$ and the subscripts $\perp$($\parallel$) refer to the components perpendicular (parallel) to the loop plane, and $g_{1(2)}$ are + or - signs. The loops are nested by the vector ${{{\bf Q}}}$. We shall refer to the NLs in Eq. (\[loop1\]) as “model-1” loops. Such NLs are protected by a local reflection[@bian] in the loop plane, ${\cal R}=(p_\perp\rightarrow -p_\perp)\otimes \tau_a$. Such a NL can be topologically characterized by a $\pi$ Berry phase along a trajectory enclosing the NL [@Burkov; @Zhang2016]. At zero chemical potential, the system is a semimetal and the FS consists of the two nested NLs. We shall also take non-degeneracy into account by considering an energy offset $\delta$ between the loops and make the replacement $H_0({{\bf k}})\rightarrow H_0({{\bf k}}) -\delta$, $H_0({{\bf k}}+\ {{\bf Q}})\rightarrow H_0({{\bf k}}+\ {{\bf Q}}) +\delta$. For positive $\delta$ and zero chemical potential, the FSs are torus-shaped, the one from $H_0({{\bf k}}+\ {{\bf Q}}) $ is in the lower (hole) band, while the FS from $H_0({{\bf k}})$ is in the upper (electron) band.

However, NLs are not necessarily protected by reflection symmetry. Here we also consider a more general model of nested NLs without reflection symmetry. The Hamiltonian reads

$$\begin{aligned}
H_0({{\bf k}})&=& \left[v_{1\parallel}\left( p_{\parallel} - p_o \right) + v_{1\perp} p_\perp \right] \tau_a + v_2 p_\perp \tau_b\,,
 \\
H_0({{\bf k}}+\ {{\bf Q}})&=& \left[g_1'v_{1\parallel}\left( p_{\parallel} - p_o\right)+  g_2'v_{1\perp}p_\perp \right]\tau_a +
g_2 v_2 p_\perp \tau_b\,,
\nonumber\\\end{aligned}$$

\[loop2\]

where $g_{1(2)}'$ are + or - signs. We refer to these as “model-2” loops. The extra $p_{\parallel}$ in the $\tau_a$ term changes the pseudospin texture near a NL, but does not affect the topological properties associated with the Berry phase. Examples of both types of NLs will be given in Section \[modelsec\]. The above two types of loops respond differently to the interactions, as shown in the following sections.

Normaly, one should expect that a perturbation arises that will lift the degeneracy between nested FSs. The perturbation may result from interactions and, in a normal system, usually takes the form of some charge or spin wave with the wave vector ${{{\bf Q}}}$. Also, superconducting pairing between fermions in different NLs will be considered.

In the rest of the paper we shall set to unity the velocity prefactors in the Hamiltonians (\[loop1\]) and (\[loop2\]), as they are not really necessary for the analysis that follows.

Density wave phases {#densec}
===================

Interaction and mean-field theory
---------------------------------

We introduce a Hubbard interaction, $$\begin{aligned}
\hat U &=& U\sum_{{\bf r}}\hat n_1({{\bf r}})\hat n_2({{\bf r}})\,,
\label{localU}\end{aligned}$$ where the indices $1,2$ refer to the orbital degree of freedom. Doing a mean field theory decoupling, the interaction takes the form $$\begin{aligned}
\hat U_{MF} = U \sum_{{\bf r}}\left[
 \langle n_1({{\bf r}}) \rangle \hat n_2({{\bf r}})  +  \hat n_1({{\bf r}})\langle n_2({{\bf r}})\rangle
- \langle n_1({{\bf r}}) \rangle \langle n_2({{\bf r}})\rangle  \right]\,.
\label{mfloop}\end{aligned}$$ A pseudospin density wave (PSDW) phase with the same nesting wave vector ${{{\bf Q}}}$ is characterized by $$\begin{aligned}
  \langle n_j({{\bf r}})\rangle  = \frac 1 2 n + \bar m (-1)^j\cos({{{\bf Q}}}\cdot{{\bf r}})\,,\end{aligned}$$ where $\bar m$ is the amplitude and $j=1,2$. Although this type of ordering describes an imbalance in orbital occupation, it is not a charge density wave (CDW) because the charge at site ${{\bf r}}$ is spatially constant, $n$. Omitting the factor $(-1)^j$, then a true CDW is obtained. We introduce the annihilation operator $\hat\psi_j({{\bf r}}) $, at point ${{\bf r}}$, with pseudospin index $j$. Then, Eq. (\[mfloop\]) can be rewritten as: $$\begin{aligned}
&&\hat U_{eff} =  \frac{Un}{2} \sum_{{\bf r}}\left( \hat\psi_1^\dagger({{\bf r}})  \hat\psi_2^\dagger({{\bf r}})    \right) \tau_0
 \left( \begin{array}{c} \hat\psi_1({{\bf r}}) \\ \hat\psi_2({{\bf r}})  \end{array}  \right)\nonumber\\ &&+\
U \bar m \sum_{{\bf r}}\cos({\boldsymbol Q}\cdot{{\bf r}}) \left( \hat\psi_1^\dagger({{\bf r}})  \hat\psi_2^\dagger({{\bf r}})    \right) \tau_3
\left( \begin{array}{c} \hat\psi_1({{\bf r}}) \\ \hat\psi_2({{\bf r}})  \end{array}  \right)\nonumber\\
&=&
 \frac{Un}{2} \sum_{{\bf k}}\left( \hat c_{{{\bf k}},1}^\dagger \hat c_{{{\bf k}},2}^\dagger    \right) \tau_0
  \left( \begin{array}{c} \hat c_{{{\bf k}},1}\\ \hat c_{{{\bf k}},2}  \end{array}  \right)
  \nonumber\\ 
&+&
\frac{U \bar m}{2} \sum_{{\bf k}}\left[ 
 \left( \hat c_{{{\bf k}}+{{\bf Q}},1}^\dagger  \hat c_{{{\bf k}}+{{\bf Q}},2}^\dagger   \right) \tau_3
  \left( \begin{array}{c} \hat c_{{{\bf k}},1}\\ \hat c_{{{\bf k}},2} \end{array}  \right) + 
   \left( \hat c_{{{\bf k}},1}^\dagger  \hat c_{{{\bf k}},2}^\dagger   \right) \tau_3
  \left( \begin{array}{c} \hat c_{{{\bf k}}+{{\bf Q}},1}\\ \hat c_{{{\bf k}}+{{\bf Q}},2} \end{array}  \right)\right]\nonumber\\
  \label{UEff}\end{aligned}$$ Replacing $\tau_3\rightarrow\tau_0$ in Eq. (\[UEff\]), we can describe a true CDW.

We write the effective Hamiltonian $H_{eff}({{\bf k}})$ matrix in operator basis $( \hat c_{{{\bf k}},1}, \hat c_{{{\bf k}},2} \hat c_{{{\bf k}}+{{\bf Q}},1}, \hat c_{{{\bf k}}+{{\bf Q}},2})
\equiv( \boldsymbol c_{{{\bf k}}} \ \boldsymbol c_{{{\bf k}}+{{\bf Q}}} )$ and introduce a factor $\frac 1 2$ to avoid double counting of momenta in the BZ: $$\begin{aligned}
\hat H_{eff}&=& 
\frac 1 2 
\sum_{{\bf k}}\left(  \boldsymbol c^\dagger_{{{\bf k}}}\ \boldsymbol c^\dagger_{{{\bf k}}+{{\bf Q}}}  \right)
H_{eff}({{\bf k}}) 
\left( \begin{array}{c}   \boldsymbol c_{{{\bf k}}} \\ \boldsymbol c_{{{\bf k}}+{{\bf Q}}}   \end{array}\right)\,,
\\
H_{eff}({{\bf k}}) &=&  \left(  \begin{array}{cc}
H_0({{\bf k}}) & U\bar m \tau_\alpha\\
 U\bar m \tau_\alpha & H_0({{\bf k}}+{{\bf Q}})  
\end{array}\right)    
\,,
\label{martizHloops} \end{aligned}$$ where $\alpha=0$ for CDW, or $\alpha=3$ for PSDW. The mean field equations for this Hamiltonian are derived in Appendix \[APMF\].

The effective Hamiltonian Eq. (\[martizHloops\]) is by no means restricted to the case of a local interaction as in Eq. (\[localU\]). A nearest-neighbor interaction, for instance, would produce an effective Hamiltonian of the same form, but where the bare interaction parameter would be multiplied by a ${{\bf Q}}$-dependent form factor which could still be denoted by “$U$”. The nesting property of the Fermi surface leads to a divergence of the susceptibilities in momentum space at the nesting wavevector[@gruner]. This always leads to a density wave with momentum $\bf Q$, described by the mean field couplings $<\hat c_{{{\bf k}}}^\dagger\hat c_{{{\bf k}}+{{\bf Q}}}>$, and the relevant interaction “$U$” would be the Fourier component of the interaction for the nesting wavevector.

model-1 loops {#model1_general}
-------------

We now introduce the Pauli matrices $t_\mu$ operating in loop space $({{\bf k}},{{\bf k}}+{{\bf Q}})$. For type-1 models the unperturbed double loop Hamiltonian has the form $$\begin{aligned}
H_{eff}^0({{\bf k}}) = \left( p_{\parallel} - p_o \right) t_i\tau_a+  p_\perp t_j\tau_b -\delta\ t_3
\label{H0eff}\end{aligned}$$ where $i,j$ can only take values 0 or 3. The effective Hamiltonian (\[martizHloops\]) can be written as $$H_{eff}({{\bf k}})= H_{eff}^0({{\bf k}}) + U\bar m\ t_1 \tau_\alpha
\label{umgap}$$ Supose that $\delta=0$, that is, degenerate loops at perfect compensation. It is clear that if the perturbing term $U\bar m t_1 \tau_\alpha$ anticommutes with only one term of $H_{eff}^0({{\bf k}}) $ then the resulting system still is a double NL semimetal. On the other hand, if $U\bar mt_1\tau_\alpha$ anticommutes with $H_{eff}^0({{\bf k}})$ then the resulting system is a gapped insulator, and if $\left[U\bar mt_1 \tau_\alpha, H_{eff}^0({{\bf k}})\right]=0$ then the original loops are shifted and a metallic phase arises with torus-shaped FSs, one of them hole-like, and the other electron-like.

Next we establish a criterion based on how a unitary transformation maps one loop onto the other. If the Hamiltonian has PT symmetry, one can always find a rotation through a Pauli matrix $\tau_\beta$ that maps one model-1 loop at ${{\bf k}}$ into the other at ${{\bf k}}+{{\bf Q}}$: $$\begin{aligned}
 H_0({{\bf k}})  =   \tau_{\beta} H_0({{\bf k}}+\ {{\bf Q}}) \tau_{\beta}   \,.
 \label{rotj}
 \end{aligned}$$ It is then convenient to apply a unitary transformation to the effective Hamiltonian in Eq. (\[martizHloops\]) according to $$\begin{aligned}
A H_{eff}({{\bf k}})  A^\dagger &=& \left( \begin{array}{cc}  
H_0({{\bf k}}) -\delta & U\bar m \tau_\alpha\tau_{\beta}\\ U \bar m \tau_{\beta}\tau_\alpha &  H_0({{\bf k}}) +\delta
 \end{array}   \right) \label{rotatedloops}\\
 A&=&\left( \begin{array}{cc}  1 & 0 \\ 0 &\tau_{\beta}  \end{array}   \right)\,.\label{Amatrix}
 \end{aligned}$$

The energy spectrum can be obtained by performing appropriate rotations on the matrix (\[rotatedloops\]), as shown in Appendix \[appa\]. We list all the four possibilities as follows.

1.  If $\tau_\alpha \tau_{\beta}=1$, the spectrum reads: $$\begin{aligned}
    E=\pm \sqrt{ \left( p_{\parallel} - p_o \right)^2+p_\perp^2}
    \pm   \sqrt{U^2\bar m^2 + \delta^2}
     \,,\label{CDWk0}\end{aligned}$$ (uncorrelated $\pm$ signs). In this case, $H_0({{\bf k}})=\tau_\alpha H_0({{\bf k}}+{{\bf Q}})\tau_\alpha$. The density wave produces a “level repulsion” effect by introducing (increasing) an energy splitting between the degenerate (non-degenerate) NLs. The density wave phase has two toroidal FSs, one hole-like and one electron-like.

2.  In the case where $\tau_\alpha \tau_{\beta}\propto\tau_a$, the energy spectrum is: $$\begin{aligned}
    E^2 &=& \left( p_{\parallel} - p_o \right)^2 + p_\perp^2+ U^2\bar m^2 +\delta^2 \nonumber\\
    &\pm& 2\sqrt{\left( p_{\parallel} - p_o \right)^2(U^2\bar m^2 + \delta^2)+
     p_\perp^2\delta^2}\,,\label{CDWka}\end{aligned}$$ which yields two NLs at $p_\perp=0$, $p_{\parallel}=p_o\pm\sqrt{U^2\bar m^2 + \delta^2}$. As $p_\parallel$ can only take a positive value, the loop with the minus sign will shrink into a point and vanish when $\sqrt{U^2\bar m^2 + \delta^2}$ becomes larger then $p_o$.

3.  For $\tau_\alpha \tau_{\beta}\propto\tau_b$ we get: $$\begin{aligned}
    E^2 &=& \left( p_{\parallel} - p_o \right)^2 + p_\perp^2+ U^2\bar m^2 +\delta^2 \nonumber\\
    &\pm& 2\sqrt{  p_\perp^2(U^2\bar m^2 + \delta^2)+
    \left( p_{\parallel} - p_o \right)^2\delta^2}\,.\label{CDWkb}\end{aligned}$$ This spectrum also gives two NLs, given by $p_\parallel = p_o$, $p_\perp=\pm \sqrt{U^2\bar m^2 + \delta^2}$. Unlike the previous case, these two loops only move along the $p_\perp$ direction when tuning $U$ or $\delta$, while their radii remain unchanged.

4.  For the case $\tau_\alpha \tau_{\beta}\propto\tau_{c(\neq a,b)}$ we have $$\begin{aligned}
    E^2 &=& U^2\bar m^2 + 
     \left[ \sqrt{ \left( p_{\parallel} - p_o \right)^2+p_\perp^2} \pm \delta \right]^2 \,,\label{CDWkc}\end{aligned}$$ which is then fully gapped. This is the case where $H_0({{\bf k}})=-\tau_\alpha H_0({{\bf k}}+{{\bf Q}})\tau_\alpha$.

At half filling (zero chemical potential), all the above energy dispersions lead to a lowering of the energy for $U\bar m\neq 0$, so, the density wave phase is energetically favorable. In a single NL, the density of states vanishes linearly at the chemical potential, so the broken symmetry phase appears only for $U$ above a finite critical value, $U_{cr}$. In any case if the phase transition is second order, then $U\rightarrow U_{cr}^+ \Rightarrow\ \bar m\rightarrow 0$. The PSDW cases occur for $U>U_{cr}>0$ but the CDW ordering requires negative $U<U_{cr}<0$, hence an attractive interaction (see Appendices \[APMF\] and \[MFloop\]).

From Eq. (\[rotj\]) we see that $\tau_\alpha H_0({{\bf k}}) \tau_\alpha =  
 \tau_\alpha\tau_{\beta} H_0({{\bf k}}+\ {{\bf Q}}) \tau_{\beta}\tau_\alpha$ and therefore, the density-wave phase is a nodal line semi-metal if $$H_0({{\bf k}})  \neq  \pm \tau_\alpha H_0({{\bf k}}+\ {{\bf Q}})\tau_\alpha \,,\label{critCDW}$$ where $\alpha=0$ for CDW, or $\alpha=3$ for PSDW.

For model-1 loops we can make the following observations regarding symmetry. The Hamiltonian (\[rotatedloops\]) for $\delta=0$ is chiral as it anti-commutes with the operator $\tau_{c\neq a,b}$, and contains the degenerate NLs. This chiral symmetry can be broken by a term of the form: [*(i)*]{} $t_\mu\tau_c$ which may fully gap the spectrum, or yield a semimetal, depending on its exact form; [*(ii)*]{} $t_\mu\tau_0$ which shifts the original loops and leads to a metallic spectrum. On the other hand, a term of the form $t_j\tau_a$ or $t_j\tau_b$ ($j=1,2,3$), preserves the chiral symmetry and yields the NL semimetal even if $\delta t_3 \neq 0$ is already present, as shown in Appendix \[appa\].

model-2 loops {#model2_general}
-------------

As we shall see, the criteria (\[critCDW\]) do not always apply to model-2 loops. The type-2 loop Hamiltonian at $\mathbf{k}$ is (omitting velocity prefactors): $$\begin{aligned}
H_0({{\bf k}})&=& \left( p_{\parallel} - p_o  + p_\perp\right)\tau_a + p_\perp \tau_b\,,\end{aligned}$$ and the loop at ${{\bf k}}+ {{\bf Q}}$ can always be related to that in ${{\bf k}}$ by either: [*(i)*]{} a rotation through a Pauli matrix if $g'_1g'_2=1$ in Eq. (\[loop2\]); or [*(ii)* ]{} a reflection in the loop plane if $g'_1g_2'=-1$. If [*(i)*]{} holds, then one can again rotate the effective Hamiltonian according to equations (\[rotatedloops\])-(\[Amatrix\]), and the resulting spectra can be obtained from Eqs (\[CDWk0\])-(\[CDWkc\]) for model-1 loops, with the replacement $p_{\parallel} - p_o\rightarrow p_{\parallel} - p_o+p_\perp$. But in case [*(ii)*]{} the two NLs can be related through a reflection in the NL’s plane,

$$\begin{aligned}
H_0({{\bf k}}) &=&  {\cal R} H_0({{\bf k}}+\ {{\bf Q}})  {\cal R}^\dagger\,,\\
{\cal R} &=&  (p_\perp\rightarrow - p_\perp) \tau_{\beta}\,, \end{aligned}$$

\[reflection\]

which corresponds to the following cases, depending on $\tau_{\beta}$: $$\begin{aligned}
\begin{aligned}
& H_0({{\bf k}}+\ {{\bf Q}})=\\
&= \left( p_{\parallel} - p_o  - p_\perp\right)\tau_a - p_\perp \tau_b\equiv H_0(-k_\perp) &\beta&=0;\\
&= \left( p_{\parallel} - p_o  - p_\perp\right)\tau_a + p_\perp \tau_b &\beta&=a;\\
&= \left[ -\left( p_{\parallel} - p_o\right)  + p_\perp\right]\tau_a - p_\perp \tau_b & \beta&=b;\\
&= \left[ -\left( p_{\parallel} - p_o\right)  + p_\perp\right]\tau_a + p_\perp \tau_b & \beta&=c\neq a,b\,.
\end{aligned}\end{aligned}$$ We apply the same rotation to the effective Hamiltonian, using Eqs. (\[rotatedloops\])-(\[Amatrix\]): $$\begin{aligned}
A H_{eff}({{\bf k}})  A^\dagger &=& \left( \begin{array}{cc}  
H_0({{\bf k}}) & U\bar m \tau_\alpha\tau_{\beta}  \\ U \bar m \tau_{\beta}\tau_\alpha &  H_0(-k_\perp)
 \end{array}   \right)\label{difficult}\,.
 \end{aligned}$$

For finite $\delta$ one cannot write the energy dispersion in closed form. We analytically deal with the degenerate case at perfect compensation, $\delta=0$, below, and show numerical results for nonzero $\delta$ in Fig.\[fig1\]. The figure shows the two inner bands of Hamiltonians (\[mode2\_case1\]), (\[mode2\_case2\]), (\[mode2\_case3\]), and (\[mode2\_case4\]) in the two-dimensional space of $(p_{\parallel} ,p_\perp)$. A NL for $U\bar m = \delta=0$ then looks like a Dirac cone at point $(p_o ,0)$. The splitting of the original NLs can be seen as the appearance of two Dirac cones in the plot. For finite $\delta$, the Dirac cone axis is tilted.

Similarly to the discussion for model-1 loops, we also list all the four possibilities.

1.  If $\tau_\alpha\tau_{\beta}=1$: $$\begin{aligned}
    A H_{eff}({{\bf k}})  A^\dagger &=& 
    \left( p_{\parallel} - p_o  \right)\tau_a +  p_\perp t_3\tau_a +  p_\perp t_3\tau_b\nonumber\\
    &+& U\bar m \ t_1 - \delta t_3\,,\label{mode2_case1}
     \end{aligned}$$ For $\delta=0$ (perfect compensation) the spectrum obeys $$\begin{aligned}
    E^2 = p_\perp^2 + \left[ \sqrt{p_\perp^2 + (U\bar m)^2} \pm \left(  p_{\parallel} - p_o\right)
    \right]^2
     \end{aligned}$$ which has two nodal lines for $p_\perp=0$ and $p_{\parallel} - p_o=\pm U\bar m$. By turning on $\delta$, the two NLs are tilted along the $p_{\perp}$ direction, and move along the $p_{\parallel}$ direction, as shown in Fig. \[fig1\](a). It can also be seen from Eq. (\[mode2\_case1\]) that if $p_\perp=0$ the dispersion relation has two nodal lines: $|p_{\parallel} - p_o|=\sqrt{U^2\bar m^2+\delta^2}$. Therefore, one of the loops will shrink into the origin as $p_{\parallel}\rightarrow 0$ when $\delta^2+U^2\overline{m}^2\rightarrow p_0^2$, and become gapped for larger $\delta$.

2.  If $\tau_\alpha\tau_{\beta}\propto\tau_a$: $$\begin{aligned}
    A H_{eff}({{\bf k}})  A^\dagger &=& 
    \left( p_{\parallel} - p_o  \right)\tau_a +  p_\perp t_3\tau_a +  p_\perp t_3\tau_b\nonumber\\
    &+& \epsilon_{k\alpha a}U\bar m \ t_2\tau_a - \delta t_3\,,\label{mode2_case2}
     \end{aligned}$$ For $\delta=0$ (perfect compensation) the spectrum obeys $$\begin{aligned}
    E^2 &=&   \left(  p_{\parallel} - p_o\right)^2 + 2 p_\perp^2 +(U\bar m)^2 \nonumber\\
    & \pm &
    2 \sqrt{ 
    \left(  p_{\parallel} - p_o\right)^2 \left(   p_\perp^2 +U^2\bar m^2\right)  + p_\perp^2 U^2\bar m^2
    }
     \end{aligned}$$ which, for $\delta=0$, has the same two nodal lines as in the previous case, and behavior of these lines with nonzero $\delta$ is also identical. However, these loops show a quadratic dispersion along the $p_\perp$ direction, as shown in Fig. \[fig1\](b).

3.  If $\tau_\alpha\tau_{\beta}\propto\tau_b$: $$\begin{aligned}
    A H_{eff}({{\bf k}})  A^\dagger &=& 
    \left( p_{\parallel} - p_o  \right)\tau_a +  p_\perp t_3\tau_a +  p_\perp t_3\tau_b\nonumber\\
    &+& \epsilon_{k\alpha b}U\bar m \ t_2\tau_b - \delta t_3\,, \label{mode2_case3}\\
    E^2 &=& p_\perp^2 + \left(
     \sqrt{  \left(  p_{\parallel} - p_o\right)^2 + U^2\bar m^2} \pm  p_\perp
     \right)^2
     \end{aligned}$$ which is a fully gapped insulator, and remains gapped with nonzero $\delta$ \[Fig. \[fig1\](c)\].

4.  If $\tau_\alpha\tau_{\beta}\propto\tau_c(\neq a,b)$: $$\begin{aligned}
    A H_{eff}({{\bf k}})  A^\dagger &=& 
    \left( p_{\parallel} - p_o  \right)\tau_a +  p_\perp t_3\tau_a +  p_\perp t_3\tau_b\nonumber\\
    &+& \epsilon_{k\alpha c}U\bar m \ t_2\tau_c - \delta t_3\,,\label{mode2_case4}\end{aligned}$$ which for $\delta=0$ has the spectrum $$\begin{aligned}
    E^2 &=&   \left(  p_{\parallel} - p_o\right)^2 + 2 p_\perp^2  + U^2\bar m^2 \nonumber\\
    &\pm&
    2p_\perp\sqrt{
     \left(  p_{\parallel} - p_o\right)^2 + 2U^2\bar m^2
    }\end{aligned}$$ with two nodal lines at $p_{\parallel} - p_o=0$ and $\sqrt{2}p_\perp=\pm U\bar m$. Unlike the nodal lines in the first two cases, these two are tilted along $p_{\parallel}$ and move away from each other along $p_{\perp}$ for $\delta\neq 0$.

![The energy dispersion of the two inner bands of the spectrum. Panels (a)-(d) correspond to the four cases in Eqs. (\[mode2\_case1\]), (\[mode2\_case2\]), (\[mode2\_case3\]), and (\[mode2\_case4\]) respectively. The energy offset $\delta$ is: (a1)-(d1) $\delta=0$; (a2)-(d2) $\delta=0.5$; and (a3)-(d3) $\delta=1$. The other parameters are: $p_o=1$, and $U\bar m =0.5$. []{data-label="fig1"}](fig1a.pdf "fig:"){width="1\linewidth"} ![The energy dispersion of the two inner bands of the spectrum. Panels (a)-(d) correspond to the four cases in Eqs. (\[mode2\_case1\]), (\[mode2\_case2\]), (\[mode2\_case3\]), and (\[mode2\_case4\]) respectively. The energy offset $\delta$ is: (a1)-(d1) $\delta=0$; (a2)-(d2) $\delta=0.5$; and (a3)-(d3) $\delta=1$. The other parameters are: $p_o=1$, and $U\bar m =0.5$. []{data-label="fig1"}](fig1b.pdf "fig:"){width="1\linewidth"} ![The energy dispersion of the two inner bands of the spectrum. Panels (a)-(d) correspond to the four cases in Eqs. (\[mode2\_case1\]), (\[mode2\_case2\]), (\[mode2\_case3\]), and (\[mode2\_case4\]) respectively. The energy offset $\delta$ is: (a1)-(d1) $\delta=0$; (a2)-(d2) $\delta=0.5$; and (a3)-(d3) $\delta=1$. The other parameters are: $p_o=1$, and $U\bar m =0.5$. []{data-label="fig1"}](fig1c.pdf "fig:"){width="1\linewidth"} ![The energy dispersion of the two inner bands of the spectrum. Panels (a)-(d) correspond to the four cases in Eqs. (\[mode2\_case1\]), (\[mode2\_case2\]), (\[mode2\_case3\]), and (\[mode2\_case4\]) respectively. The energy offset $\delta$ is: (a1)-(d1) $\delta=0$; (a2)-(d2) $\delta=0.5$; and (a3)-(d3) $\delta=1$. The other parameters are: $p_o=1$, and $U\bar m =0.5$. []{data-label="fig1"}](fig1d.pdf "fig:"){width="1\linewidth"}

Example models {#examples}
==============

In this section we investigate the density wave phases in several specific lattice models which have been widely studied in the literature. The examples in Secs. \[model1\_example\] and \[model2\_example\] fall into the two types of NLs discussed in previous sections. The remainder of this section will be devoted to model beyond the simplest two-band cases.

A model-1 loop example {#model1_example}
----------------------

A lattice model with two “model-1" loops can be described by the Hamiltonian, $$\begin{aligned}
H_0({{\bf k}})&=& \left( \cos k_x + \cos k_y - m \right)\tau_a + \cos{k_z} \tau_b\label{model1}\,.\end{aligned}$$ When $|m|<2$, the system has two NLs at $\cos{k_x}+\cos{k_y}-m=0$, in two parallel planes $k_z=\pm \pi/2$. They are nested by the vector $\mathbf{Q}=(0,0,\pi)$, and can be mapped to each other as $$\begin{aligned}
H_0({{\bf k}})=\tau_aH_0({{\bf k}}+{\bf Q})\tau_a\,.\end{aligned}$$ Introducing a Hubbard interaction, the effective Hamiltonian can be written as $$\begin{aligned}
H_{eff}({{\bf k}}) &=&
 \left( \cos k_x + \cos k_y - m \right)t_0\tau_a +  \cos{k_z} t_3\tau_b \nonumber\\ 
 &+& U\bar m\  t_1\tau_{\alpha}\,.\end{aligned}$$

We analyze first the CDW case ($\tau_\alpha=\tau_0$). Following the discussion in the previous Section, this condition corresponds to the second case in Sec. \[model1\_general\], where each NL gets split into two. Indeed, from the commuting relations between different terms, we can obtain the energy dispersion as $$E_{CDW}= \pm \sqrt{\cos^2 {k_z} + \left[
\cos k_x + \cos k_y - m \pm U\bar m
\right]^2
}\,,$$ and the split NLs are given by $$\begin{aligned}
\cos{k_z}  = 0\,,\mbox{and}\
\cos k_x + \cos k_y = m \pm U\bar m.\end{aligned}$$

In the PSDW case ($\tau_\alpha=\tau_3$), different choices of $a$ and $b$ in Eq. (\[model1\]) will lead to different phases of the system. In such case, the effective Hamiltonian reads $$\begin{aligned}
H_{eff}({{\bf k}}) &=&
 \left( \cos k_x + \cos k_y - m \right)t_0\tau_a +  \cos{k_z} t_3\tau_b \nonumber\\ 
 &+& U\bar m\  t_1\tau_3,\end{aligned}$$ and the possibilities are summarized in Table \[tab1\], and listed explicitly as follows.

   a$\backslash$ b     1       2       3
  ----------------- ------- ------- --------
          1            X     loops   gapped
          2          loops     X     gapped
          3          metal   metal     X

  : Possibilites for $\tau_a$, $\tau_b$ matrices in the loop model (\[model1\]). And the resulting different outcomes of the PSDW phase.[]{data-label="tab1"}

In the case of $\tau_a=\tau_3$, the energy spectrum reads $$\begin{aligned}
E_{PSDW}&=&\pm\sqrt{ 
\left( \cos k_x + \cos k_y - m \right)^2 +  \cos {k_z}  ^2 }\pm U\bar m 
 \nonumber\\\end{aligned}$$ with uncorrelated $\pm$ signs. This is the first case in Sec. \[model1\_general\], and the ordered phase is metallic.

If $a,b\neq 3$, the spectrum reads $$\begin{aligned}
E_{PSDW}^2&=&  
\left( \cos k_x + \cos k_y - m \right)^2 + \left(    \cos{k_z}    \pm U\bar m  \right)^2\,,\end{aligned}$$ where each original NL splits into two with different $k_z$, as in the third cases in Sec. \[model1\_general\].

Finally, when $\tau_b=\tau_3$, the spectrum takes the form, $$\begin{aligned}
E^2&=& 
\left( \cos k_x + \cos k_y - m \right)^2 +  \cos {k_z} ^2 + U^2\bar m^2\,.\end{aligned}$$ Thus the system is fully gapped by a nonzero $U \bar m$ and becomes an insulator, as in the fourth case in Sec. \[model1\_general\].

A model-2 loop example {#model2_example}
----------------------

By including $\cos{k_z}$ in the $\tau_a$ term in Eq. (\[model1\]), we obtain a system with “model-2" loops, described by $$\begin{aligned}
H_0({{\bf k}})= \left( \cos k_x + \cos k_y + \cos k_z- m \right)\tau_a + \cos k_z  \tau_b.\label{model2}\end{aligned}$$ This system has two parallel NLs with $k_z=\pm\pi/2$ when $-2<m<2$. The continuous approximation for these two loops, as in Eqs. (\[loop2\]), satisfies $g'_1=1$, and $g_2=g'_2=-1$. Thus, the effective Hamiltonian including the Hubbard interaction is given by $$\begin{aligned}
H_{eff} &=&  \left( \cos k_x + \cos k_y - m \right)t_0\tau_a +   \cos k_z  t_3\tau_a\nonumber\\ &+&
\cos k_z  t_3\tau_b + U\bar m t_1\tau_{\alpha}\,.
\label{16}\end{aligned}$$

For the CDW case ($\tau_{\alpha}=\tau_0$), the spectrum reads $$\begin{aligned}
E^2_{CDW}&=& \left[
\cos k_x + \cos k_y  - m \pm \sqrt{U^2\bar m^2+\cos^2{k_z}}
\right]^2\nonumber\\
&&+\cos {k_z}^2 
\,,\end{aligned}$$ thus each of the original NLs splits into two. The new NLs are given by $$\begin{aligned}
\cos{k_z}  = 0\,, \hspace{0.2cm}
\cos k_x + \cos k_y = m \pm U\bar m.\end{aligned}$$

In the case of PSDW ordering ($\tau_{\alpha}=\tau_3$), there are three possible outcomes as summarized in Table \[tab2\]. In the case of $\tau_a=\tau_3$, the energy spectrum reads $$E^2_{PSDW}=(\cos k_x+\cos k_y-m)^2+2\cos^2 k_z+U^2\bar m^2$$ $$\begin{aligned}
\pm2\sqrt{(\cos k_x+\cos k_y-m)^2\left[  \cos^2 k_z + U^2\bar m^2  \right] + (U\bar m\cos k_z)^2 }\,,\nonumber\\\end{aligned}$$ where each NL splits into two as in the second case in Sec. \[model2\_general\]. The condition of the NLs is the same as for the CDW case; however, as discussed in the previous section \[Fig. \[fig1\](b)\], these NLs have a quadratic dispersion along $k_z$.

   a$\backslash$ b     1       2       3
  ----------------- ------- ------- --------
          1            X     loops   gapped
          2          loops     X     gapped
          3          loops   loops     X

  : Possibilites for $\tau_a$, $\tau_b$ matrices in the loop model (\[model2\]). And the resulting different outcomes of the PSDW phase.[]{data-label="tab2"}

If $\tau_b=\tau_3$ then the spectrum reads $$\begin{aligned}
E_{PSDW}^2&=& \left[
\cos k_z\pm \sqrt{ U^2\bar m^2+ \left(   \cos k_x + \cos k_y -m\right)^2}
\right]^2\nonumber\\
 &+& \cos^2k_z \,,\end{aligned}$$ which is an insulating phase as in the third case in Sec. \[model2\_general\].

The remaining possibility, $a,b\neq 3$, where the spectrum reads $$\begin{aligned}
E^2_{PSDW}&=&(\cos k_x+\cos k_y-m)^2+2\cos^2 k_z+U^2\bar m^2\nonumber\\
&&\pm2\cos k_z\sqrt{
(\cos k_x+\cos k_y-m)^2
+2U^2\bar m^2
}\,,\nonumber\\\end{aligned}$$ yields two NLs given by $$\begin{aligned}
\cos{k_z}  = \pm U\bar m /\sqrt{2}\,, \
\cos k_x + \cos k_y = m\,.\end{aligned}$$ This is the splitting along $k_z$, as in the fourth case in Sec. \[model2\_general\]. Finally, if we add an extra term $H_{\delta}=\delta\sin{k_z}\tau_0$ to the original Hamiltonian of Eq. (\[model2\]), we can induce a energy offset $\delta$ of the two original NLs, and tilt the resulting NLs after the Hubbard interaction is introduced.

Nested $\mathbb{Z}_2$ NLs  {#sec:Z2}
--------------------------

We have hitherto considered examples of two-band NLSMs, which verify our results for general two-band models. If an extra degree of freedom, say, a (pseudo)spin-1/2 subspace, is introduced, the Hilbert space is enlarged and the possibility of nodal lines carrying a $\mathbb{Z}_2$ monopole charge arises, which must be created in pairs[@Z2_loop]. It is beyond the scope of this paper to extend the general analysis of Sec III. Instead, we explicitly consider a recent four-band model for $Z_2$ NLs [@Z2_loop] and study density-wave order due to nesting. The model reads $$\begin{aligned}
H_0({{\bf k}}) &=&  \sin k_x \tau_0\sigma_1+ \sin k_y \tau_2\sigma_2+ 
\sin k_z   \tau_0\sigma_3 + m \tau_1\sigma_1\,,\nonumber\\ 
\label{Z2_H}\end{aligned}$$ and has the spectrum $$\begin{aligned}
E_0^2({{\bf k}})&=&
\left(   \sqrt{\sin^2k_x + \sin^2k_y}\pm m  \right)^2 + \sin^2k_z\,.\label{spectrum_Z2}\end{aligned}$$ There are eight NLs, centered at momenta ${{\bf k}}$ with Cartesian components $k_i=0,\pi$ for $i=x,y,z$. We introduce a Hubbard interaction in four-dimensional space assuming that the repulsion exists between the two orbitals in subspace of $\boldsymbol\sigma$: $$\begin{aligned}
\hat U&=&  U  \sum_{{{\bf r}}, \nu}n_{\nu 1}({{\bf r}}) n_{\nu 2}({{\bf r}}) = \frac{U}{2}  \sum_{{{\bf r}}}  n({{\bf r}})\tau_0\sigma_1 n({{\bf r}})\,,\end{aligned}$$ so that the remaining index $\nu=1,2$ for the two components in the subspace of $\boldsymbol\tau$ produces a two-fold degeneracy. In the mean-field approximation, the interaction reads (apart from unimportant constants) $$\begin{aligned}
\hat U_{eff}&=&  U\bar m \sum_{{{\bf k}}} \hat c_{{{\bf k}}+{{\bf Q}}}^\dagger \tau_0\sigma_z \hat c_{{{\bf k}}} \,,\end{aligned}$$ where we defined the PSDW as: $$\langle   n_{\nu j}({{\bf r}})\rangle = \frac{n}{2} + \bar m (-1)^j\cos({{\bf Q}}\cdot{{\bf r}}) \,,$$ The effective 8x8 Hamiltonian has the same form as that in Eq. (\[martizHloops\]), but the anti-diagonal blocks are now written as $U\bar m \tau_0\sigma_3$. Two of the original $\mathbb{Z}_2$ loops are now coupled by the interaction. Such coupling provides an extra pseudospin-1/2 subspace, and the interaction term breaks the SU(2) symmetry in this space. As a result, the pair of NLs may either survive, shrink to point nodes, or be gapped out by the PSDW.

Depending on vector ${{\bf Q}}$, the effective Hamiltonian takes the form of $$\begin{aligned}
H_{eff} &=& \sin k_x t_{\alpha_x} \tau_0\sigma_1 + \sin k_y t_{\alpha_y} \tau_2\sigma_2 + \sin k_z t_{\alpha_z} \tau_0\sigma_3 \nonumber\\ 
&+&m t_0 \tau_1\sigma_1 + U\bar m t_1 \tau_0\sigma_3\,,\end{aligned}$$ where the Pauli matrix $t_{\alpha_i}$ equals $t_0$ ($t_3$) for $Q_i=0$ ($Q_i=\pi$). Thus, different choices of $Q_i$ determine whether each term commutes or anticommutes with the interaction term $U\bar m t_1 \tau_0\sigma_3$. The possibilities with different choices of ${{\bf Q}}$ are summarized as follows.

![The effect of a PSDW on the $\mathbb{Z}_2$ loops for different nesting vectors ${{\bf Q}}$. (a) to (d) show four different cases of ${{\bf Q}}$ where the resulting system is still a semimetal. The red lines are the original NLs of Hamiltonian $H_0({{\bf k}})$ in Eq. (\[Z2\_H\]), whereas the blue lines and stars are the nodal lines or points in the PSDW phase. The parameters are: $m=\sqrt{2}/2$ and $U\bar m=0.5$. We only plot for $k_i\in [-\pi/2,\pi/2]$, as the graph repeats itself with a period of $\pi$.\[fig2\]](fig2.pdf){width="1\linewidth"}

1.  ${{\bf Q}}=(1,1,1)\pi$: $$\begin{aligned}
    H_{eff} &=& \sin k_x t_3 \tau_0\sigma_1 + \sin k_y t_3 \tau_2\sigma_2 + \sin k_z t_3 \tau_0\sigma_3 \nonumber\\ &+&
    m t_0 \tau_1\sigma_1 + U\bar m t_1 \tau_0\sigma_3\,, 
    \\
    E^2 &=& \left(  \sqrt{\sin^2 k_x + \sin^2k_y} \pm  \sqrt{m^2  +U^2\bar m^2}\right)^2 + \sin^2k_z\,. \nonumber\\ 
    \label{E111}\end{aligned}$$ The spectrum is composed of four two-fold degenerate bands, and has the same NLs as the original $H_0$, albeit with enlarged radius \[Fig. \[fig2\](a)\]: $$\begin{aligned}
    k_z=0,\pi; \   \sqrt{\sin^2 k_x + \sin^2k_y} =  \sqrt{m^2  +U^2\bar m^2}\,.
    \label{eq56}\end{aligned}$$ Since the energy dispersions in Eq. (\[E111\]) are two-fold degenerate, the NLs are four-fold degenerate.

2.  ${{\bf Q}}=(1,1,0)\pi$: $$\begin{aligned}
    H_{eff} &=& \sin k_x t_3 \tau_0\sigma_1 + \sin k_y t_3 \tau_2\sigma_2 + \sin k_z t_0 \tau_0\sigma_3  \nonumber\\ &+&
    m t_0 \tau_1\sigma_1 + U\bar m t_1 \tau_0\sigma_3\,,\end{aligned}$$ which has the same NLs as in Eq. (\[eq56\]), and as shown in FIg. \[fig2\](b), despite some minor variation of the spectrum near the NLs.

3.  ${{\bf Q}}=(1,0,1)\pi$: $$\begin{aligned}
    H_{eff} &=& \sin k_x t_3 \tau_0\sigma_1 + \sin k_y t_0 \tau_2\sigma_2 + \sin k_z t_3 \tau_0\sigma_3 \nonumber\\ &+&
    m t_0 \tau_1\sigma_1 + U\bar m t_1 \tau_0\sigma_3\,.\end{aligned}$$ The full spectrum also has the two-fold degeneracy. while the NLs are gapped in most regions, leaving only pairs of four-fold degenerate points at $$\begin{aligned}
    \sin k_y = \sin k_z=0; \   \sin k_x  = \pm \sqrt{m^2  +U^2\bar m^2}\,,\end{aligned}$$ as shown in Fig. \[fig2\](c).

    We note that the case ${{\bf Q}}=(0,1,1)\pi$ has a similar spectrum, which can be obtained from the above by performing the substitution $k_x\leftrightarrow k_y$.

4.  ${{\bf Q}}=(1,0,0)\pi$: $$\begin{aligned}
    H_{eff} &=& \sin k_x t_3 \tau_0\sigma_1 + \sin k_y t_0 \tau_2\sigma_2 + \sin k_z t_0 \tau_0\sigma_3 \nonumber\\ &+&
    m t_0 \tau_1\sigma_1 + U\bar m t_1 \tau_0\sigma_3\,.\end{aligned}$$ Energy zeros are obtained if two conditions are simultaneously satisfied:

    $$\begin{aligned}
    \sum_{i=1}^3 \sin^2k_i =  U^2\bar m^2 + m^2\\
    m\sin k_z=\pm  U\bar m \sin k_y\,.\end{aligned}$$

    Geometrically, one can think of the first condition as defining a spherical surface, for small ${{\bf k}}$, and the second one as two planes. The intersection between the surface and the planes yields two NLs, obtained by rotating the original loops (in the XY plane) around the $x$ axis in opposite directions. These two NLs are given by different pairs of energy bands, and thus form a nodal chain \[Fig. \[fig2\](d)\].

    We note that spectrum for the case ${{\bf Q}}=(0,1,0)\pi$ is similar to the case ${{\bf Q}}=(\pi,0,0)$, and differs only in the interchange $k_x\leftrightarrow k_y$. The nodal chain is obtained from the two original NLs by a rotation around the $y$ axis.

5.  ${{\bf Q}}=(0,0,1)\pi$ $$\begin{aligned}
    H_{eff} &=& \sin k_x t_0 \tau_0\sigma_1 + \sin k_y t_0 \tau_2\sigma_2 + \sin k_z t_3 \tau_0\sigma_3 \nonumber\\ &+&
    m t_0 \tau_1\sigma_1 + U\bar m t_1 \tau_0\sigma_3\,.\end{aligned}$$ In this case the system is a fully gapped insulator.

Nested NLs in Dirac systems {#Diracsec}
===========================

In this Section we study density wave phases in four-band spinfull Hamiltonians with NLs. The spin degree of freedom allows us to distinguish two types of ordered phases: true and hidden spin density waves (SDWs). In subsection \[Diracperturb\] we consider NLs obtained from a perturbed[@Burkov; @mitschell] Dirac Hamiltonian.

Spin degenerate loops
---------------------

The simplest way to go from a Weyl to a Dirac loop is to introduce spin degeneracy, $$H_0({{\bf k}})\  \rightarrow\ H_0({{\bf k}}) \sigma_0$$ where $\sigma_\mu$ acts in spin space.

We assume Hubbard repulsion between two fermions having opposite spins in the same orbital (labeled by the index $j$): $$\begin{aligned}
\hat U = \sum_{{{\bf r}}, j=1,2} \hat n_{j{\uparrow}}({{\bf r}})\hat n_{j{\downarrow}}({{\bf r}})\,.
\label{Usigma}\end{aligned}$$ In principle one could have ordered phases with ferromagnetic (FM) or antiferromagnetic (or SDW) configurations of the spin: $$\begin{aligned}
\langle n_{j\sigma}\rangle &=& 
  \frac 1 4 n + \bar m  \sigma  \qquad
  \sigma=\pm 1\,,\qquad\mbox{Stoner FM}\\
  \langle n_{j\sigma}\rangle &=& 
  \frac 1 4 n + \bar m \sigma (-1)^j   \qquad\mbox{hidden Stoner FM}\\
 \langle n_{j\sigma}\rangle &=& 
  \frac 1 4 n + \bar m  \sigma \cos({{\bf Q}}\cdot{{\bf r}}) \qquad\mbox{true SDW}\\
  \langle n_{j\sigma}\rangle &=& 
  \frac 1 4 n + \bar m  \sigma (-1)^j \cos({{\bf Q}}\cdot{{\bf r}})   \qquad\mbox{hidden  SDW} \label{hiddenSDW}
 \end{aligned}$$ Because a NL’s density of states vanishes linearly with energy, the Stoner criterion precludes the FM orderings for weak interactions[@roy], and they are not related to the nesting ${{\bf Q}}$. In the following, we shall concentrate on SDW phases. Considering the hidden SDW, Eq. (\[hiddenSDW\]), the effective interaction reads: $$\begin{aligned}
\hat U_{eff} &=&
-U\sum_{{{\bf r}},j} \left[ \left( \frac n 4  \right)^2  -   \bar m^2  \cos^2({{\bf Q}}\cdot{{\bf r}}) 
- \frac n 4 \sum_{\sigma} \hat\psi^\dagger_{j\sigma}({{\bf r}}) \hat\psi_{j\sigma}({{\bf r}}) 
\right]\nonumber\\
&+&U\bar m\sum_{{{\bf r}},j,j',\sigma\sigma'} \hat\psi^\dagger_{j\sigma}({{\bf r}})\
  \tau^{jj'}_3\sigma^{\sigma\sigma'}_3  \cos({{\bf Q}}\cdot{{\bf r}})\   \hat\psi_{j'\sigma'}({{\bf r}}) \,,\label{tSDW}\end{aligned}$$ where the field operator $\hat\psi_{j\sigma}({{\bf r}})$ now includes the spin index $\sigma$. Similarly, the field operator in momentum space $\boldsymbol c_{{\bf k}}$ now denotes $ c_{{{\bf k}},j,\sigma}$ for all $j,\sigma$. The Hamiltonian matrix in $\left(  \boldsymbol c^\dagger_{{{\bf k}}}\ \boldsymbol c^\dagger_{{{\bf k}}+{{\bf Q}}}  \right)$ space reads (apart from unimportant constants), $$\begin{aligned}
H_{eff}({{\bf k}}) &=&  \left(  \begin{array}{cc}
H_0({{\bf k}}) & U\bar m \tau_\alpha\sigma_3 \\
 U\bar m \tau_\alpha\sigma_3  & H_0({{\bf k}}+{{\bf Q}})  
\end{array}\right)\,,\label{SDWDirac}\end{aligned}$$ where $\alpha=3$ describes a hidden SDW \[as in Eq. (\[tSDW\])\], and $\alpha=0$ describes a true SDW. The off-diagonal block $U\bar m t_1\tau_3\sigma_3$ has exactly the same (anti)commutation relations with the other Hamiltonian terms, as in the Weyl case of Sections \[model1\_general\] and \[model2\_general\] . All the criteria and spectra established for the Weyl case still hold, if one just replaces $\bar m\rightarrow\bar m \sigma_3$. Since this term appears as $(U\bar m)^2$ in the dispersion relations, there is no spin splitting in the spectra.

We note that single antiferromagnetic NLs have been discussed in the literature[@jingwang].

Two perturbed Dirac points {#Diracperturb}
--------------------------

A nodal line Dirac semimetal can be obtained starting from a pristine 3D Dirac semimetal[@mitschell] of the form $H_D({{\bf k}})  = -\tau_3{{\bf p}}\cdot\boldsymbol\sigma$ and perturbing it with terms of the form $a_\mu\tau_\mu\otimes b_\nu\sigma_\nu$. Suppose, for instance, $$H_0({{\bf k}})  = -\tau_3{{\bf p}}\cdot\boldsymbol\sigma + \tau_1\boldsymbol b\cdot\boldsymbol\sigma$$ Without loss of generality assume $\boldsymbol b \parallel \hat z$. The term $p_z\tau_3\sigma_3$ anticommutes with the others, so: $$E^2 =p_z^2 + \left( \sqrt{p_x^2+p_y^2}\pm b   
 \right)^2 \,.$$ We note that this dispersion relation is very similar to that in Eq. (\[spectrum\_Z2\]) for the $\mathbb{Z}_2$ loops. However, the effect of the Hubbard interaction is different, as these two systems couple the two sets of (pseudo)spin-1/2 subspace in different ways. On the other hand, Dirac points described by $H_D$ above do not exist alone if an additional symmetry, such as time-reversal or inversion, is present. For instance, a time-reversal symmetry (TRS) relates two Dirac points at $-\frac 1 2 {{\bf Q}}$ and $+\frac 1 2 {{\bf Q}}$ in such a way that: $$\begin{aligned}
\sigma_2H_0^*\left(\frac {{\bf Q}}2 -{{\bf k}}\right)\sigma_2 = H_0\left(-\frac {{\bf Q}}2 +{{\bf k}}\right)\,,\end{aligned}$$ it then follows that, for $\boldsymbol b=0$, $H_0\left(- {{\bf Q}}/2 +{{\bf k}}\right) =H_0\left( {{\bf Q}}/2 +{{\bf k}}\right) =H_D({{\bf k}})$. Therefore, the two unperturbed Dirac points have the same $k\cdot p$ Hamiltonian. Including the $\tau_1\boldsymbol b\cdot\boldsymbol\sigma$ term, which breaks TRS, we obtain the model,

$$\begin{aligned}
 H_0\left(-\frac {{\bf Q}}2 +{{\bf k}}\right) &=& -\tau_3{{\bf p}}\cdot\boldsymbol\sigma +\tau_1\boldsymbol b\cdot\boldsymbol\sigma\,,\\
 H_0\left(\frac {{\bf Q}}2 +{{\bf k}}\right) &=& -\tau_3{{\bf p}}\cdot\boldsymbol\sigma +\tau_1\boldsymbol b\cdot\boldsymbol\sigma\,,\end{aligned}$$

\[Tpartial\]

and the effective Hamiltonian has equal diagonal blocks.

A different version of the above model that would preserve TRS symmetry reads:

$$\begin{aligned}
 H_0\left(-\frac {{\bf Q}}2 +{{\bf k}}\right) &=& -\tau_3{{\bf p}}\cdot\boldsymbol\sigma +\tau_1\boldsymbol b\cdot\boldsymbol\sigma\,,\\
 H_0\left(\frac {{\bf Q}}2 +{{\bf k}}\right) &=& -\tau_3{{\bf p}}\cdot\boldsymbol\sigma -\tau_1\boldsymbol b\cdot\boldsymbol\sigma\,.\end{aligned}$$

\[Tcase\]

If we now consider the role of inversion symmetry ${{\bf k}}\rightarrow -{{\bf k}}$, the two Dirac points are related by

$$\begin{aligned}
H_0\left(\frac {{\bf Q}}2 -{{\bf k}}\right) = H_0\left(-\frac {{\bf Q}}2 +{{\bf k}}\right) &=& -\tau_3{{\bf p}}\cdot\boldsymbol\sigma + \tau_1\boldsymbol b\cdot\boldsymbol\sigma\,,\hspace{1cm} \\ 
\Rightarrow
H_0\left(\frac {{\bf Q}}2 +{{\bf k}}\right) &=& \tau_3{{\bf p}}\cdot\boldsymbol\sigma +\tau_1\boldsymbol b\cdot\boldsymbol\sigma\,,
\label{Pcase}\end{aligned}$$

therefore, the two Dirac points have different $k\cdot p$ Hamiltonians.

Next we study the effects of a hidden SDW and a true SDW for these different cases, still assuming $\boldsymbol b \parallel \hat z$. For a hidden SDW, the effective Hamiltonian for the TRS breaking model in Eq. (\[Tpartial\]) is then $$\begin{aligned}
\hat H_{eff} =  t_0 \left[-\tau_3{{\bf p}}\cdot\boldsymbol\sigma + b\tau_1 \sigma_3\right] +  U\bar m t_1\tau_3\sigma_3\,,\end{aligned}$$ which, by inspection produces the eight band spectrum: $$E^2 = \left( - p_z\pm U\bar m \right)^2 + \left(  \sqrt{p_x^2+p_y^2}\pm b   \right)^2$$ where the $\pm$ signs are uncorrelated. This corresponds to splitting each loop along the $k_z$ direction.

If one considers, instead, a true SDW phase, $$\begin{aligned}
\hat H_{eff} =  t_0 \left[-\tau_3{{\bf p}}\cdot\boldsymbol\sigma + b\tau_1\sigma_3\right]+  U\bar m t_1\tau_0\sigma_3\,,\end{aligned}$$ then there are four doubly degenerate bands: $$\begin{aligned}
E^2&=& b^2 +  U^2\bar m^2+ {{\bf p}}^2\nonumber\\ &\pm& 2\sqrt{ b^2\left(  U^2\bar m^2 + p_x^2 + p_y^2  \right) + U^2\bar m^2p_z^2
} \label{mathspec}\end{aligned}$$ with nodal lines given by $ p_z=0\,,  p_x^2 + p_y^2  = b^2 - U^2\bar m^2$. So, the initial two loops still exist but their radius shrinks.

In the TRS model, Eq. (\[Tcase\]), the hidden SDW phase is described by the effective Hamiltonian: $$\begin{aligned}
\hat H_{eff}&=&  -t_0\tau_3{{\bf p}}\cdot\boldsymbol\sigma + bt_3\tau_1\sigma_3 +  U\bar m t_1\tau_3\sigma_3\,.\end{aligned}$$ The spectrum is the same as in Eq. (\[mathspec\]). So, the initial two loops still exist but their radius shrinks. And a true SDW phase is described by the effective Hamiltonian: $$\begin{aligned}
\hat H_{eff} &=&  -t_0 \tau_3( p_x\sigma_1 + p_y\sigma_2) - p_z t_0\tau_3\sigma_3 + b\   t_3\tau_1\sigma_3 \nonumber\\ &+& U\bar m t_1\tau_0\sigma_3\,,\end{aligned}$$ which produces the spectrum with eight bands: $$E^2 =\left(    \sqrt{p_x^2+p_y^2} \pm b  \right)^2 + \left( p_z \pm     U\bar m  \right)^2\,,$$ where the $\pm$ signs are uncorrelated. This corresponds to splitting each loop along $p_z=\pm U\bar m$.

For the case with inversion symmetry, in Eq. (\[Pcase\]), a hidden SDW phase is described by the effective Hamiltonian: $$\begin{aligned}
H_{eff} ({{\bf k}})=  -t_3 \tau_3{{\bf p}}\cdot\boldsymbol\sigma + b \ t_0\tau_1\sigma_3 +  U\bar m t_1\tau_3\sigma_3\,,
\label{PhiddenSDW}\end{aligned}$$ which produces the eight band spectrum: $$E^2 =p_z^2 +  \left(    \pm \sqrt{b^2 + U^2\bar m^2} \pm   \sqrt{p_x^2+p_y^2}
  \right)^2\,,$$ (uncorrelated$\pm$ signs). This corresponds to splitting each nodal line by changing its radius. A true SDW is obtained by changing $\tau_3\rightarrow\tau_0$ in the last term of Eq. (\[PhiddenSDW\]). The resulting spectrum, $$E^2 =p_z^2 +  \left(   \sqrt{p_x^2+p_y^2} \pm b \pm  U\bar m   \right)^2\,,$$ (with uncorrelated $\pm$ signs) also has NLs given by $p_z=0$, $\sqrt{p_x^2+p_y^2} =|b\pm  U\bar m|$.

In the remaining case, where $H_0(-{{\bf Q}}/2 + {{\bf k}})=-H_0({{\bf Q}}/2  + {{\bf k}})$, a hidden SDW phase is described by the effective Hamiltonian: $$\begin{aligned}
H_{eff} ({{\bf k}})=  t_3 \left[ -\tau_3{{\bf p}}\cdot\boldsymbol\sigma + b \tau_1\sigma_3 \right]+  U\bar m t_1\tau_3\sigma_3\,,\end{aligned}$$ which, by inspection produces the eight band spectrum: $$E^2 =p_z^2 +  \left(    \sqrt{p_x^2+p_y^2} \pm b \pm  U\bar m \right)^2\,,$$ with uncorrelated$\pm$ signs. Therefore, each loop splits in the radial direction. The true SDW is described by the effective Hamiltonian: $$\begin{aligned}
\hat H_{eff} = t_3 \left[ -\tau_3{{\bf p}}\cdot\boldsymbol\sigma + b \tau_1\sigma_3 \right]+  U\bar m t_1\tau_0\sigma_3\,,\end{aligned}$$ which, by inspection produces the eight band spectrum: $$E^2 = p_z^2 +  \left(  \pm  \sqrt{p_x^2+p_y^2} \pm \sqrt{ b^2 + U^2\bar m^2}  \right)^2 \,,$$ where the $\pm$ signs are uncorrelated. Again, this corresponds to splitting each nodal loop in the radial direction.

Superconductivity {#supersec}
=================

When considering a single Weyl NL, the pairing block of the Bogolyubov-deGennes (BdG) matrix in the particle-hole basis[@sacramento1; @sacramento2; @beri], takes the form $$\hat \Delta({{\bf k}}) = 
 \left[ 
d_0({{\bf k}})\tau_0+ \boldsymbol d({{\bf k}})\cdot\boldsymbol\tau
 \right] i\tau_2
\,,$$ and fermionic statistics imposes that $\hat \Delta({{\bf k}}) = \hat \Delta^T(-{{\bf k}})$. Close to the nodal lines the 3D momentum, ${{\bf p}}=\hbar{{\bf k}}$, can be parametrized as

$$\begin{aligned}
  p_x &=& (p_0 + \tilde p\cos\phi) \cos\theta\\ 
  p_y &=& (p_0 + \tilde p\cos\phi) \sin\theta \,,\\
p_z&=&p_\perp=\tilde p\sin\phi\,, \end{aligned}$$

\[pcoordinates\]

which is to be inserted in the $k\cdot p$ loops models. Here, $\theta$ is the azimuthal angle along the loop, $\tilde p$ is the radius of a torus involving the NL, and the angle $\phi$ wraps around the latter[@Nandkishore1]. Note that, according to Eq. (\[pcoordinates\]), momentum inversion ${{\bf p}}\rightarrow -{{\bf p}}$ is equivalent to $\theta\rightarrow\theta+\pi$ and $\phi\rightarrow -\phi$, while reflection in the loop’s plane, $p_z\rightarrow -p_z$, is equivalent to $\phi\rightarrow -\phi$.

In the semimetal case (undoped, or compensated case) the FS reduces to the NL and ${{\bf p}}$ reduces to the angle $\theta$ on the loop. In the doped case, any point on the torus shaped FS can be labeled by two angles, $\theta,\phi$. The functions $ d_0({{\bf k}})$ and $\boldsymbol d({{\bf k}})$ describe (pseudo-spin) singlet and triplet pairing, respectively. One can expand the singlet pairing function quite generally as $$d_0({{\bf k}})  =  \sum_{l_1,l_2} e^{il_1\theta}\left[  \Delta_{l_1l_2}\cos(l_2\phi) + \tilde \Delta_{l_1l_2}\sin(l_2\phi)
\right]\,.$$ An analogous expansion can be written for $\boldsymbol d({{\bf k}})$.

If there are two nested Weyl loops, then an additional loop label must be introduced and the Pauli matrix $t_\mu$ operates in the two-dimensional loop space. For a two-loop system then, we write the pairing matrix as $$\hat \Delta({{\bf k}}) = 
 \left[ 
d_0({{\bf k}})\tau_0+ \boldsymbol d({{\bf k}})\cdot\boldsymbol\tau
 \right] i\tau_2
 t_\mu\,.$$ The BdG Hamiltonian matrix in the particle-hole basis has the form $$H({{\bf k}})= \left( \begin{array}{cc}
\hat\Xi ({{\bf k}})& \hat \Delta ({{\bf k}})\\
\hat\Delta^\dagger ({{\bf k}}) & -\hat\Xi^T(-{{\bf k}})
 \end{array}\right) 
 \label{BdGmatrix}$$ with $\hat\Xi= diag(H_1, H_2)$ . The Hamiltonians $H_{1(2)}$ are the $k\cdot p$ Weyl NL models. The total Hamiltonian is then $$\begin{aligned}
\hat H=\frac 1 2  \sum_{{{\bf k}}}{\boldsymbol c}^\dagger H({{\bf k}})
 {\boldsymbol c}\,,\end{aligned}$$ where $ {\boldsymbol c}=(\hat {\boldsymbol c}_{{{\bf k}},1},\   \hat  {\boldsymbol c}_{{{\bf k}},2},\ 
 \hat  {\boldsymbol c}_{-{{\bf k}},1}^\dagger,\   {\boldsymbol c}_{-{{\bf k}},2}^\dagger
)^T$. If the two NLs are centered at BZ points $\pm{{\bf Q}}/2$ respectively, then the inter-NL pairing is the “usual” pairing between opposite momenta, and we shall take this to be the case. If not, then the Cooper pair would have a finite quasi-momentum (a Fulde-Ferrel-Larkin-Ovchinnikov state[@LOFF1; @LOFF2; @LOFF3]).

The cases $\mu=0,1,3$, are different from the case $\mu=2$ regarding the parity of the functions $ d_0({{\bf k}})$ and $\boldsymbol d({{\bf k}})$. In the cases $\mu=1,2$, electrons on different loops are being paired: an electron $({{\bf k}},1)$ is being paired with another $(-{{\bf k}},2)$. The cases $\mu=0,3$ describe intra-NL pairing, where the scattering of two particles from one NL into the other may be included, and $id_0\tau_2t_3$ would describe sign-reversed s-wave pairing, analogous to that in pnictide superconductors[@BangChoi].

Inter-NL pairing with $\mu=1$ (interloop triplet pairing), requires $d_0$ to be even and ${\boldsymbol d}$ to be odd function of ${{\bf k}}$; if $\mu=2$ (interloop singlet), then $d_0$ and ${\boldsymbol d}$ have the opposite parities. The BdG matrix decouples into two blocks each associated with the vector spaces $(\hat c_{{{\bf k}},1}, \hat c_{-{{\bf k}},2}^\dagger)^T$ and $(\hat c_{{{\bf k}},2}, \hat c_{-{{\bf k}},1}^\dagger)^T$, respectively. Since we expect a fully gapped excitation spectrum to have higher condensation energy than a nodal spectrum, we shall examine the cases where $d_0$ and $\boldsymbol d$ are constant on the FS (in the $\mu=1$ and $\mu=2$ cases, respectively). If TRS holds, then these order parameters must also be real.

Model-1 loops {#model-1-loops}
-------------

Assuming a positive energy offset, $\delta$, the interband pairing occurs between the electronic toroidal FS from the $H_1-\delta$ loop, and the hole-like FS from the $H_2+\delta$ loop. As in previous literature, this is best done by considering projective form factors[@wangye; @Nandkishorepro] onto the conduction or valence band. Let $U_{1(2)}$ be the unitary matrices which diagonalize $H_{1(2)}$, so that $U_sH_sU_s^\dagger=\sqrt{\left( |{{\bf p}}_\parallel| - p_0 \right)^2 + p_\perp^2}\ \tau_3\equiv  \tilde p\tau_3$ for $s=1,2$. The positive and the negative branches are the conduction and valence bands, respectively. Because for model-1 loops there is always a Pauli matrix $\tau_\beta$ such that $H_1=\tau_\beta H_2\tau_\beta$, it then follows that $U_2=U_1\tau_\beta$. We can apply this same unitary transformation to the BdG matrix in Eq. (\[BdGmatrix\]) as: $$\begin{aligned}
\left( \begin{array}{cc}\Lambda & 0\\ 0 & \Lambda^*(-{{\bf k}})\end{array}\right)H({{\bf k}}) 
\left( \begin{array}{cc}\Lambda^\dagger & 0\\ 0 &\Lambda^T(-{{\bf k}}) \end{array}\right)
\nonumber\\ =
\left( \begin{array}{cc} \tilde p\tau_3t_0 -\delta\tau_0 t_3 & \Lambda \hat\Delta \Lambda^T(-{{\bf k}}) \\ 
\Lambda^*(-{{\bf k}})  \hat\Delta^\dagger \Lambda^\dagger& -\tilde p\tau_3t_0+\delta \tau_0 t_3 \end{array}\right)
\label{UBdGU}\end{aligned}$$ where $\Lambda=diag(U_1, U_2)$. The off-diagonal pairing block is then $\Lambda({{\bf k}}) \hat\Delta({{\bf k}}) \Lambda^T(-{{\bf k}})$.

For $\delta>0$, only the pairing between the conduction band of $H_1$ and the valence band of $H_2$ is considered. From the BdG matrix in Eq. (\[UBdGU\]) we obtain the submatrix operating in this two-fold subspace as: $$\begin{aligned}
H^{FS}=\left( \begin{array}{cc} \tilde p -\delta & 
\Delta_{FS}({{\bf k}}) \\ 
\Delta_{FS}^*({{\bf k}})
& \tilde p -\delta  \end{array}\right)
\label{condval1}\end{aligned}$$ where $\Delta_{FS}({{\bf k}})$ is the pairing function on the FS which, from Eq. (\[UBdGU\]) and for $\mu=1$ reads: $$\begin{aligned}
\Delta_{FS}({{\bf k}}) = 
 \left[
 U_1({{\bf k}})  \left(d_0+{\boldsymbol d}\cdot{\boldsymbol\tau}\right) i\tau_2U_2^T(-{{\bf k}})\right]_{12}
 \label{pairingfunction}\end{aligned}$$ It is then clear from Eq. (\[condval1\]) that the spectrum is $E=\tilde p -\delta \pm |\Delta_{FS}({{\bf k}})|$, and is gapless. At finite doping, no gapped state is to be expected from FS interloop pairing between non-degenerate model-1 Weyl loops.

The situation is different for the degenerate ($\delta=0$) case, however, where the FS is composed of two nodal lines. From Eq. (\[UBdGU\]) and $t_\mu=t_1$ we obtain a BdG matrix restricted to the subspace $\left(  U_1({{\bf k}})  \hat{\boldsymbol c}_{{{\bf k}},1} \,,    
U_2^*(-{{\bf k}}) \hat{\boldsymbol c}_{-{{\bf k}}, 2}^\dagger     \right)$, as $$\begin{aligned}
H_{12}'=\left( \begin{array}{cc}
\tilde p\tau_3  &  U_1(d_0 + {\boldsymbol d}\cdot{\boldsymbol \tau})i\tau_2  U_2^T  \\
-iU_2^*\tau_2 (d_0 + {\boldsymbol d}\cdot{\boldsymbol \tau}) U_1^\dagger  & -\tilde p\tau_3 
\end{array}\right)
\label{UBU2}\end{aligned}$$ For the sake of definiteness we consider the NL models with $\tau_a=\tau_1,\tau_b=\tau_2$, so that $$\begin{aligned}
  H_1(\phi) &=& \tilde p\left( \cos\phi\tau_1 + \sin\phi\tau_2\right)\,,\label{h1phi}\\
U_1({{\bf k}})&=& \frac{1}{\sqrt{2}}\left( \begin{array}{cc}
1 & e^{-i\phi}\\
-1 & e^{-i\phi}
\end{array}\right)\,.\label{u1phi}\\end{aligned}$$ We note that all the other $(\tau_a,\tau_b)$ cases can be related to this through a suitable rotation in pseudo-spin space. From Eq. (\[u1phi\]), one can see that $U_1^\dagger({{\bf k}}) =U_1^T(-{{\bf k}}) $. Interloop pairing is described by the off-diagonal block in Eq (\[UBU2\]): $$\begin{aligned}
 && U_1({{\bf k}})  \left(d_0+{\boldsymbol d}\cdot{\boldsymbol\tau}\right) i\tau_2\tau_\beta^TU_1^T(-{{\bf k}})=\nonumber\\
  &=& 
 \left(\begin{array}{cc}
 id_0\sin\phi + id_2+d_3\cos\phi &  d_0\cos\phi + d_1 + id_3\sin\phi\\
  -d_0\cos\phi + d_1 - id_3\sin\phi &  - id_0\sin\phi + id_2-d_3\cos\phi 
\end{array}\right) \,, \nonumber \\
 &=&  \left(\begin{array}{cc}
d_3 +id_2\cos\phi- id_1\sin\phi  & -d_0 - d_1\cos\phi- d_2\sin\phi  \\ 
-d_0 + d_1\cos\phi + d_2\sin\phi  & d_3 -id_2\cos\phi+  id_1\sin\phi \end{array}\right) \,, \nonumber\\
 &=& \left(\begin{array}{cc} 
 -id_0 -id_1\cos\phi-id_2\sin\phi   &   d_1\sin\phi - d_2\cos\phi +id_3\\
 - d_1\sin\phi + d_2\cos\phi +id_3 &  -id_0 +id_1\cos\phi+id_2\sin\phi   \end{array}\right)\,,  \nonumber\\
 &=&
  \left(\begin{array}{cc} 
-d_0\cos\phi  - d_1-id_3\sin\phi &   -id_0\sin\phi - id_2 -d_3\cos\phi\\
  id_0\sin\phi - id_2 +d_3\cos\phi &  d_0\cos\phi  - d_1+id_3\sin\phi  \end{array}\right)\,,\nonumber\\
  \label{udelu}\end{aligned}$$ for the cases $\beta=0,1,2,3$, respectively. Note that for the case $t_\mu=t_2$, we simply have to multiply both sides of Eq. (\[udelu\]) by $-i$.

For $\mu=1$ a fully gapped FS can only happen for constant $d_0$ because $\boldsymbol d$ is an odd function and must have nodes on the NLs. In this case, only for $\beta=2$ a gapped spectrum is obtained: $E^2 = \tilde p^2 + d_0^2$.

For $\mu=2$ (interloop singlet) and constant real $\boldsymbol d$ there are more possibilities. If $\beta=0$ a fully gapped spectrum $E^2 = \tilde p^2 + d_2^2$; if $\beta=1$ a fully gapped spectrum $E^2 = \tilde p^2 + d_3^2$; for $\beta=3$ the fully gapped spectrum $E^2 = \tilde p^2 + d_1^2$. Gapped spectra result from intraband pairing. Interband pairing leads to nodal spectra for the same reason as in the $\delta>0$ case.

Model-2 loops {#model-2-loops}
-------------

In a model-2 loop we replace Eqs. (\[h1phi\])-(\[u1phi\]) with $$\begin{aligned}
  H_1(\phi) &=& \tilde p\left[ (\cos\phi+ \sin\phi)\tau_1 + \sin\phi\tau_2\right]\,,\label{h21phi}\\
U_1({{\bf k}})&=& \frac{1}{\sqrt{2}}\left( \begin{array}{cc}
e^{i \omega} & 1\\
e^{-i \omega}& -1
\end{array}\right)\,,\label{u21phi}\end{aligned}$$ where $\omega=arg(e^{i\phi}+\sin\phi)$. If $H_1 = \tau_\beta H_2\tau_\beta$, then the conclusions are the same as for model-1 loops, with the replacement $\tilde p\rightarrow \tilde p|e^{i\phi}+\sin\phi|$.

We now consider the case where the two loops are related through the reflection operation in Eq. (\[reflection\]). Because the reflection implies $\phi\rightarrow -\phi$, the energy dispersions are different for $H_1$ and $H_2$. In the non-degenerate case ($\delta>0$), $H^{FS}$ now takes the form $$\begin{aligned}
H^{FS}=\left( \begin{array}{cc} \tilde p |e^{i\phi}+\sin\phi|-\delta & 
\Delta_{FS}({{\bf k}}) \\ 
\Delta_{FS}^*({{\bf k}})
& \tilde p|e^{i\phi}-\sin\phi| -\delta  \end{array}\right)\end{aligned}$$ and the resulting spectrum allows gapless excitations, as was the case for model-1 loops.

In the degenerate case, we find it more convenient not to perform the rotation in Eq. (\[UBdGU\]), and diagonalize the original BdG matrix restricted to the subspace $\left(  \hat{\boldsymbol c}_{{{\bf k}},1} \,,     \hat{\boldsymbol c}_{-{{\bf k}}, 2}^\dagger     \right)$, instead. In this subspace, the two diagonal blocks of the BdG matrix, which follow from Eq. (\[BdGmatrix\]), are $ H_1(\phi)$ and $$\begin{aligned}
  -H_2^T(-\phi)&=&- \tau_\beta H_1^T(\phi)\tau_\beta\,,
  \label{H1H2phi}\end{aligned}$$ which follows from (\[h21phi\]) and the reflection operation that relates both loops: $H_2(\phi)=\tau_\beta H_1(-\phi)\tau_\beta$. For $t_\mu=t_1$ (interloop triplet) the BdG reads: $$\begin{aligned}
H_{12}'=\left( \begin{array}{cc}
H_1\phi)  & (d_0 + {\boldsymbol d}\cdot{\boldsymbol \tau})i\tau_2  \\
-i\tau_2 (d_0 + {\boldsymbol d}\cdot{\boldsymbol \tau})  & - \tau_\beta H_1^T(\phi)\tau_\beta 
\end{array}\right)
\label{UBU3}\end{aligned}$$ We identify the TRS cases where the excitation spectrum is fully gapped. For constant $d_0$ and ${\boldsymbol d}=0$, the gapped spectra are obtained for $\beta=1$ and $\beta=3$, respectively: $$\begin{aligned}
\beta=1:  E^2 &=&d_0^2 + {\tilde p}^2 \sin^2\phi + \left(  d_0\pm \tilde p| \sin\phi+\cos\phi |  \right)^2\,,\nonumber
\\ \\
\beta=3: E^2 &=&d_0^2 + {\tilde p}^2 \left[ \sin^2\phi + (\sin\phi+\cos\phi )^2\right]\,.\end{aligned}$$ In the case of interloop singlet $t_\mu=t_2$, we consider $d_0=0$ and constant $\boldsymbol d$. Gapped spectra exist for: $\beta=0$ and nonzero $d_3$; $\beta=1$ and nonzero $d_2$; $\beta=2$ and nonzero $d_1$. All these cases have similar spectra:

$$\begin{aligned}
\beta=0:  E^2 &=&d_3^2 +  {\tilde p}^2 \left[ \sin^2\phi + (\sin\phi+\cos\phi )^2\right] \,,\ \ \\
\beta=1: E^2 &=&d_2^2 +  {\tilde p}^2 \left[ \sin^2\phi + (\sin\phi+\cos\phi )^2\right]\,,\ \ \\
\beta=2: E^2 &=&d_1^2 +  {\tilde p}^2 \left[ \sin^2\phi + (\sin\phi+\cos\phi )^2\right]\,.\ \\end{aligned}$$

Pairing between Dirac loops
---------------------------

Including the spin degree of freedom, we may discuss the pairing between spin degenerate loops $H_1\otimes\sigma_0$ and $H_2\otimes\sigma_0$ described by the BdG matrix: $$\begin{aligned}
H_{12,s}'=\left( \begin{array}{cc}
H_1\otimes\sigma_0 & (d_0 + {\boldsymbol d}\cdot{\boldsymbol \tau})i\tau_2 (t_\mu)_{12}\sigma_s  \\
-i\sigma_s\tau_2 (d_0 + {\boldsymbol d}\cdot{\boldsymbol \tau})(t_\mu)_{21}  & -H_2^T \otimes\sigma_0 
\end{array}\right)\nonumber\\\end{aligned}$$ Whatever the choice for $s=1,2,3$, $H_{12,s}'$ decouples in subblocks for which the results obtained above for Weyl systems can be applied. The (anti)symmetric property of the matrices $t_\mu$, $\sigma_s$ will determine whether the functions $d_0$,$\boldsymbol d$ should be odd of even: if for instance, $s=1,3$ then the parity of $d_0$,$\boldsymbol d$ is as in the Weyl case; if, however, $s=2$ (spin singlet), then the parities should be reversed.

Summary and Conclusions {#concsec}
=======================

We have described broken symmetry phases of nested Weyl and Dirac NLs that are induced by a short range interaction. We made a systematic analysis for two-band Hamiltonians with PT symmetry, where the two nested Weyl NLs can be mapped onto each other through a rotation or reflection operator. Charge and (pseudo)spin density waves always lower the energy and the broken symmetry phase can be metallic, semimetallic or insulating, depending on the operator that maps the the initial NLs onto each other, and on whether they enjoy a local reflection symmetry in the loop plane. This outcome does not depend on whether the initial system is semimetallic or metallic (when the initial FS is composed of two toroidal FSs, one hole- and one electron-like). If the initial system is semimetallic, spontaneous symmetry breaking requires a finite interaction which is attractive for CDWs and repulsive for PSDWs.

We have also studied specific four-band models, including the $\mathbb{Z}_2$ NLs, spin degenerate Dirac systems, and NL’s derived from perturbed spinful Dirac nodal points. The PSDW phases from $\mathbb{Z}_2$ NLs include nodal point and nodal chain semimetals.

Fully gapped superconducting phases from electron pairing in different NLs (interloop pairing), with TRS, have been found. They include all possibilities of triplet and singlet pairing in loop space and spin space.

There has recently been an intensive search for topological semimetal materials. Given that point nodes tend to appear in pairs for symmetry reasons, it is conceivable that suitable engineering can produce double NLs. Indeed, a recent proposal for realizing point nodes (Dirac or Weyl), and pairs of NLs by strain engineering in SnTe and GeTe is relevant here[@LauOrtix]. Another recent proposal concerns layered ferromagnetic rare-earth-metal monohalides [*LnX*]{} ([*Ln*]{}=La, Gd; [*X*]{}=Cl, Br) and a pair of mirror-symmetry protected nodal lines in La[*X*]{} and Gd[*X*]{}[@Nie]. Also, splitting of Dirac rings into pairs of Weyl rings by spin-orbit interaction in InNbS$_2$ has been proposed[@Du]. Two groups of Dirac nodal rings have been experimentally detected[@Lou] in ZrB$_2$. However, the detection of pairs of NLs at the Fermi level is presently still lacking.

We have not addressed the competition between different orderings or interactions, but such an extension of our work might be relevant to real materials. We have also neglected the effect of the long-ranged tail of the Coulomb interaction which could be present if the starting system is a NL semimetal with the screening radius diverging near the Fermi level. In this respect, the study in Ref\[[@roy]\] for a single NL suggests that the critical interaction strength for orderings where a fully gapped spectrum arises could be lowered.

Ordered phases[@gruner], such as orbital and/or spin density waves, can be detected through neutron scattering, or resonant soft x-ray scattering[@xray]. The band structure itself may be studied with angle-resolved photoemission spectroscopy.

Acknowledgments {#acknowledgments .unnumbered}
===============

M.A.N.A. acknowledges partial support from Fundação para a Ciência e Tecnologia (Portugal) through Grant No. UID/CTM/04540/2013, and the hospitality of Computational Science Research Center, Beijing, China, where this work was initiated. M.A.N.A. would like to thank Vítor R. Vieira, Bruno Mera, and Tilen Cadez for a discussion.

Energy spectrum for Hamiltonian (\[rotatedloops\]) {#appa}
==================================================

Taking $\alpha=0$, (CDW case), for instance, and $(a,b)=(1,2)$, the rotated effective Hamiltonian in Eq. (\[rotatedloops\]) then reads: $$A H_{eff}({{\bf k}})  A^\dagger  = 
\left( p_{\parallel} - p_o \right) t_0\tau_1+  p_\perp t_0\tau_2 +U\bar m\ t_1 \tau_{\beta}- \delta\ t_3\tau_0\,.
\label{AHAdelta}$$ We perform a suitable rotation on the Hamiltonian $AH_{eff}A^\dagger$ in Eq. (\[AHAdelta\]) so that its energy spectrum can be written down by inspecting the (anti)commutation relations among its terms. If $\beta=1$ or 2, we introduce a SU(2) rotation in $t_\mu$ space so that in the end, only the matrix $t_3$ appears. The required rotation is $$W=\cos\frac \theta 2 -it_2\tau_{\beta} \sin\frac \theta 2 \,,$$ with the rotation angle $\theta$ given by

$$\begin{aligned}
\sin\theta &=& \frac{U\bar m}{\sqrt{U^2\bar m^2 + \delta^2}}\,,\\
\cos\theta &=& \frac{\delta}{\sqrt{U^2\bar m^2 + \delta^2}}\,,\end{aligned}$$

so that the rotated Hamiltonian for $\beta=1$ reads $$\begin{aligned}
WAH_{eff}A^\dagger W^\dagger &=& 
 \left( p_{\parallel} - p_o \right)\tau_1 + p_\perp\cos\theta\ \tau_2 \nonumber\\
  &+& p_\perp\sin\theta\ t_2\tau_3\nonumber\\
 &-&  \sqrt{U^2\bar m^2 + \delta^2}\ t_3\,,\end{aligned}$$ and the energy spectrum obeys $$\begin{aligned}
E^2&=& \left[ 
\sqrt{  \left( p_{\parallel} - p_o \right)^2 + p_\perp^2\cos^2\theta}\pm
\sqrt{U^2\bar m^2 + \delta^2}
\right]^2 \nonumber\\
 &+&  p_\perp^2\sin^2\theta\,,\end{aligned}$$ which is equivalent to Eq. (\[CDWka\]).

For $\beta=2$ we have $$\begin{aligned}
WAH_{eff}A^\dagger W^\dagger &=&
 \left( p_{\parallel} - p_o \right)\cos\theta\ \tau_1  + 
 p_\perp\tau_2 \nonumber\\ &-&  \left( p_{\parallel} - p_o \right)\sin\theta\  t_2\tau_3\nonumber\\
 &-&  \sqrt{U^2\bar m^2 + \delta^2}\ t_3\,,\end{aligned}$$ and the energy spectrum obeys $$\begin{aligned}
E^2&=& \left[ 
\sqrt{  \left( p_{\parallel} - p_o \right)^2 \cos^2\theta+ p_\perp^2}\pm
\sqrt{U^2\bar m^2 + \delta^2}
\right]^2 \nonumber\\
&+&  \left( p_{\parallel} - p_o \right)^2\sin^2\theta\,,\end{aligned}$$ which is equivalent to Eq. (\[CDWkb\]).

In the case $\beta=3$, it is preferable to perform a SU(2) rotation in $\tau$ space in order to eliminate one of the Pauli matrices $\tau$. To this aim, we introduce $$R_c=\cos\frac \theta 2 -i\tau_3\sin\frac \theta 2 \,,$$ with

$$\begin{aligned}
\sin\theta &=& \frac{ p_{\parallel} - p_o }{\sqrt{\left( p_{\parallel} - p_o \right)^2 + \delta^2}}
\,,\\
\cos\theta &=& \frac{p_\perp}{\sqrt{ \left( p_{\parallel} - p_o \right)^2 + \delta^2}}\,,\end{aligned}$$

so that the rotation of the Hamiltonian now works out as, $$\begin{aligned}
R_cAH_{eff}A^\dagger R_c^\dagger &=&
{\sqrt{ \left( p_{\parallel} - p_o \right)^2 + \delta^2}}\ \tau_2 + U\bar m \ t_1\tau_3-\delta \ t_3\,.\nonumber\\\end{aligned}$$

More generally, if the product $\tau_\alpha\tau_{\beta}=i\epsilon_{\alpha\beta j}\tau_j$, then the above results for the energy still hold, because the appearance of the factor $i$ in the $U\bar m$ term would lead to the replacement $t_1\rightarrow t_2$ in (\[AHAdelta\]), which does not change the (anti)commutation relations among the Hamiltonian terms.

Mean field treatment of PSDW/CDW {#APMF}
================================

Given the order parameter for a PSDW, $ \langle n_s({{\bf r}})\rangle  = \frac 1 2 n + \bar m (-1)^s\cos({{{\bf Q}}}\cdot{{\bf r}})$, one may transform to Fourier space as $$\begin{aligned}
\langle n_s({{\bf r}})\rangle  =\frac 1 N\sum_{{{\bf q}}} \langle  c_{{{\bf q}}s}^\dagger c_{{{\bf q}}s} \rangle +
e^{i{{\bf Q}}\cdot{{\bf r}}}\langle  c_{{{\bf q}}s}^\dagger c_{{{\bf q}}+{{\bf Q}}s} \rangle + 
e^{-i{{\bf Q}}\cdot{{\bf r}}}\langle  c_{{{\bf q}}+{{\bf Q}}s}^\dagger c_{{{\bf q}}s} \rangle\,.\nonumber\\\end{aligned}$$ Using $\hat\psi_s({{\bf r}}) =\sum_{{{\bf q}}} e^{i{{\bf q}}\cdot{{\bf r}}} c_{{{\bf q}}s}/\sqrt{N}$, where $N$ is the number of momentum values in the summation, we see that the above $\langle n_s({{\bf r}})\rangle$ is obtained if $$\begin{aligned}
\frac 1 N\sum_{{{\bf q}}} \langle  c_{{{\bf q}}s}^\dagger c_{{{\bf q}}s} \rangle&=& \frac 1 2 n\,,\\
\frac 1 N\sum_{{{\bf q}}} \langle  c_{{{\bf q}}+{{\bf Q}}s}^\dagger c_{{{\bf q}}s} \rangle =
\frac 1 N\sum_{{{\bf q}}} \langle  c_{{{\bf q}}s}^\dagger c_{{{\bf q}}+{{\bf Q}}s} \rangle &=&
 \frac 1 2  \bar m (-1)^s\,, \label{MFQ}\end{aligned}$$ for $\alpha=3$ (PSDW). If $\alpha=0$ (CDW) then the factor $(-1)^s$ should be omitted. The Hamiltonian is given by $$\begin{aligned}
\hat H_{eff} = \frac 1 2\sum_{{{\bf k}}}
(\hat {\boldsymbol  c}_{{{\bf k}}}^\dagger \hat {\boldsymbol  c}_{{{\bf k}}+{{\bf Q}}}^\dagger) H_{eff}({{\bf k}})
  \left( \begin{array}{c}\hat {\boldsymbol  c}_{{{\bf k}}}\\ \hat {\boldsymbol  c}_{{{\bf k}}+{{\bf Q}}} \end{array}  \right)\end{aligned}$$ where $\hat{\boldsymbol c}_{{{\bf k}}} = (\hat c_{{{\bf k}},1} \hat c_{{{\bf k}},2})$. We assume that $H_{eff}({{\bf k}})$ is diagonalized by a unitary matrix, $S$, so that $SH_{eff}({{\bf k}})S^\dagger$ is the diagonal matrix composed of the eigenenergies. Then, the operators $\hat{\boldsymbol\gamma}$ which destroy the elementary excitations are given by $$\begin{aligned}
\hat{\boldsymbol\gamma} = S 
 \left( \begin{array}{c}\hat {\boldsymbol  c}_{{{\bf k}}}\\ \hat {\boldsymbol  c}_{{{\bf k}}+{{\bf Q}}} \end{array}  \right)\,.
 \label{gammaS}\end{aligned}$$ Following Eq. (\[MFQ\]) we can see that $$\begin{aligned}
\frac 1 N\sum_{{{\bf k}}}  \langle
(\hat {\boldsymbol  c}_{{{\bf k}}}^\dagger \hat {\boldsymbol  c}_{{{\bf k}}+{{\bf Q}}}^\dagger) 
 \left( \begin{array}{cc} 
 0 & \tau_\alpha \\  \tau_\alpha  & 0 \end{array}  \right)
  \left( \begin{array}{c}\hat {\boldsymbol  c}_{{{\bf k}}}\\ \hat {\boldsymbol  c}_{{{\bf k}}+{{\bf Q}}} \end{array}  \right) 
  \rangle
   = \mp 2\bar m \left\{
   \begin{array}{c} \alpha=3\\ \alpha=0\end{array}
  \right.\end{aligned}$$ or, in the eigenbasis using (\[gammaS\]), $$\begin{aligned}
&&\frac 1 N\sum_{{{\bf k}}} \langle \hat{\boldsymbol\gamma}^\dagger S
\left( \begin{array}{cc} 
 0 & \tau_\alpha \\  \tau_\alpha  & 0 \end{array}  \right) S^\dagger\hat{\boldsymbol\gamma} \rangle\ 
 =\mp 2\bar m \left\{
   \begin{array}{c} \alpha=3\\ \alpha=0\end{array}
  \right. \nonumber\\
&=&
\frac 1 N\sum_{{{\bf k}}} 
\sum_j \left[ 
S
\left( \begin{array}{cc} 
 0 & \tau_\alpha \\  \tau_\alpha  & 0 \end{array}  \right) S^\dagger
\right]_{jj} f(E_j({{\bf k}})) 
   \,,
  \label{SalphaS}\end{aligned}$$ where $f(x)=1/(1+e^{x/T})$ denotes the Fermi-Dirac distribution function, and $j=1,...,4$ denotes a band index. The energy dispersions, $E_j({{\bf k}})$, are given in the main text. However, it is more convenient to work with the transformed Hamiltonian $A H_{eff}A^\dagger$ as in the main text, which implies that all operators are similarly rotated and $S\rightarrow SA^\dagger$ above. Then, Eq. (\[SalphaS\]) can be written in the form, $$\begin{aligned}
&&\frac 1 N\sum_{{{\bf k}}} 
\sum_j \left[ 
SA^\dagger
\left( \begin{array}{cc} 
 0 & \tau_\alpha\tau_{\beta} \\  \tau_{\beta}\tau_\alpha  & 0 \end{array}  \right)
AS^\dagger
\right]_{jj} f(E_j({{\bf k}}))\nonumber\\
&=& \mp 2\bar m \left\{
   \begin{array}{c} \alpha=3\\ \alpha=0\end{array}
  \right. \,.\nonumber\\
  \label{APMFeq}\end{aligned}$$

Critical interaction $U_{cr}$ for degenerate NLs {#MFloop}
================================================

We consider a circular nodal line and use the momentum parametrization in Eq. (\[pcoordinates\]).

We linearize the theory in a toroidal region surrounding the NL up to a momentum cutoff: $0<\tilde p<\tilde p_c$, $0<\theta,\phi<2\pi$. The volume element is $d^3p=(p_0 + \tilde p\cos\phi)\tilde p\cdot d\tilde p  d\theta d\phi$. The number of ${{\bf k}}$ terms in the toroidal region around the NL is then given by $$\begin{aligned}
N=\frac{1}{(2\pi\hbar)^3}\int d^3p &=&   \frac{p_0\tilde p_c^2}{4\pi\hbar^3}\,.\label{N}\end{aligned}$$

In order to simplify the calculations, it is assumed that the dispersion relation has the same velocity, $v$, in the NL plane and perpendicular to it. Then, using $p_\parallel - p_0 = \tilde p\cos\phi$, the model-1 NL Eq. (\[loop1\].a) reads $$\begin{aligned}
H_0({{\bf k}}) = v\tilde p \left( \cos\phi\tau_1 + \sin\phi\tau_2\right)\,,\end{aligned}$$ where $v=v_1=v_2$ and we shall take $\delta=0$. As before, we proceed considering the velocity $v=1$ and but shall restore it in the final result for $U_{cr}$. In the cases where $\tau_\alpha \tau_{\beta}\propto\tau_1$ and $\tau_\alpha\tau_{\beta}\propto\tau_2$ the ordered phases are semimetalic and yield similar mean field equations. In the case where $\tau_\alpha\tau_{\beta}\propto\tau_2$, for instance, the negative energy bands are $$\begin{aligned}
E_{1(2)} &=&-\sqrt{ \tilde p^2\pm 2\tilde pU\bar m\sin\phi + U^2\bar m^2}\end{aligned}$$ and the l.h.s. of Eq. (\[APMFeq\]) takes the form, $$\begin{aligned}
&&\frac 1 N\sum_{{{\bf k}}} 
\sum_j \left[ 
SA^\dagger
\left( \begin{array}{cc} 
 0 & \tau_\alpha\tau_{\beta} \\  \tau_{\beta}\tau_\alpha  & 0 \end{array}  \right)
AS^\dagger
\right]_{jj} f(E_j({{\bf k}})) \nonumber\\ &=&  \frac 1 N\
\int \frac{2(\tilde p\sin\phi-U\bar m) (p_0 + \tilde p\cos\phi)
d^3p }{(2\pi\hbar)^3\sqrt{ \tilde p^2- 2\tilde pU\bar m\sin\phi + U^2\bar m^2}}\,.\end{aligned}$$ In the limit $U\bar m\rightarrow 0$ one can Taylor expand the integrand to first order. The mean field equations (\[APMFeq\]) yield $$\begin{aligned}
\frac{p_0\tilde p_c}{4\pi N\hbar^3} =  \pm \frac{1}{U_{cr}}
\Longrightarrow\
U_{cr}=\pm v \tilde p_c
 \left\{
   \begin{array}{c} \alpha=3\\ \alpha=0\end{array}
  \right. \,,\end{aligned}$$ where we used (\[N\]) and the velocity $v$ has been restored. The finiteness of $U_{cr}$ stems from the linear form of the density of states near the Fermi level. The case where $\tau_\alpha\tau_{\beta}\propto\tau_{c\neq a,b}$, where the PSDW phase is insulating, yields a similar result modified by a prefactor of $1/2$: $U_{cr}=\pm v \tilde p_c/2$. This is valid also for the case $\tau_\alpha\tau_{\beta}=1$, where the density wave phase is metallic. We see then that the CDW phase requires an attractive interaction.

[99]{}

K. Mullen, B. Uchoa, and D. T. Glatzhofer, Phys. Rev. Lett. [**115**]{}, 026403 (2015).

J. L. Lu, X. Y. Li, S. Q. Yang, J. X. Cao, X. G. Gong, and H. J. Xiang, Chin. Phys. Lett. [**34**]{}, 057302 (2017).

A. A. Burkov, M. D. Hook, and L. Balents, Phys. Rev. B [**84**]{}, 235126 (2011).

S. A. Yang, H. Pan, and F. Zhang, Phys. Rev. Lett. [**113**]{}, 046401 (2014).

Y. Chen, Y. Xie, S. A. Yang, H. Pan, F. Zhang, M. L. Cohen, and S. Zhang, Nano Letters [**15**]{}, 6974 (2015).

H. Weng, Y. Liang, Q. Xu, R. Yu, Z. Fang, X. Dai, and Y. Kawazoe, Phys. Rev. B [**92**]{}, 045108 (2015).

D.-W. Zhang, Y. X. Zhao, R.-B. Liu, Z.-Y. Xue, S.-L. Zhu, and Z. D. Wang, Phys. Rev. A [**93**]{}, 043617 (2016).

C. Zhong, Y. Chen, Z.-M. Yu, Y. Xie, H. Wang, S. A. Yang, and S. Zhang, Nat. Comm. [**8**]{}, 15641 (2017).

M. Ezawa, Phys. Rev. B [**96**]{}, 041202 (2017).

Z. Yan, R. Bi, H. Shen, L. Lu, S.-C. Zhang, and Z. Wang, Phys. Rev. B [**96**]{}, 041103 (2017).

W. Chen, H.-Z. Lu, and J.-M. Hou, Phys. Rev. B [**96**]{}, 041102 (2017).

Y. Zhou, F. Xiong, X. Wan, and J. An, Phys. Rev. B [**97**]{}, 155140 (2018).

L. Li, C. H. Lee, and J. Gong, Phys. Rev. Lett. [**121**]{}, 036401 (2018).

Y. Kim, B. J. Wieder, C. L. Kane, and A. M. Rappe, Phys. Rev. Lett. [**115**]{}, 036806 (2015).

R. Yu, H. Weng, Z. Fang, X. Dai, and X. Hu, Phys. Rev. Lett. [**115**]{}, 036807 (2015).

A. Narayan, Phys. Rev. B [**94**]{}, 041409 (2016).

K. Taguchi, D.-H. Xu, A. Yamakage, and K. T. Law, Phys. Rev. B [**94**]{}, 155206 (2016).

C.-K. Chan, Y.-T. Oh, J. H. Han, and P. A. Lee, Phys. Rev. B [**94**]{}, 121106 (2016).

X.-X. Zhang, T. T. Ong, and N. Nagaosa, Phys. Rev. B [**94**]{}, 235137 (2016).

Z. Yan and Z. Wang, Phys. Rev. Lett. [**117**]{}, 087402 (2016).

Z. Yan and Z. Wang, Phys. Rev. B [**96**]{}, 041206 (2017).

L. Li, C. Yin, S. Chen, and M. A. N. Araújo, Phys. Rev. B [**95**]{}, 121107 (2017).

L. Li, H. H. Yap, M. A. N. Araújo, and J. Gong, Phys. Rev. B [**96**]{}, 235424 (2017).

S. Sur and R. Nandkishore, New J. Phys. [**18**]{} 115006 (2016).

R. Nandkishore, Phys. Rev B [**93**]{}, 020506 (2016).

H. Shapourian, Y. Wang, and S. Ryu, Phys. Rev B [**97**]{}, 094508 (2018).

B. Roy, Phys. Rev B [**96**]{}, 041113 (2017).

Y. Wang and P. Ye Phys. Rev B [**94**]{}, 075115 (2016).

P.-Y. Chang, and C.-H. Yee, Phys. Rev. B [**96**]{}, 081114 (2017).

C. Fang, Y. Chen, H.-Y. Kee, and L. Fu, Phys. Rev. B [**92**]{}, 081201 (2015).

G. Bian et. at., Phys. Rev B [**93**]{}, 121113 (2016).

George Grüner, [*Density waves in solids* ]{} (Addison-Wesley Publishing Company, 1994).

A. K. Mitchell and L. Fritz, Phys. Rev B [**93**]{}, 035137 (2016).

J. Wang, Phys. Rev B [**96**]{}, 081107 (2017).

P. D. Sacramento, M. A. N. Araújo, and E. V. Castro, Europhys. Lett. [**105**]{}, 37011 (2014).

P. D. Sacramento, M. A. N. Araújo, V. R. Vieira, V. K. Dugaev, and J. Barnaś, Phys. Rev B [**85**]{}, 014518 (2012).

B. Béri, Phys. Rev B [**81**]{}, 134515 (2010).

Y. Bang, and H.-Y. Choi, Phys. Rev B [**78**]{}, 134523 (2008).

P. Fulde, and A. Ferrel, Phys. Rev. [**135**]{}, A550 (1964).

A. I. Larkin, and Y. N. Ovchinnikov, Sov. Phys. JETP-USSR [**20**]{}, 762 (1965).

H. Shimahara, Phys. Rev B [**50**]{}, 12760 (1994).

Y. Wang, and R. M. Nandkishore, Phys. Rev. [**95**]{}, 060506 (2017).

A. Lau, and C. Ortix, arXiv:1804.09574.

S. Nie, H. Weng, and F. B. Prinz, arXiv:1803.08486.

Y. Du, X. Bo, D. Wang, E.-J. Kan, C.-G. Duan, S. Y. Savrasov, and X. Wan, Phys. Rev. B [**96**]{}, 235152 (2017).

R. Lou, P. Guo, M. Li, Q. Wang, Z. Liu, S. Sun, C. Li, X. Wu, Z. Wang, Z. Sun, D. Shen, Y. Huang, K. Liu, Z-Y Lu, H. Lei, H. Ding, and S. Wang, npj Quantum Materials [**3**]{}, 43 (2018).

B. Pan, H. Jang, J.-S. Lee, R. Sutarto, F. He, J. F. Zeng, Y. Liu, X. Zhang, Y. Feng, Y. Hao, J. Zhao, H. C. Xu, Z. H. Chen, J. Hu, and D. Feng, arXiv:1808.08562.
