# Thermal conductivity of nonunitary triplet superconductors: application to
UTe2

Vivek Mishra  Department of Physics, University of Florida, Gainesville,
Florida 32611  Ge Wang  Department of Physics, University of Florida,
Gainesville, Florida 32611  P. J. Hirschfeld  Department of Physics,
University of Florida, Gainesville, Florida 32611

###### Abstract

There is considerable evidence that the heavy fermion material UTe2 is a spin-
triplet superconductor, possibly manifesting time-reversal symmetry breaking,
as measured by Kerr effect below the critical temperature, in some samples.
Such signals can arise due to a chiral orbital state, or possible nonunitary
pairing. Although experiments at low temperatures appear to be consistent with
point nodes in the spectral gap, the detailed form of the order parameter and
even the nodal positions are not yet determined. Thermal conductivity
measurements can extend to quite low temperatures, and by varying the heat
current direction should be able to provide information on the order parameter
structure. Here we derive a general expression for the thermal conductivity of
a spin triplet superconductor, and use it to compare the low-temperature
behavior of various states proposed for UTe2.

##  I Introduction

The uranium-based superconductor $\mathrm{UTe}_{2}$ has stimulated a large
number of experimental and theoretical studies, initially because of its
apparent role as a paramagnetic end-point of a family of ferromagnetic
superconductors(Ran  _et al._ , 2019a; Aoki  _et al._ , 2019a, 2022a), and
later as evidence for spin-triplet superconductivity accumulated. The nuclear
magnetic resonance (NMR) Knight shift measurements on the earlier samples did
not show any change in the superconducting stateFujibayashi  _et al._ (2022),
although recent Knight shift measurements on high quality samples show a small
reduction along all three axesMatsumura  _et al._ (2023). Both measurements
support spin triplet pairing, however, the spin structure of Cooper pairs
remains unclear. Another piece of evidence that indicate spin triplet pairing
is the size of the upper critical field $H_{c2}$ that exceeds the Pauli limit
for all field directionsRan  _et al._ (2019a); Aoki  _et al._ (2019b).
Measured power law temperature dependence in NMR relaxation, specific heat(Ran
_et al._ , 2019a) and thermal conductivity(Metz  _et al._ , 2019; Hayes  _et
al._ , 2024) were found to be consistent with point nodes, expected for a
triplet superconductor in a system with strong spin-orbit coupling(Volovik and
Gorkov, 1985; Ueda and Rice, 1985; Blount, 1985). Finally, a reentrant
superconducting phase was shown to be stabilized at high magnetic fields(Ran
_et al._ , 2019b).

A second set of measurements relevant to the nature of the superconducting
state purports to exhibit evidence for time reversal symmetry breaking (TRSB)
below $T_{c}$, suggesting that UTe2 may support the long-sought chiral
$p$-wave state that may serve as a quantum computing platform(Wang  _et al._ ,
2020; Kitaev, 2001; Wilczek, 2009; Ivanov, 2001; Chung  _et al._ , 2011).
Initially, polar Kerr effect(Hayes  _et al._ , 2021) experiments suggested
TRSB occurred in the superconductor, implying the existence of a
multicomponent spin triplet order parameter. According to group theoretical
classifications of the 1-dimensional irreducible representations (irreps)
allowed in orthorhombic symmetry, order parameters corresponding to single
irreps must be unitary triplet states, meaning that any TRSB must arise from a
nonunitary multicomponent state. Such combinations of 1D representations were
discussed intensively, particularly because a double specific heat transition
was sometimes seen in early samples, recalling the specific heat experiments
in multicomponent UPt3.

More recently, measurements of a new generation of high quality UTe2 crystals
grown in molten salt flux have challenged this characterization of UTe2 as a
chiral triplet state breaking time reversal symmetry. The Kerr effect was
observed in a sample showing two specific heat jumps, but as the quality of
the samples improved, only a single transition was observed(Rosa  _et al._ ,
2022; Sakai  _et al._ , 2022). A recent Kerr effect investigation of both the
old and new generation $\mathrm{UTe}_{2}$ single crystals displaying a single
specific heat jump found no evidence for TRSB superconductivity(Ajeesh  _et
al._ , 2023). Similarly, muon spin relaxation ($\mu$SR) measurements of the
molten-salt flux grown samples found no evidence of TRSB(Azari  _et al._ ,
2023). Finally, sound velocity changes across $T_{c}$(Theuss  _et al._ ,
2024), as well as recent NMR Knight shift measurements on similar
samples(Matsumura  _et al._ , 2023) both point to a single-component, odd
parity order parameter, i.e. are inconsistent with the previous hypothesis of
nonunitary pairing.

The thermal conductivity $\kappa(T)$ is an important probe of the gap
structure of unconventional superconductors, reflecting the ability of the
superconductor to carry heat current in various directions. The theory of
thermal conductivity in unitary triplet superconductors is quite similar to
the well-known theory developed for singlet superconductors(Ambegaokar and
Griffin, 1965; Kadanoff and Martin, 1963). Most of the popular model triplet
states in the literature, including the 3He-A phase, are in this class. In
that case, the triplet quasiparticle energies are $E_{\bf k}=\sqrt{\xi_{\bf
k}^{2}+|{\bf d}({\bf k})|^{2}}$, where ${\bf d}({\bf k})$ is the triplet order
parameter vector defining its structure in spin space via
$\Delta_{\sigma\sigma^{\prime}}=\left[{\bf d}({\bf
k})\cdot\vec{\sigma}(i\sigma_{2})\right]_{\sigma\sigma^{\prime}}$. Since the
thermal current response depends only on the quasiparticle energies, and the
same expressions can be used for triplet superconductors(Schmitt-Rink  _et
al._ , 1986; Hirschfeld  _et al._ , 1986, 1988) with $|\Delta({\bf k})|^{2}$
replaced by $|{\bf d}({\bf k})|^{2}$. As we show below, however, in the
nonunitary case, additional terms involving the spin moment ${\bf q}\equiv
i{\bf d}({\bf k})\times{\bf d}^{*}({\bf k})$ carried by quasiparticles of
momentum $\bf k$ occur in both the quasiparticle energies and the weights of
scattering processes. Furthermore, in nonunitary triplet superconductors, the
zeros of $|{\bf d}({\bf k})|^{2}$ differ from those of the zeros of $E_{\bf
k}$ even when ${\bf k}$ is on the Fermi surface (“spectral nodes”). This
distinction may be important; it was suggested by Ishihara et al.(Ishihara
_et al._ , 2023) that in UTe2, complex linear combinations of 1D irreducible
representations could support spectral nodes pointing in generic directions in
the orthorhombic Brillouin zone and thereby explain early experiments
exhibiting TRSB and relative isotropy of the low-temperature penetration
depth. On the other hand, order parameters corresponding to a single 1D irrep
must be real, with nodes along high symmetry axes.

In this paper, we derive a general form of the thermal conductivity of a
triplet superconductor in the presence of nonmagnetic disorder, and evaluate
it for various types of triplet states that have been proposed for UTe2. The
aim is to see if there are qualitative distinctions between the thermal
conductivity temperature and heat current direction dependence of unitary and
nonunitary states, and whether or not predictions of low-temperature behavior
can be used, by comparison with experiments, to identify the ground state of
UTe2.

##  II Model & Formalism

$\Gamma$ | Gap function ${\bf d}({\bf k})$ | Nodes
---|---|---
$A_{1u}$ | ($p_{1}k_{x},p_{2}k_{y},p_{3}k_{z}$) | accidental
$B_{1u}$ | ($p_{1}k_{y},p_{2}k_{x},p_{3}k_{x}k_{y}k_{z}$) |  $z$ axis
$B_{2u}$ | ($p_{1}k_{z},p_{2}k_{x}k_{y}k_{z},p_{3}k_{x}$) |  $y$ axis
$B_{3u}$ | ($p_{1}k_{x}k_{y}k_{z},p_{2}k_{z},p_{3}k_{y}$) |  $x$ axis
Table 1: List of possible spin-triplet superconducting states for an orthorhombic crystal with strong spin orbit coupling. Here $p_{i=1,..,3}$ are constants and $\forall p_{i}\in\mathbf{R}$. Classification | ${\bf d}({\bf k})$ | Cylindrical FS | Spherical FS
---|---|---|---
|  | $r$ | Nnodes | $r$ | Nnodes
Antiferromagnetic non-unitary | $A_{1u}+irB_{1u}$ | $r<1$ | 0 | $r<1$ | 0
$r=1$ | 4 ($xy$ plane) | $r=1$ | 4 ($xy$ plane)
$1<r\leq\sqrt{2}$ | 8 ($xz$, $yz$ planes) | $r>1$ | 8 ($xz$, $yz$ planes)
$r>\sqrt{2}$ | 0 |  |
$A_{1u}+irB_{2u}$ | $r<1$ | 0 | $r<1$ | 0
$r=1$ | 2 ($x$ axis) | $r=1$ | 4 ($xz$ plane)
$1<r\leq\sqrt{2}$ | 4 ($xy$ plane) | $r>1$ | 8 ($xy$, $yz$ planes)
$r>\sqrt{2}$ | 8 ($xy$, $yz$ planes) |  |
$A_{1u}+irB_{3u}$ | $r<1$ | 0 | $r<1$ | 0
$r=1$ | 2 ($y$ axis) | $r=1$ | 4 ($yz$ plane)
$1<r\leq\sqrt{2}$ | 4 ($xy$ plane) | $r>1$ | 8 ($xy$, $xz$ planes)
$r>\sqrt{2}$ | 8 ($xy$, $xz$ planes) |  |
Ferromagnetic non-unitary |  $B_{1u}+irB_{2u}$ ( $\hat{x}$ spin moment ) | $r<\frac{1}{\sqrt{2}}$ | 0 | $r<1$ | 4 ($xz$ plane)
$\frac{1}{\sqrt{2}}\leq r<1$ | 4 ($xz$ plane) | $r=1$ | 2 ($x$ axis)
$r=1$ | 2 $x$ axis | $r>1$ | 4 ($xy$ plane)
$r>1$ | 4 ($xy$ plane) |  |
$B_{1u}+irB_{3u}$ ( $\hat{y}$ spin moment ) | $r<\frac{1}{\sqrt{2}}$ | 0 | $r<1$ | 4 ($yz$ plane)
$\frac{1}{\sqrt{2}}\leq r<1$ | 4 ($yz$ plane) | $r=1$ | 2 ($y$ axis)
$r=1$ | 2 ($y$ axis) | $r>1$ | 4 ($xy$ plane)
$r>1$ | 4 ($xy$ plane) |  |
$B_{2u}+irB_{3u}$ ($\hat{z}$ spin moment ) | $r\leq\frac{1}{\sqrt{2}}$ | 4 ($yz$ plane) | $r<1$ | 4 ($yz$ plane)
$\frac{1}{\sqrt{2}}<r<\sqrt{2}$ | 0 | $r=1$ | 2 ($z$ axis)
$r\geq\sqrt{2}$ | 4 ($xz$ plane) | $r>1$ | 4 ($xz$ plane)
Table 2: Six possible mixed IR states

###  II.1 Superconducting State

In a general triplet superconductor, the spin-structure of the superconducting
order parameter is constrained by the underlying crystal symmetries. The
structure of the UTe2 crystals corresponds to the orthorhombic point group
$D_{2h}$, and the symmetry of the odd-parity pairing states can be deduced
according to the irreducible representations of the $D_{2h}$ point
group(Marchenko, 1986; Sigrist and Rice, 1987; Annett, 1990). Table 1 shows
the odd parity triplet superconducting states that we consider in this
article. Here, we do not consider the even parity states for the $D_{2h}$
point group, and we further assume strong spin-orbit coupling (SOC) due to
heavy atoms like U and Te. In the weak SOC limit, the odd-parity states for
the $D_{2h}$ point group come with line nodes that are not consistent with
experimental measurements.

The ${\bf d}$-vector is real for the superconducting states listed in Table 1,
hence these states preserve time reversal symmetry (TRS). These states are
unitary triplet states, i.e.
$\hat{\Delta}^{\dagger}\hat{\Delta}\propto\hat{\mathbb{1}}$. We denote a
2$\times$2 matrix in the spin space with $\hat{\square}$ and a 4$\times$4
matrix in the Nambu-spin space with $\check{\square}$. A TRSB state is not
possible with a single component order parameter, noting that the $D_{2h}$
group has only one dimensional irreducible representations. We construct the
TRSB superconducting state with a combination of two irreducible
representations, and all such possible combinations are shown in Table 2. In
principle, a combination of two different irreducible representations involves
six real constants, however, we introduce a single parameter model for the
TRSB states. The effective ${\bf d}$-vector is $({\bf d}_{1}+ir{\bf
d}_{2})/\sqrt{1+r^{2}}$, where $r$ is the mixing parameter that determines the
relative strength of the individual order parameter. The individual ${\bf
d}$-vectors are listed in Table 1 with all the coefficients $p_{i=1,2,3}$ are
set to unity.

In the Nambu-spin basis, the mean field Hamiltonian reads,

$\check{\mathbf{H}}=\begin{pmatrix}\xi_{\bf
k}\sigma_{0}&i\Delta_{0}(\mathbf{d}\cdot\boldsymbol{\sigma})\sigma_{y}\\\
-i\Delta_{0}\sigma_{y}\left(\mathbf{d}^{*}\cdot\boldsymbol{\sigma}\right)&-\xi_{\bf
k}\sigma_{0}\end{pmatrix},$ (1)

where $\xi_{\bf k}$ is the electronic dispersion, $\Delta_{0}$ is the
superconducting gap energy scale, $\boldsymbol{\sigma}$ is the Pauli vector in
the spin space spanned by the Pauli matrices, and $\sigma_{0}$ is the identity
matrix in the spin space. We adopt a model where the electronic dispersion
reads,

$\displaystyle\xi_{\bf
k}=\frac{\hbar^{2}k_{x}^{2}}{2m_{a}}+\frac{\hbar^{2}k_{y}^{2}}{2m_{b}}-\mu-2t_{\perp}\cos
k_{z},$ (2)

where $m_{a/b}$ is the effective masses along $\hat{x}/\hat{y}$ directions,
$t_{\perp}$ is the hopping energy that controls the $\hat{z}$ velocity, and
$\mu$ is the chemical potential. We further assume $t_{\perp}\ll\mu$. An
alternate dispersion with closed Fermi surface has been considered in the
appendix. The bare Green’s function is,

$\check{G}_{0}=(\check{\mathbb{1}}\omega-\check{\mathbf{H}})^{-1}=\begin{pmatrix}\hat{G}_{11}&\hat{G}_{12}\\\
\hat{G}_{21}&\hat{G}_{22}\end{pmatrix}.$ (3)

Here, $\omega$ is the quasi-particle energy. The Matsubara Green’s function
can be obtained by $\omega\rightarrow i\omega_{n}$. The 2$\times$2 matrices in
the spin space are,

$\displaystyle\hat{G}_{11}$ $\displaystyle=$
$\displaystyle\frac{(\omega+\xi)}{{D}}\left[(\omega^{2}-\xi^{2}-\Delta_{0}^{2}|{\bf
d}|^{2})\sigma_{0}+\Delta_{0}^{2}\mathbf{q}\cdot\boldsymbol{\sigma}\right],$
(4) $\displaystyle\hat{G}_{12}$ $\displaystyle=$
$\displaystyle\left[(\omega^{2}-\xi^{2}-\Delta_{0}^{2}|{\bf
d}|^{2})\sigma_{0}+\Delta_{0}^{2}\mathbf{q}\cdot\boldsymbol{\sigma}\right]\frac{i\Delta_{0}({\bf
d}\cdot\boldsymbol{\sigma})\sigma_{y}}{{D}},$ (5) $\displaystyle\hat{G}_{21}$
$\displaystyle=$ $\displaystyle-\left[(\omega^{2}-\xi^{2}-\Delta_{0}^{2}|{\bf
d}|^{2})\sigma_{0}+\Delta_{0}^{2}\mathbf{q}\cdot\boldsymbol{\sigma}^{T}\right]\frac{\Delta_{0}i\sigma_{y}({\bf
d}^{*}\cdot\boldsymbol{\sigma})}{{D}},$ (6) $\displaystyle\hat{G}_{22}$
$\displaystyle=$
$\displaystyle\frac{(\omega-\xi)}{{D}}\left[(\omega^{2}-\xi^{2}-\Delta_{0}^{2}|{\bf
d}|^{2})\sigma_{0}+\Delta_{0}^{2}\mathbf{q}\cdot\boldsymbol{\sigma}^{T}\right].$
(7)

Here ${\bf q}=i({\bf d}\times{\bf d}^{*})$ and the denominator ${D}$ is,

$\displaystyle{D}$ $\displaystyle=$
$\displaystyle\left(\omega^{2}-\xi^{2}-\Delta_{0}^{2}|{\bf
d}|^{2}\right)^{2}-\Delta_{0}^{4}|\mathbf{q}|^{2}=(\xi^{2}-\omega^{2}+\Delta^{2}_{+})(\xi^{2}-\omega^{2}+\Delta^{2}_{-}).$
(8)

Here we introduce $\Delta^{2}_{\pm}=\Delta_{0}^{2}(|{\bf d}|^{2}\pm|{\bf
q}|)$. For the unitary case, ${\bf q}=0$, therefore, there is only a single
energy scale. In contrast, for the TRSB non-unitary states, ${\bf q}\neq 0$,
which leads to non-degenerate excitation energies. For single component order
parameters, ${\bf q}$ vanishes. However, for a mixture of multiple irreducible
representations ${\bf q}$ remains finite, and it can be interpreted as spin
moment of the Cooper pairs. The average of ${\bf q}$ over the Fermi surface
may or may not vanish. The non-unitary states can therefore be further divided
in anti-ferromagnetic (AF) and ferromagnetic (FM) states, where the average of
${\bf q}$ vanishes over the Fermi surface for the former, while remains finite
for the latter(Volovik and Gorkov, 1985). The six possible nonunitary states
are shown in Table 2, with the possibility of nodes on a spherical or a
cylindrical Fermi surface open along the $\hat{z}$ axis. For cylindrical Fermi
surface, we adopt cylindrical coordinates with $k_{z}$ dependence of the gap
functions replaced with $\sin(k_{z}d/2)$, where $d$ is the $z$-axis. The
factor of half is added to ensure only a single pair of point nodes in the
first Brillouin zone in the unitary limit. However, this does not have any
qualitative effect on our results.

###  II.2 Impurity Scattering & Thermal Transport

In order to calculate the thermal conductivity, we need to include the effect
of impurity scattering that dominates all other relaxation mechanism at low
temperatures. We consider elastic impurity scattering due to pointlike
defects, and include its effect through a disorder-averaged self-energy. The
impurity self-energy is calculated within the self-consistent T-matrix
approximation. The momentum integrated Green’s function ${\hat{g}}$ has both
vector and scalar components for the normal Green’s functions,

$\displaystyle\hat{g}_{11}$ $\displaystyle=$ $\displaystyle\pi
N_{0}\left(g_{0}+\mathbf{g}\cdot\boldsymbol{\sigma}\right),$ (9)
$\displaystyle\hat{g}_{22}$ $\displaystyle=$ $\displaystyle\pi
N_{0}\left(g_{0}+\mathbf{g}\cdot\boldsymbol{\sigma}^{T}\right),$ (10)

and the anomalous Green’s functions $\hat{g}_{12}$ and $\hat{g}_{21}$ vanish
because the odd-parity order parameter averages to zero. In Eq. (9) and (10),
$\mathbf{g}$ is directly related to the Fermi surface average of the spin
moment ${\bf q}$ and remains finite for the chiral non-unitary states, only.
Using these integrated Green’s function, we can write the T-matrix self-energy
for the non-magnetic impurities as,

$\displaystyle\check{\Sigma}$ $\displaystyle=$
$\displaystyle\check{\tau}_{3}\cdot\begin{pmatrix}\Sigma_{3}+\mathbf{\Sigma}_{3}\cdot\boldsymbol{\sigma}&0\\\
0&\Sigma_{3}+\mathbf{\Sigma}_{3}\cdot\boldsymbol{\sigma}^{T}\end{pmatrix}+\check{\tau}_{0}\cdot\begin{pmatrix}\Sigma_{0}+\mathbf{\Sigma}\cdot\boldsymbol{\sigma}&0\\\
0&\Sigma_{0}+\mathbf{\Sigma}\cdot\boldsymbol{\sigma}^{T}\end{pmatrix},$ (11)
$\displaystyle\Sigma_{3}$ $\displaystyle=$
$\displaystyle\Gamma_{u}\frac{\cot\delta_{s}\left[\cot^{2}\delta_{s}-\left(g_{0}^{2}+\mathbf{g}\cdot\mathbf{g}\right)\right]}{{D}_{imp}},$
(12) $\displaystyle\mathbf{\Sigma}_{3}$ $\displaystyle=$
$\displaystyle\Gamma_{u}\frac{2\cot\delta_{s}g_{0}\mathbf{g}}{{D}_{imp}},$
(13) $\displaystyle\Sigma_{0}$ $\displaystyle=$
$\displaystyle\Gamma_{u}\frac{g_{0}\left[\cot^{2}\delta_{s}-\left(g_{0}^{2}-\mathbf{g}\cdot\mathbf{g}\right)\right]}{{D}_{imp}},$
(14) $\displaystyle\mathbf{\Sigma}$ $\displaystyle=$
$\displaystyle\Gamma_{u}\frac{\mathbf{g}\left[\cot^{2}\delta_{s}+\left(g_{0}^{2}-\mathbf{g}\cdot\mathbf{g}\right)\right]}{{D}_{imp}},$
(15) $\displaystyle{D}_{imp}$ $\displaystyle=$
$\displaystyle\cot^{4}\delta_{s}-2\cot^{2}\delta_{s}\left(g_{0}^{2}+\mathbf{g}\cdot\mathbf{g}\right)+\left(g_{0}^{2}-\mathbf{g}\cdot\mathbf{g}\right)^{2},$
(16)

where $\delta_{s}\equiv\tan^{-1}(\pi N_{0}V_{imp})$ is the $s$-wave scattering
phase shift, $\Gamma_{u}=n_{imp}/(\pi N_{0})$, $n_{imp}$ and $V_{imp}$ are the
impurity concentration and impurity potential, respectively. . Here
$\check{\tau}_{3}$ is the Pauli matrix in the Nambu space. The $\tau_{3}$
component of the self-energy that renormalizes the electronic dispersion is
ignored. It can be absorbed in the chemical potential. The impurity dressed
Green’s function reads,

$\displaystyle\check{\mathbf{G}}^{-1}$ $\displaystyle=$
$\displaystyle\check{G}_{0}^{-1}-\check{\Sigma}=\begin{pmatrix}\tilde{\omega}-\xi\sigma_{0}-\mathbf{\Sigma}\cdot\boldsymbol{\sigma}&-\hat{\Delta}\\\
-\hat{\Delta}^{\dagger}&\tilde{\omega}+\xi\sigma_{0}-\mathbf{\Sigma}\cdot\boldsymbol{\sigma}^{T}\end{pmatrix},$
$\displaystyle\check{\mathbf{G}}$ $\displaystyle=$
$\displaystyle\begin{pmatrix}\hat{\mathbf{G}}_{11}&\hat{\mathbf{G}}_{12}\\\
\hat{\mathbf{G}}_{21}&\hat{\mathbf{G}}_{22}\end{pmatrix}.$ (17)

Here impurity renormalized $\tilde{\omega}=\omega-\Sigma_{0}$, which is
obtained self-consistently. Unlike unitary superconductors, the impurity
dressed non-unitary Green’s function acquires a different structure than the
bare Green’s function, in particular the structure in spin space for the
normal component. The individual component of the Green’s function
$\check{\mathbf{G}}$ are,

$\displaystyle\hat{\mathbf{G}}_{11}$ $\displaystyle=$
$\displaystyle\frac{L_{0}+\mathbf{L}_{1}\cdot\boldsymbol{\sigma}}{\widetilde{\mathbf{D}}},$
(18) $\displaystyle\hat{\mathbf{G}}_{22}$ $\displaystyle=$
$\displaystyle\frac{L_{0}(\xi\rightarrow-\xi)+\mathbf{L}_{1}(\xi\rightarrow-\xi)\cdot\boldsymbol{\sigma}^{T}}{\widetilde{\mathbf{D}}},$
(19) $\displaystyle\hat{\mathbf{G}}_{12}$ $\displaystyle=$
$\displaystyle\left[2\xi(\boldsymbol{\Sigma}\cdot{\bf
d})+(\tilde{\omega}^{2}-\xi^{2}-\Delta_{0}^{2}|{\bf
d}|^{2}+\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}){\bf
d}\cdot\boldsymbol{\sigma}+i\Delta_{0}^{2}({\bf q}\times{\bf
d})\cdot\boldsymbol{\sigma}-2(\boldsymbol{\Sigma}\cdot{\bf
d})\boldsymbol{\Sigma}\cdot\boldsymbol{\sigma}\right.$ (20)
$\displaystyle\left.+2i\tilde{\omega}(\boldsymbol{\Sigma}\times{\bf
d})\cdot\boldsymbol{\sigma}\right]\frac{i\sigma_{y}\Delta_{0}}{\widetilde{\mathbf{D}}},$
$\displaystyle\hat{\mathbf{G}}_{21}$ $\displaystyle=$
$\displaystyle\left[-2\xi\boldsymbol{\Sigma}\cdot{\bf
d}^{\ast}+(\tilde{\omega}^{2}-\xi^{2}-\Delta_{0}^{2}|{\bf
d}|^{2}+\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}){\bf
d}^{\ast}\cdot\boldsymbol{\sigma}^{T}-i\Delta_{0}^{2}({\bf q}\times{\bf
d}^{\ast})\cdot\boldsymbol{\sigma}^{T}-2\boldsymbol{\Sigma}\cdot{\bf
d}^{\ast}\boldsymbol{\Sigma}\cdot\boldsymbol{\sigma}^{T}\right.$ (21)
$\displaystyle\left.-2i\tilde{\omega}(\boldsymbol{\Sigma}\times{\bf
d}^{\ast})\cdot\boldsymbol{\sigma}^{T}\right]\frac{i\sigma_{y}\Delta_{0}}{\widetilde{\mathbf{D}}},$

where

$\displaystyle L_{0}$ $\displaystyle=$ $\displaystyle
X_{0}a_{0}-\Delta_{0}^{2}b_{0}|{\bf d}|^{2}-\Delta_{0}^{2}{\bf
q}\cdot\boldsymbol{\Sigma},$ (22) $\displaystyle\mathbf{L}_{1}$
$\displaystyle=$ $\displaystyle(X_{0}+\Delta_{0}^{2}|{\bf
d}|^{2})\boldsymbol{\Sigma}+\Delta_{0}^{2}b_{0}{\bf
q}-\Delta_{0}^{2}\boldsymbol{\Sigma}\cdot{\bf d}{\bf
d}^{\ast}-\Delta_{0}^{2}\boldsymbol{\Sigma}\cdot{\bf d}^{\ast}{\bf d},$ (23)

Here $a_{0}/b_{0}=\tilde{\omega}\mp\xi$,
$X_{0}=b_{0}^{2}-\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}$ and the
denominator
$\widetilde{\mathbf{D}}=(\xi^{2}+\tilde{Q}_{+}^{2})(\xi^{2}+\tilde{Q}_{-}^{2})$,
where $\tilde{Q}_{\pm}$ is,

$\displaystyle\tilde{Q}_{\pm}^{2}$ $\displaystyle=$
$\displaystyle\Delta_{0}^{2}|{\bf
d}|^{2}-\tilde{\omega}^{2}-\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}\pm\sqrt{\Delta_{0}^{4}{\bf
q}\cdot{\bf
q}+4\tilde{\omega}^{2}\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}+4\Delta_{0}^{2}\tilde{\omega}{\bf
q}\cdot\boldsymbol{\Sigma}-4\Delta_{0}^{2}(\boldsymbol{\Sigma}\cdot{\bf
d})(\boldsymbol{\Sigma}\cdot{\bf d}^{\ast})}$ (24)

For the unitary states, the nodes are symmetry imposed, and for the non-
unitary states, nodes may shift away from the high symmetry directions, and
their positions remain protected against disorder as long as the
$\boldsymbol{\Sigma}$ component of the impurity self-energy vanishes,
$\boldsymbol{\Sigma}$ can be interpreted as impurity induced magnetization.
For the chiral states, this term remains finite and give rise to non-
degenerate quasiparticle spin density, and in principle, can change the nodal
positions. Here the nodes do not refer to the zeros of gap or the order
parameter, instead, they are the zeros in the quasi-particle spectrum on the
Fermi surface. In the unitary states, the gap nodes and the spectral nodes are
same unlike the nonunitary states. There are some additional triplet terms
$\boldsymbol{\Sigma}\cdot\boldsymbol{\sigma}$ and
$(\boldsymbol{\Sigma}\times{\bf d})\cdot\boldsymbol{\sigma}$ in Eq. (40),
these terms reflect impurity induced modification of spin structure of the
Cooper pairs. It is worth mentioning that there is an impurity induced odd-
frequency pairing for the chiral non-unitary states, which is spin singlet and
odd parity in nature.

After obtaining the impurity dressed Green’s function, we calculate the
electronic thermal conductivity $\kappa$ using the Kubo formula that relates
the thermal conductivity to the heat-current responseAmbegaokar and Griffin
(1965). We ignore the vertex corrections and restrict ourselves to the bare
thermal current response function. The vertex corrections are small in the
strong scattering limit that we focus on this article(Fledderjohann and
Hirschfeld, 1995). The diagonal thermal conductivity for a general triplet
superconductor reads,

$\frac{\kappa_{ii}}{T}=\int_{-\infty}^{\infty}d\omega\frac{\omega^{2}}{T^{2}}\left(-\frac{dn_{F}(\omega)}{d\omega}\right)\left\langle
N_{0}v_{Fi}^{2}\frac{(-c_{1}c_{4}+c_{2}c_{3})+c_{3}b_{1}+c_{1}b_{2}+b_{3}(-c_{3}+c_{1}c_{2})/c_{4}}{(c_{3}^{2}+c_{1}^{2}c_{4}-c_{1}c_{2}c_{3})}\right\rangle_{FS}.$
(25)

Here $c_{1}=-2\mathrm{Re}(\tilde{Q}_{+}+\tilde{Q}_{-})$,
$c_{2}=|\tilde{Q}_{+}|^{2}+|\tilde{Q}_{-}|^{2}+4\mathrm{Re}(\tilde{Q}_{+})\mathrm{Re}(\tilde{Q}_{-})$,
$c_{3}=-2|\tilde{Q}_{+}|^{2}\mathrm{Re}(\tilde{Q}_{-})-2|\tilde{Q}_{-}|^{2}\mathrm{Re}(\tilde{Q}_{+})$,
$c_{4}=|\tilde{Q}_{+}|^{2}|\tilde{Q}_{-}|^{2}$and the coefficients
$b_{i=1,..,3}$ are,

$\displaystyle b_{1}$ $\displaystyle=$
$\displaystyle(|\tilde{\omega}|^{2}-\Delta_{0}^{2}|{\bf
d}|^{2}+|\boldsymbol{\Sigma}|^{2})+\mathrm{Re}\left[\tilde{Q}_{+}^{2}+\tilde{Q}^{2}_{-}\right]$
(26) $\displaystyle b_{2}$ $\displaystyle=$
$\displaystyle\frac{1}{4}|\tilde{Q}_{+}^{2}+\tilde{Q}^{2}_{-}|^{2}+(|\tilde{\omega}|^{2}+|\boldsymbol{\Sigma}|^{2}-\Delta_{0}^{2}|{\bf
d}|^{2})\mathrm{Re}\left[\tilde{Q}_{+}^{2}+\tilde{Q}^{2}_{-}\right]+3\Delta_{0}^{4}{\bf
q}\cdot{\bf
q}+4|\boldsymbol{\Sigma}|^{2}(|\tilde{\omega}|^{2}+\mathrm{Re}[\tilde{\omega}^{2}])$
(27) $\displaystyle-4\Delta_{0}^{2}\left(|\boldsymbol{\Sigma}\cdot{\bf
d}|^{2}+|\boldsymbol{\Sigma}\cdot{\bf
d}^{\ast}|^{2}+\mathrm{Re}[(\boldsymbol{\Sigma}\cdot{\bf
d})(\boldsymbol{\Sigma}\cdot{\bf d}^{\ast})]\right)$
$\displaystyle+4|\tilde{\omega}|^{2}\mathrm{Re}[\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}]+4\Delta_{0}^{2}\mathrm{Re}[\tilde{\omega}{\bf
q}\cdot\boldsymbol{\Sigma}^{\ast}+2\tilde{\omega}{\bf
q}\cdot\boldsymbol{\Sigma}]$ $\displaystyle b_{3}$ $\displaystyle=$
$\displaystyle(|\tilde{\omega}|^{2}-\Delta_{0}^{2}|{\bf
d}|^{2})\left[|\alpha_{+}|^{2}+\Delta_{0}^{4}{\bf q}\cdot{\bf
q}\right]-2\Delta_{0}^{4}{\bf q}\cdot{\bf
q}\mathrm{Re}[\alpha_{+}]+\mathcal{Y}(\boldsymbol{\Sigma})$ (28)
$\displaystyle{\mathcal{Y}}(\boldsymbol{\Sigma})$ $\displaystyle=$
$\displaystyle\left\\{|(\tilde{\omega}^{2}+\Delta_{0}^{2}|{\bf
d}|^{2}-\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma})|^{2}-\Delta_{0}^{4}{\bf
q}\cdot{\bf
q}\right\\}|\boldsymbol{\Sigma}|^{2}-2\Delta_{0}^{2}\mathrm{Re}[(\alpha_{-}-\alpha_{-}^{\ast})\tilde{\omega}{\bf
q}\cdot\boldsymbol{\Sigma}^{\ast}]$ (29)
$\displaystyle-4\Delta_{0}^{2}(|\tilde{\omega}|^{2}+|\boldsymbol{\Sigma}|^{2}-\Delta_{0}^{2}|{\bf
d}|^{2})\mathrm{Re}[\tilde{\omega}{\bf
q}\cdot\boldsymbol{\Sigma}^{\ast}]+2\Delta_{0}^{4}\mathrm{Re}\left[|{\bf
q}\cdot\boldsymbol{\Sigma}|^{2}-({\bf
q}\cdot\boldsymbol{\Sigma})^{2}+\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}{\bf
q}\cdot{\bf q}\right]$
$\displaystyle+2\Delta_{0}^{2}\mathrm{Re}\left[(|\tilde{\omega}|^{2}-\tilde{\omega}^{2})-(|\boldsymbol{\Sigma}|^{2}-\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma})\right](|\boldsymbol{\Sigma}\cdot{\bf
d}|^{2}+|\boldsymbol{\Sigma}\cdot{\bf d}^{\ast}|^{2})-4\Delta_{0}^{2}|{\bf
d}|^{2}|\tilde{\omega}|^{2}|\boldsymbol{\Sigma}|^{2}$
$\displaystyle+4\Delta_{0}^{2}\mathrm{Re}\left[\alpha_{+}(\boldsymbol{\Sigma}^{\ast}\cdot{\bf
d}^{\ast})(\boldsymbol{\Sigma}^{\ast}\cdot{\bf
d})\right]-4\Delta_{0}^{2}\mathrm{Re}[(\tilde{\omega}^{2}-\Delta_{0}^{2}|{\bf
d}|^{2})\boldsymbol{\Sigma}^{\ast}\cdot\boldsymbol{\Sigma}^{\ast}],$

where $\alpha_{\pm}=\tilde{\omega}^{2}-\Delta_{0}^{2}|{\bf
d}|^{2}\pm\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}$. The derivation of
thermal conductivity is provided in the supplementary materials. We calculate
the full temperature dependence of thermal conductivity using self-
consistently determined superconducting gap using an effective paring
potential to give a single transition temperature (see appendix ).

##  III Results & Discussion

Figure 1: Thermal conductivity for the single component unitary states allowed
by the $D_{2h}$ point group: The thermal conductivity normalized to its value
at $T_{c}$ for the four irreducible representations of the $D_{2h}$ point
group shown column-wise for the $A_{1u}$, $B_{1u}$, $B_{2u}$ and $B_{3u}$
representations from panel (a) to (d), respectively. The first row shows the
residual thermal conductivity in the zero temperature limit as a function of
relative reduction in the transition temperature $\delta t_{c}$. Panels (e) to
(h) in the second row show $\kappa T_{c}/\kappa_{N}T$ for the weak scatterers
with $\tan\theta_{s}=2$ and for the intermediate strength scatterers with
$\tan\theta_{s}=5$ from panels (i) to(l) in the third row. For the temperature
dependence, the $T_{c}$ is reduced by $2\%(\delta t_{c}=0.02)$ w.r.t. the
clean limit.

###  III.1 $T=0$ limit of density of states and thermal conductivity

We start with the single component states based on four irreducible
representations of the $D_{2h}$ point group symmetry. The basis functions for
these four states are listed in Table 1, where the $A_{1u}$ state remains
gapped unless the coefficient of one of the basis functions is set to zero. We
exclude that possibility and choose the same coefficients for all three basis
functions, and this choice of coefficients is adopted for the other states,
too. For qualitative understanding of the low energy properties, this is a
reasonable choice. In principle, it is also possible to generate line nodes
with appropriate choice of basis function coefficients, but those
possibilities are excluded considering the recent experimental results on
UTe2. Apart from the $A_{1u}$ state, the $B_{1u}$ state also remains gapped
because the open Fermi surface along the $\hat{z}$-axis forbids the nodes for
this state. For the $B_{2u}$ and $B_{3u}$ states, there is a pair of point
nodes that exists along the $\hat{y}$-axis and $\hat{x}$-axis, respectively.
Fig. 1(a) to (d) show the disorder dependence of the thermal conductivity
normalized to the normal state value at the transition as a function of
relative reduction in the transition temperature $\delta
t_{c}\equiv=1-T_{c}/T_{c0}$, where $T_{c0}$ is the transition temperature in
the clean limit. For a clean system, $\delta t_{c}=0$ and $\delta t_{c}$
reaches unity as the impurity scattering kill superconductivity. The
normalized residual thermal conductivity $\kappa T_{c}/\kappa_{N}T$ in the
zero temperature limit remains zero up to a threshold value of disorder
$\Gamma_{th}$ for all four states, this $\Gamma_{th}$ corresponds to a
threshold level of $T_{c}$ suppression $\delta t_{c}^{th}$. This threshold
value of disorder depends on the superconducting gap structure and the
strength of the impurity potentialHirschfeld  _et al._ (1986). For the gapped
states $A_{1u}$ and $B_{1u}$, the residual thermal conductivity remains zero
for slightly higher values of disorder compared to the other two states
$B_{2u}$ and $B_{3u}$ as expected due to the presence of impurity induced
quasiparticle states near the nodes. The normalized thermal conductivity
remains very isotropic for the $A_{1u}$ and $B_{1u}$ states. However, for the
nodal states, the residual thermal conductivity is enhanced for thermal
current along the nodal directions. This trend in anisotropy in thermal
conductivity continues even at finite temperatures, as shown in the
temperature evolution of $\kappa T_{c}/\kappa_{N}T$ for the four single
component states in Fig. 1(e) to (h) for $\tan\theta_{s}=2$ and in panels (i)
to (l) for $\tan\theta_{s}=5$, which represents the stronger scatterers. For
both impurity strengths, $\delta t_{c}$ is $0.02$. In the presence of the
point nodes, the thermal conductivity shows a weak maximum at very low
temperature along the nodal direction for the weak scatterers, which
disappears as the scattering rate increases. This is a known behavior for the
superconducting states with point nodesHirschfeld  _et al._ (1988), which is
not present for the stronger impurity potentials.

It is useful to examine the density of states (DOS) and the structure of the
low energy quasiparticle states before discussing the thermal transport for
the nonunitary states. We first report the average density of states per spin
for the nonunitary states on a cylindrical Fermi surface that is open along
the $\hat{z}$-axis. Fig. 2(a) shows the DOS for the $A_{1u}+irB_{2u}$ state,
which is a chiral state, and Fig. 2(b), (c) and (d) show $\Delta_{-}$ for this
state. In the gapped phase ($r<1$), this state has minima near the
$\hat{x}$-axis and a small gap is visible in the DOS. For $r=1$, a pair of
point nodes appears along the $\hat{x}$-axis and the low energy DOS shows
$\omega^{2}$ behavior that is expected for linear point nodes. Here and in
subsequent discussion, “node” refers to spectral nodes, not the gap nodes.
However, in contrast to a unitary state, the low energy quadratic behavior
remains confined to very low energy scale as compared to the unitary $B_{3u}$
state. For $r>1$, the nodes move away from the $\hat{x}$-axis and the
positions of four nodes are determined by $\tan\phi=\pm\sqrt{r^{2}-1}$, where
$\phi$ is the polar angle on the cylindrical Fermi surface. As the mixing
parameter $r$ increases, additional pairs of nodes appear in the $y{}z$ plane
at $\sin(k_{z}/2)=\pm 1/\sqrt{r^{2}-1}$. In the $r\rightarrow\infty$ limit,
only two point nodes along $\hat{y}$-axis survive, as expected for a pure
$B_{2u}$ state. The low energy DOS remains quadratic in all these cases. For
$1\leq r<\sqrt{2}$, the nodes remain closer to the $\hat{x}$-axis, for
$r\geq\sqrt{2}$, the nodes move closer to $\hat{y}$-axis. The
$A_{1u}+irB_{3u}$ state also shows similar DOS as $A_{1u}+irB_{2u}$ state, but
it has different nodal structure. It has gap minima in quasiparticle spectrum
near the $\hat{y}$-axis, and nodes appear along the $\hat{y}$-axis. For $r>1$,
a set of four nodes appear near $\hat{y}$-axis at
$\phi=\pm\cot^{-1}\sqrt{r^{2}-1}$ and move towards the $\hat{x}$-axis in the
$x{}y$ plane as the value of $r$ increases. For $r>\sqrt{2}$, four more nodes
appear at $k_{z}=\pm 2\sin^{-1}(1/\sqrt{r^{2}-1})$ in the $x{}z$ plane. Both
these states are chiral and show finite quasiparticle spin density along
$\hat{y}$ and $\hat{x}$ directions in the spin space. The last AF non-unitary
state is $A_{1u}+irB_{1u}$, which is not a chiral state. It has nodes along
the $\hat{x}$ and $\hat{y}$ directions for $r=1$, a quadruple pair of point
nodes in the $xz$-plane, and another set of four point nodes $yz$ plane, where
the $k_{z}$ for the nodal position is determined by
$\sin(k_{z}/2)=\pm\sqrt{r^{2}-1}$. The DOS shows quadratic behavior at low
energies (see appendix).

Figure 2: Density of states: Panel (a) shows the DOS per spin for the
$A_{1u}+irB_{2u}$ state for different values of the mixing parameter $r$. The
inset shows the energy dependence of DOS at low energies. A dashed line
illustrating the DOS for the unitary $B_{3u}$ state has been shown for
comparison. Panel (b), (c) and (d) show the variation of $\Delta_{-}$ for
different values of $r$ in the $xy$, $xz$ and $yz$ planes, respectively. Panel
(f) shows the DOS per spin for the $B_{1u}+irB_{2u}$ state and $\Delta_{-}$ in
the $x{}y$, $x{}z$ and $y{}z$ planes are shown in panels (f), (g) and (h),
respectively.

Next, we consider the ferromagnetic nonunitary states on the cylindrical Fermi
surface, which are chiral states with finite Cooper pair spin moment. Fig.
2(d) shows the DOS for the $B_{1u}+irB_{2u}$ state along with the $\Delta_{-}$
in the panels (e) and (f). This state is gapped for $r<1/\sqrt{2}$, and a
quadruple pair of point nodes appears in the $y{}z$ plane for $1/\sqrt{2}\leq
r<1$, whose positions are determined by $\sin(k_{z}/2)=r/\sqrt{1-r^{2}}$. For
$r>1$, a quadruple pair of point nodes appear in the $x{}y$ plane at
$\phi=\pm\tan^{-1}\sqrt{r^{2}-1}$ close to the $\hat{x}$-axis and move towards
$\hat{y}$-axis as $r$ increases. These states show quadratic DOS at low
energies. For $r=1$, this state shows point nodes along the $\hat{x}$-axis,
but these are the quadratic point nodes. The first derivative of the
$\Delta_{-}$ vanishes at the nodes for a quadratic or second order point
nodes. This leads to linear DOS at low energies, as shown in the Fig. 2(d).
Similarly, for the $B_{1u}+irB_{3u}$ a twin pair of quadratic point nodes
appear along the $\hat{y}$-axis and shows linear DOS at low energies for
$r=1$. For $1/\sqrt{2}\leq r<1$ four point nodes appear at
$\sin(k_{z}/2)=\pm\sqrt{1-r^{2}}/r$ in the $y{}z$ plane, and for $r>1$ there
is a quadruple pair of point nodes in the $x{}y$ plane at
$\cot\phi=\pm\sqrt{r^{2}-1}$, closer to the $\hat{y}$-axis for $r\geq 1$,
which move towards the $\hat{x}$-axis for $r\gg 1$. The DOS remains quadratic,
as expected. For the $B_{2u}+irB_{3u}$ state, four point nodes are either
located in the $y{}z$ plane at $\sin(k_{z}/2)=\pm r/\sqrt{1-r^{2}}$ for $r\leq
1/\sqrt{2}$ or in the $x{}z$ plane at $\sin(k_{z}/2)=\pm 1/\sqrt{r^{2}-1}$ for
$r\geq\sqrt{2}$. These states show $\omega^{2}$ behavior in the low energy
DOS. For $1/\sqrt{2}<r<\sqrt{2}$, there exists a gap in the quasiparticle
spectrum due to lack of nodes (see supplementary materials). In contrast to
the AF nonunitary states on the cylindrical Fermi surface, the FM nonunitary
states can have at most four nodes and are expected to be more anisotropic.

Figure 3: Residual thermal conductivity: Zero temperature limit $\kappa
T_{c}/\kappa_{N}T$ for the $A_{1u}+irB_{1u}$ in panels (a) to (d) and for the
$B_{2u}+irB_{3u}$ state in panels (e) to (g) for various values of mixing
parameter $r$. Panel (h) shows $\kappa T_{c}/\kappa_{N}T|_{T\rightarrow 0}$
for the unitary $B_{3u}$ state. The $s$-wave scattering phase shift is
$\theta_{s}=\tan^{-1}(2)$ for all the panels. The number of nodes $N_{n}$ is
indicated for each case.

We now discuss the zero temperature limit of the thermal conductivity, which
is very sensitive to the gap structure. We compare the normalized $\kappa
T_{c}/T\kappa_{N}$ along three principle directions, where the thermal
conductivity along a particular direction is normalized to its normal state
value at $T_{c}$ along that direction. This suppresses the intrinsic
anisotropy present in the electronic structure and highlights the effect of
order parameter anisotropy. Panels (a) to (d) in Fig. 3 show the $\kappa
T_{c}/T\kappa_{N}$ in the zero temperature limit for the $A_{1u}+irB_{1u}$
state as a function of relative reduction in the transition temperature
$\delta T_{c}$. For a weakly disordered system, $\kappa/T$ vanishes in the
zero temperature limit, but as the disorder level increases and crosses
$\Gamma_{th}$, $\kappa/T|_{T\rightarrow 0}$ becomes finite and reaches the
normal state value as the superconductivity vanishes. For the
$A_{1u}+irB_{1u}$ state, the $\kappa T_{c}/\kappa_{N}T$ shows isotropic
behavior in the $x{}y$-plane. Note, for this state, the nodes are always along
the $\hat{x}$ and $\hat{y}$ axes at the same $k_{z}$ value. In contrast, the
$B_{2u}+irB_{3u}$ state shows relatively weaker level anisotropy than a
unitary state. Panels (e) to (g) in Fig. 3 show the $\kappa T_{c}/\kappa_{N}T$
for the $B_{2u}+irB_{3u}$ state, while Fig. 3(h) shows $\kappa
T_{c}/\kappa_{N}T$ for the unitary $B_{3u}$ state with nodes along $\hat{x}$
axis. As shown in Fig. 3(e), with nodes in the $y{}z$-plane $\kappa_{yy}$ and
$\kappa_{zz}$ exceed $\kappa_{xx}$, and this behavior reverses as the nodes
move to the $x{}z$-plane for the $B_{2u}+irB_{3u}$ state, as depicted in Fig.
3(h). As the value of $r$ increases, the anisotropy also reduces, and for
$r=1$ the superconducting state becomes fully isotropic, as shown in Fig.
3(f). This is a special case, which has a four-fold symmetric $\Delta_{-}$,
leads to a very isotropic normalized thermal conductivity along the three
principal directions. As $r$ becomes larger than unity, the $\kappa_{xx}$
starts to dominate.

Figure 4: Residual thermal conductivity: Panels (a) to (h) show $\kappa/T$ in
the zero temperature limit normalized to $\kappa_{N}/T$ at $T_{c}$ as a
function of relative $T_{c}$ reduction $\delta t_{c}$ for the
$A_{1u}+irB_{2u}$ state for different values of the mixing parameter $r$. The
$s$-wave scattering phase for the panels (a) to (d) is
$\theta_{s}=\tan^{-1}(2)$ and for the panels (e) to (h) is
$\theta_{s}=\tan^{-1}(5)$. For the $B_{2u}+irB_{3u}$ state, $\kappa
T_{c}/\kappa_{N}T|_{T\rightarrow 0}$ is shown for various values of parameter
$r$ for $s$-wave scattering phase shift $\theta_{s}=\tan^{-1}(2)$ in the
panels (i) to (l) and with $\theta_{s}=\tan^{-1}(5)$ in the panels (m) to (p).
The number of nodes $N_{n}$ is indicated for each case.

Next, we look at the zero temperature limit thermal conductivity for the
$A_{1u}+irB_{2u}$ state for two different impurity potential strengths. Fig. 4
(a) to (d) show $\kappa T_{c}/\kappa_{N}T$ for $s$-wave scattering phase shift
$\tan^{-1}(2)$ and panels (e) to (h) show $\kappa T_{c}/\kappa_{N}T$ for
$\theta_{s}=\tan^{-1}(5)$. The $\kappa/T|_{T\rightarrow 0}$ becomes finite
above a threshold disorder level, like earlier cases. This threshold
scattering rate is smaller for the stronger impurity potentials. This state
has minimum and maximum along the $xz$-plane, the energy gap is small in the
$yz$-plane and there is a weak maximum along the $\hat{y}$-axis in the
$xy$-plane. This excitation energy spectrum is reflected in the thermal
conductivity. As the $\kappa/T$ becomes finite, the in-plane anisotropy is
very weak with slightly larger value along the $\hat{x}$-axis, as long as
there is no node in the $yz$-plane. As the impurity induced quasiparticles
overcome the minima along the $\hat{y}$ axis, the thermal conductivity along
the $\hat{y}$ axis starts to dominate. In case of eight nodes, with four in
the $xy$-plane and another four $yz$-plane, the thermal conductivity is always
larger along the $\hat{y}$-axis. The $\hat{z}$-axis thermal conductivity
remains close to the $\hat{x}$-axis thermal conductivity. For the
$A_{1u}+irB_{3u}$ state, the thermal conductivity along the $\hat{x}$-axis and
$\hat{y}$-axis show the same behavior as the $\kappa_{yy}/T$ and
$\kappa_{xx}/T$ in the $A_{1u}+irB_{2u}$ state, respectively.

For the $B_{1u}+irB_{3u}$ state, in the gapped phase i.e. $r<1/\sqrt{2}$, the
spectral gap on the Fermi surface is small in the $yz$-plane, and near the
$\hat{y}$ axis in the $xy$-plane, which leads to larger
$\kappa/T|_{T\rightarrow 0}$ along the $\hat{y}$ axis followed by the
$\hat{z}$ direction, as shown in Fig. 4(i) and Fig. 4(m), for
$\tan\theta_{s}=2$ and $\tan\theta_{s}=5$, respectively. As a pair of
quadruple nodes appears in the $yz$-plane, the relative anisotropy remains the
same, as depicted in Fig.4(j) and Fig.4(n). This state also shows an elusive
quadratic node for $r=1$ along the $\hat{y}$-axis. Unlike the other states
with linear point nodes, for this case the $\kappa/T|_{T\rightarrow 0}$ term
remain finite along the nodes and for other directions residual thermal
conductivity remains zero below the threshold disorder level. This state shows
linear DOS at low energies, like superconductors with line nodes. Finally, the
nodes appear in the $xy$-plane for this state as $r$ goes beyond unity and
remains closer to the $\hat{y}$-axis, and thermal conductivity along the
$\hat{y}$-axis becomes dominant, while other two directions show very similar
$\kappa T_{c}/\kappa_{N}T$. As $r$ increases, the in-plane anisotropy
decreases and the $\hat{x}$-axis $\kappa$ increases and becomes stronger along
the $\hat{y}$ direction in $r\gg 1$ limit. For the $B_{1u}+irB_{2u}$ state,
the in-plane anisotropy found for the $B_{1u}+irB_{3u}$ state gets
interchanged.

Figure 5: $A_{1u}+irB_{1u}$: Thermal conductivity normalized to its normal
state value at $T_{c}$ for the $A_{1u}+irB_{1u}$ shown as a function of
temperature normalized to $T_{c}$ for scatterers with $\tan\theta_{s}=5$ for
$\delta t_{c}=0.02$ from panels (a) to (d) and for $\delta t_{c}=0.15$ in
panels (e) to (h) for various values of mixing parameter $r$.

###  III.2 Finite $T$ electronic thermal conductivity

Now, we look at the temperature evolution of the normalized thermal
conductivity for the nonunitary states. At very low temperatures, the elastic
scattering by the impurities is the main mechanism of relaxation. However, as
the temperature increases the inelastic scattering also becomes important,
which we discuss in the subsequent section. Apart from the electronic
contribution to thermal conductivity, phonon thermal conductivity can also
become significant. Here, we focus on electronic thermal conductivity and the
effect of impurity scattering on it, and the effect of underlying spectral
nodes on the anisotropy in the thermal conductivity. As shown in the Fig. 1,
the thermal conductivity shows a weak maximum as function of temperature at
very low temperatures along the nodal directions. There is no evidence for
such a feature in the experimental measurementsMetz  _et al._ (2019); Suetsugu
_et al._ (2024); therefore, we have set $\theta_{s}=\tan^{-1}(5)$ for the rest
of our discussion (see supplementary materials for the weaker impurity
scatterers). We first look at the $A_{1u}+irB_{1u}$ state, which shows very
isotropic residual thermal conductivity, and the thermal conductivity remains
isotropic as a function of temperature, as shown in Fig. 5(a) to (d). For the
clean system with $\delta t_{c}=0.02$, at very low temperatures, the thermal
conductivity remains vanishingly small and as the temperature increases,
$\kappa/T$ increases $T^{2}$ at very low temperatures. This behavior is
observed for both nodal and gapped systems. One should keep in mind that in
superconductors with point nodes, there are very few states available at the
Fermi energy, as shown in Fig. 2. For dirtier systems such as $\delta
t_{c}=0.15$, there are sufficient quasiparticle states at the Fermi level to
give nonzero thermal conductivity, which is depicted in Fig. 5(e) to (h). At
low temperatures, $\kappa/T$ remains independent of temperature, and it begins
to increase, once more quasiparticles become relevant for transport as the
temperature increases.

Figure 6: $A_{1u}+irB_{2u}$: Thermal conductivity normalized to its normal
state value at $T_{c}$ for the $A_{1u}+irB_{2u}$ shown as a function of
temperature normalized to $T_{c}$ for a scatterers with $\tan\theta_{s}=5$ for
$\delta t_{c}=0.02$ from panels (a) to (d) and for $\delta t_{c}=0.15$ in
panels (e) to (h) for various values of mixing parameter $r$. The dashed lines
show the normalized thermal conductivities for the $B_{2u}$ state with same
impurity parameters.  Figure 7: $B_{2u}+irB_{3u}$: Thermal conductivity
normalized to its normal state value at $T_{c}$ for the $B_{2u}+irB_{3u}$
shown as a function of temperature normalized to $T_{c}$ for a scatterers with
$\tan\theta_{s}=5$ for $\delta t_{c}=0.02$ from panels (a) to (c) and for
$\delta t_{c}=0.15$ in panels (d) to (f) for various values of mixing
parameter $r$. The dashed lines show the normalized thermal conductivities for
the $B_{2u}$ state in panels (a) and (d) and for the $B_{3u}$ state in panels
(c) and (f). The impurity parameters are the same for the single component
superconducting states.

Next, we look at the $A_{1u}+irB_{2u}$ state, which shows enhanced thermal
conductivity along the $\hat{y}$-axis. The zero temperature trends in
anisotropy continue as the temperature increases, as shown in Fig. 6. For
$\delta t_{c}=0.02$, $\kappa/T$ is close to zero and rises as temperature
increases beyond a critical value. For the gapped case ($r<1$) and for the
$\hat{x}$-axis nodes ($r=1$), $\kappa/T$ is slightly larger than other two
directions, but as the temperature increases, the other two directions begin
to increase and dominate because the temperature overcomes the gap minima in
the $yz$-plane. For $r>1$, there are four point nodes which move towards the
$\hat{y}$ axis, and it is reflected as a larger $\kappa/T$ along the
$\hat{y}$-axis. For $r>\sqrt{2}$, four more nodes in $yz$-plane make
$\kappa/T$ along $\hat{y}$ and $\hat{z}$ directions larger than the $\hat{x}$
direction. For dirtier systems, thermal conductivity becomes finite and
$T$-independent at very low temperature like the earlier case, but the
anisotropy remains similar to its zero temperature limit. For
$A_{1u}+iB_{3u}$, the behavior of $\hat{x}$ and $\hat{y}$ directions get
interchanged (see supplementary materials).

Next, we consider the $B_{2u}+irB_{3u}$ state, which is one of the FM
nonunitary chiral states. For this state, there are four point nodes in
$yz$-plane for $r<1/\sqrt{2}$ and for $r>\sqrt{2}$, there is a quadruple pair
of nodes in the $xz$-plane. Fig. 7 shows the thermal conductivity as a
function of temperature, and as expected, $\kappa/T$ along the nodal
directions dominate. For $r<1$, the thermal conductivity is enhanced along the
$\hat{y}$ and $\hat{z}$ directions, and as $r$ increases, the system becomes
gapped with enhanced thermal conductivity in the $\hat{y}$ and $\hat{z}$
directions, but the anisotropy reduces. For $r=1$, this state shows a
completely isotropic behavior and as $r$ goes beyond unity, this state show
$\hat{x}$-axis dominated thermal response, as the nodes reappear in the
$xz$-plane. The anisotropy remains qualitatively same as the temperature
increases.

Figure 8: $B_{1u}+irB_{3u}$: Thermal conductivity normalized to its normal
state value at $T_{c}$ for the $B_{1u}+irB_{3u}$ shown as a function of
temperature normalized to $T_{c}$ for a scatterers with $\tan\theta_{s}=5$ for
$\delta t_{c}=0.02$ from panels (a) to (d) and for $\delta t_{c}=0.15$ in
panels (e) to (h) for various values of mixing parameter $r$. The dashed lines
show the normalized thermal conductivities for the $B_{2u}$ state with same
impurity parameters.

Finally, we look at the $B_{1u}+irB_{3u}$ state, which is another possible FM
nonunitary state that shows enhanced thermal conductivity along the $\hat{y}$
directions. Here, the anisotropy changes significantly as the temperature
increases. At low temperatures, $\hat{y}$ axis dominates due to its vicinity
to the point nodes, however, as the temperature increases, we find enhancement
of $\kappa T_{c}/\kappa_{N}T$ along the $\hat{x}$ directions. This happens
because at the lower temperatures, the lower energy branch of the
quasiparticle excitations $\sqrt{\xi_{\bf k}^{2}+\Delta_{-}^{2}}$ dominates,
which has more quasi-particle states along the $\hat{y}$ directions, but at
higher temperatures $\sqrt{\xi_{\bf k}^{2}+\Delta_{+}^{2}}$ branch of the
quasiparticle excitations begin to contribute, which for this state has minima
along the $\hat{x}$ axis. The overall anisotropy for this state is relatively
less compare to other states, except those that show fully isotropic behavior
as a function of temperature or disorder. For the $B_{1u}+irB_{2u}$ state,
$\kappa/T$ along the $\hat{x}$ and $\hat{y}$ directions gets interchanged.

###  III.3 Inelastic scattering effects

As we mentioned in the previous section, at very low temperatures i.e. $T\ll
T_{c}$, the elastic scattering from impurities is the only mechanism that
determines the scattering rate. However, at higher temperature inelastic
scattering from a bosonic mode is possible. We consider a simple scenario
where there is a dispersionless bosonic mod that couples with the fermions
with an effective coupling constant $g_{fb}$. We further assume that the
coupling does not depend on the fermion’s spin degree of freedom. The lowest
order self-energy for the fermions reads,

$\displaystyle\Sigma_{in}(i\omega_{n},\mathbf{k})$ $\displaystyle=$
$\displaystyle
g^{2}_{fb}T\sum_{m,\mathbf{q}}G(i\omega_{n}-i\Omega_{m},\mathbf{k}-\mathbf{q})D(i\Omega_{m},\mathbf{q}).$
(30)

Here $D$ is the bosonic Green’s function, $\omega_{n}$ and $\Omega_{m}$ are
the fermionic and bosonic Matsubara frequencies. After performing the
Matsubara summation, we get,

$\displaystyle\Sigma_{in}(\omega)$ $\displaystyle=$
$\displaystyle-\frac{g^{2}_{fb}}{2\pi}\int_{-\infty}^{\infty}dx\int_{-\infty}^{\infty}dy\frac{N(y)D^{\prime\prime}(x)}{x+y-\omega-i\eta}\left(\coth(\beta
x/2)+\tanh(\beta y/2)\right).$ (31)

Here the real part of the self-energy contributes to mass renormalization and
the imaginary part modifies the scattering rate, which is a function of
quasiparticle energy and temperature. Noting that the UTe2 has very high
effective mass, we ignore the real part of the self energy. The imaginary part
of the self-energy is,

$\Sigma_{in}^{\prime\prime}(\omega)=-\frac{g^{2}_{fb}}{2}\int_{-\infty}^{\infty}dxN(\omega-x)D^{\prime\prime}(x)\left(\coth(\beta
x/2)+\tanh(\beta(\omega-x)/2)\right),$ (32)

where $D^{\prime\prime}(x)\equiv x/(x^{2}+\Omega_{0}^{2})$ is the bosonic DOS,
and $\Omega_{0}$ is the characteristic energy scale associated to the bosonic
mode. In the context of UTe2, we expect $\Omega_{0}\gg T_{c}$, therefore, we
approximate the bosonic DOS as $D^{\prime\prime}(x)\approx x/\Omega_{0}^{2}$.
In the zero temperature limit,
$\Sigma_{in}^{\prime\prime}\propto\omega^{n+2}$, where the DOS for the
superconducting state behaves like $\omega^{n}$ at low energies. For linear
point nodes, $\Sigma_{in}^{\prime\prime}\propto\omega^{4}$ and for quadratic
point nodes or line nodes, it behaves like $\omega^{3}$. Similarly, in the
static limit ($\omega\rightarrow 0$), the imaginary part of the inelastic
self-energy reads,

$\Sigma_{in}^{\prime\prime}(\omega=0,T)=-{2g^{2}_{fb}}\int_{0}^{\infty}dx\frac{N(x)D^{\prime\prime}(x)}{\sinh\beta
x}.$ (33)

Since the thermal response integrand is peaked at $\omega=0$, we only retain
the temperature dependence of the inelastic self-energy at $\omega=0$. At very
low temperatures, the inelastic scattering rate behaves likes $T^{4}$ for the
point nodes and $T^{3}$ for the line nodes in the static limit. This is
sufficient to understand the qualitative effect of the inelastic scattering.
The prefactor in Eq. (33) is fixed by the value of the inelastic scattering
rate at $T_{c}$.

Figure 9: Effect of inelastic scattering scattering: Panels (a) to (d) show
temperature dependence of normalized thermal conductivity for the
$A_{1u}+irB_{2u}$ state for $r=1.2$. Panels (e) to (h) show the temperature
dependence of the normalized thermal conductivity for the $B_{1u}+irB_{3u}$
state. The imaginary part of the the inelastic self energy at $T=T_{c0}$ and
relative reduction in $T_{c}$ are indicated in each panel. For brevity, we
denote $\Sigma^{\prime\prime}_{in}(T=T_{c0})$ as
$\Sigma^{\prime\prime}_{T=T_{c0}}$. The $s$-wave scattering phase is
$\tan^{-1}(5)$.

Fig. 9 shows the effect of the inelastic scattering on the electronic thermal
conductivity. We show the results for two cases, as the rest of the cases are
qualitatively similar. The first row of the Fig. 9 shows the normalized
thermal conductivity for the $A_{1u}+irB_{2u}$ state, which show enhanced
thermal conductivity along the $\hat{y}$-axis. Note, the normal state thermal
conductivity includes the inelastic scattering. We have fixed the value of $r$
at $1.2$, which leads to four point nodes in the $xy$-plane. The primary
effect of inelastic scattering is formation of a peak below $T_{c}$. This peak
appears first for the directions that have larger thermal conductivity, as
shown in Fig. 9(a) for $\delta=0.02$. As the system becomes dirtier, the peaks
get smeared as depicted in Fig. 9(b) for $\delta t_{c}=0.15$. If the strength
of inelastic scattering is increased, the peaks also strengthen, as shown in
Fig. 9(c) and it may survive in dirtier systems, as illustrated in Fig. 9(d).
The similar trends continue for the other states also. In panels 9(e) to (h),
we show the normalized thermal conductivity for the $B_{1u}+irB_{3u}$ state,
which shows qualitatively same behavior as a function of disorder and
inelastic scattering strength. However, this state shows a nonmonotonic
variation of anisotropy, and the peaks for different directions remain quite
close to each other.

##  IV Summary & Concluding Remarks

We studied low energy quasiparticle excitations and the thermal transport for
the single component and two components pairing states allowed by the
irreducible representations of the $D_{2h}$ point group symmetry, which is
relevant for the orthorhombic UTe2 crystals. Since $D_{2h}$ has only one
dimensional representations, the four single component pairing states
corresponding to the irreducible representations can not break the time
reversal symmetry, and describe unitary triplet pairing states. Therefore, we
also considered the pairing states that are combinations of two of the four
irreducible representations using a single mixing parameter $r$. We examined
all six two component superconducting states as a function of $r$ on a
cylindrical Fermi surface, which describe either gapped states or state with
spectral point nodes depending on the value of the mixing parameter. The
spectral point nodes are not necessarily the zeros of order parameters, but
these are the points on the Fermi surface hosting quasiparticle excitations.
Of course, the spectral nodes are identical to the gap nodes in the case of
single component or unitary states. These six states can be divided into AF or
FM categories depending on the Fermi surface average of the Cooper pair spin-
moment, which vanishes for the AF states and remains finite for the FM states.
Except for the $A_{1u}+irB_{1u}$ state, all other states are chiral on the
cylindrical Fermi surface, as the average angular momentum of the Cooper pairs
remains finite.

After introducing the one parameter model for the two component states, we
calculated the effect of impurity scattering within the self-consistent
t-matrix approximation. One of the new findings is the spin-dependent impurity
scattering rate for the chiral states. This happens due to finite
quasiparticle spin densities for the chiral states. This can be interpreted as
accumulation of magnetization near the impurity sites and this leads to
qualitative changes in the quasiparticle excitation spectrum. For the two
component states, the nodes are accidental not symmetry imposed like the
single component states, and the spin-dependent self-energy or the impurity
pinned magnetization can change the position of the spectral nodes. In
principle, the removal of spectral point nodes by impurity scattering is
possible in the chiral superconducting states, but we have not found such
effect for the cases that we consider. Next, we calculate the thermal
conductivity using the Kubo formula for the thermal current response function.
Due to spin-dependent impurity self-energies, the thermal conductivity
significantly differs from the thermal conductivity reported for the unitary
states. We examine the thermal transport for all the single and double
component states that are possible for the $D_{2h}$ point group.

We have considered a single band with cylindrical Fermi surface in our
theoretical calculations. The quantum oscillation experiment reports two
cylindrical Fermi surfaces, where one of them is an electron-like and the
other one is a hole-like Fermi surface with comparable effective masses(Aoki
_et al._ , 2022b). However, we expect our analysis to be valid for two band
system too, because the interband scattering is always pair-breaking due to
odd-parity order parameters. The impurity will renormalize the quasiparticle
energies, there will not be any off-diagonal impurity self-energy. Therefore,
multiple bands will lead to higher impurity scattering rate, this should not
affect the anisotropy of the thermal transport. There are also some
speculations about a closed Fermi surface near the $\mathbf{Z}$ point,
therefore, we also considered a spherical Fermi surfaces (see appendix). The
key qualitative difference is the possibility of nodes along the $\hat{z}$
axis and strong thermal conductivity along that directions, however, there is
no experimental data available to support that scenario.

Based on our thermal transport study and some recent experimental data, we can
identify some states that could possibly describe the gap structure in UTe2.
Definitive conclusions are not possible at this time due to a lack of
sufficient direction-dependent data on the newer samples, but we can make some
qualitative statements and rule out some states. For the $A_{1u}+irB_{1u}$
state on a cylindrical Fermi surface, the normalized thermal conductivity
shows isotropic behavior as a function of impurity scattering in the zero
temperature limit and in its temperature dependence for a fixed disorder
level. The limited data that are available for the thermal conductivity
indicates weak in-plane anisotropy, but not absolute isotropic behaviorMetz
_et al._ (2019). The thermal conductivity measurements by Suetsugu et
al.Suetsugu  _et al._ (2024) claims a fully gapped superconducting state. In
contrast, another independent thermal conductivity measurements by Hayes et
al. find evidence for point nodes without finding any residual thermal
conductivity in the zero temperature limitHayes  _et al._ (2024). Therefore,
absence of finite zero temperature limit thermal conductivity in high quality
samples is not sufficient to rule out point nodes. Other probes such as the
field dependence of specific heat suggest a superconducting state with nodes
closer to $\hat{y}$-axisLee  _et al._ (2023). The superfluid density
measurement indicate stronger low energy quasiparticle excitations along the
$\hat{y}$ directions and weakest along the $\hat{x}$ axisIshihara  _et al._
(2023). The relative anisotropy in the penetration depth measurements is
weaker compare to the single component states. We find that the
$A_{1u}+irB_{2u}$ state, the $B_{2u}+irB_{3u}$ state with dominant $B_{2u}$
component and the $B_{1u}+irB_{3u}$ state show stronger quasi-particle
excitations along the $\hat{y}$-axis. The $B_{1u}+irB_{3u}$ state shows change
of anisotropy as a function of temperature, which can be used to distinguish
it from other two states. This state also shows a quadratic point node for
$r=1$, which shows linear DOS at low energies, hence, it can be also ruled
out. The two ferromagnetic states show higher thermal conductivity along the
$\hat{y}$ direction than the $B_{2u}$ state, which has point nodes along the
$\hat{y}$-axis.

As we mentioned earlier that the phonon thermal conductivity could be
significant, especially in the low $T_{c}$ samples. At very low temperatures,
the scattering of phonons from defects dictates the phonon mean free path,
hence, the phonon thermal conductivityKlemens (1955); Bardeen  _et al._
(1959); Uher (1990). Phonon thermal conductivity is expected to be
insignificant in the high quality samples due to low concentration of defects.
Therefore, a systematic measurement of thermal conductivity along all three
directions in the samples with high residual resistivity ratios is highly
desirable. One common feature among all these superconducting phases is zero
$\kappa/T$ in the zero temperature limit in clean samples. For sufficient
disorder, even point nodal states acquire a very small residual $kappa/T$,
which would require very careful low-T measurements to detect.

## Acknowledgments

The authors are grateful to S. Anlage, I. Hayes, T. Metz, J. P. Paglione, and
T. Shibauchi for useful discussions. P.J.H. was supported by the NSF under
DMR-2231821.

##  Appendix A Thermal Conductivity

The heat current vertex $\check{\mathbf{j}}_{Q_{i}}$ along the ith
($i=\hat{x},\hat{y},\hat{z}$) reads,

$\displaystyle\check{\mathbf{j}}_{Q_{i}}$ $\displaystyle=$
$\displaystyle(\omega+\Omega/2)\left[v_{Fi}\tau_{3}\otimes\sigma_{0}\right],$
(34)

where $i$ denotes the direction of the heat current in the real space along
$\hat{x}/\hat{y}/\hat{z}$. The thermal conductivity $\kappa$ is obtained from
the current-current correlation function for the thermal current(Ambegaokar
and Griffin, 1965). After doing the Matsubara summation, thermal conductivity
becomes,

$\frac{\kappa}{T}=\frac{1}{\pi}\int_{-\infty}^{\infty}\frac{d\omega}{T}\frac{\omega^{2}}{T^{2}}\left(-\frac{dn_{F}(\omega)}{d\omega}\right)\frac{\mathrm{Tr}\left[\check{\mathcal{A}}({\bf
k},\omega)\check{\mathbf{j}}_{Q_{i}}\check{\mathcal{A}}({\bf
k},\omega)\check{\mathbf{j}}_{Q_{i}}\right]}{2}.$ (35)

Here a factor of half accounts for the double counting of the spins while
taking the trace, $n_{F}$ is the Fermi function and $\check{\mathcal{A}}$ is
the spectral function, which reads,

$\displaystyle\check{\mathcal{A}}$ $\displaystyle=$
$\displaystyle\frac{\left[\check{\mathbf{G}}({\bf
k},\omega+i\eta)-\check{\mathbf{G}}({\bf
k},\omega-i\eta)\right]}{2i}=\frac{(\check{\mathbf{G}}^{R}-\check{\mathbf{G}}^{A})}{2i},$
(36)

where the impurity dressed Green’s function $\check{\mathbf{G}}$ is,

$\displaystyle\check{\mathbf{G}}$ $\displaystyle=$
$\displaystyle\begin{pmatrix}\hat{\mathbf{G}}_{11}&\hat{\mathbf{G}}_{12}\\\
\hat{\mathbf{G}}_{21}&\hat{\mathbf{G}}_{22}\end{pmatrix},$ (37)

where the block Green’s functions in the spin-space are,

$\displaystyle\hat{\mathbf{G}}_{11}$ $\displaystyle=$
$\displaystyle\frac{L_{0}+\mathbf{L}_{1}\cdot\boldsymbol{\sigma}}{\widetilde{\mathbf{D}}},$
(38) $\displaystyle\hat{\mathbf{G}}_{22}$ $\displaystyle=$
$\displaystyle\frac{L_{0}(\xi\rightarrow-\xi)+\mathbf{L}_{1}(\xi\rightarrow-\xi)\cdot\boldsymbol{\sigma}^{T}}{\widetilde{\mathbf{D}}},$
(39) $\displaystyle\hat{\mathbf{G}}_{12}$ $\displaystyle=$
$\displaystyle\left[2\xi(\boldsymbol{\Sigma}\cdot{\bf
d})+(\tilde{\omega}^{2}-\xi^{2}-\Delta_{0}^{2}|{\bf
d}|^{2}+\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}){\bf
d}\cdot\boldsymbol{\sigma}+i\Delta_{0}^{2}({\bf q}\times{\bf
d})\cdot\boldsymbol{\sigma}-2(\boldsymbol{\Sigma}\cdot{\bf
d})\boldsymbol{\Sigma}\cdot\boldsymbol{\sigma}\right.$ (40)
$\displaystyle\left.+2i\tilde{\omega}(\boldsymbol{\Sigma}\times{\bf
d})\cdot\boldsymbol{\sigma}\right]\frac{i\sigma_{y}\Delta_{0}}{\widetilde{\mathbf{D}}},$
$\displaystyle\hat{\mathbf{G}}_{21}$ $\displaystyle=$
$\displaystyle\left[-2\xi\boldsymbol{\Sigma}\cdot{\bf
d}^{\ast}+(\tilde{\omega}^{2}-\xi^{2}-\Delta_{0}^{2}|{\bf
d}|^{2}+\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}){\bf
d}^{\ast}\cdot\boldsymbol{\sigma}^{T}-i\Delta_{0}^{2}({\bf q}\times{\bf
d}^{\ast})\cdot\boldsymbol{\sigma}^{T}-2\boldsymbol{\Sigma}\cdot{\bf
d}^{\ast}\boldsymbol{\Sigma}\cdot\boldsymbol{\sigma}^{T}\right.$ (41)
$\displaystyle\left.-2i\tilde{\omega}(\boldsymbol{\Sigma}\times{\bf
d}^{\ast})\cdot\boldsymbol{\sigma}^{T}\right]\frac{i\sigma_{y}\Delta_{0}}{\widetilde{\mathbf{D}}},$

where

$\displaystyle L_{0}$ $\displaystyle=$ $\displaystyle
X_{0}a_{0}-\Delta_{0}^{2}b_{0}|{\bf d}|^{2}-\Delta_{0}^{2}{\bf
q}\cdot\boldsymbol{\Sigma},$ (42) $\displaystyle\mathbf{L}_{1}$
$\displaystyle=$ $\displaystyle(X_{0}+\Delta_{0}^{2}|{\bf
d}|^{2})\boldsymbol{\Sigma}+\Delta_{0}^{2}b_{0}{\bf
q}-\Delta_{0}^{2}\boldsymbol{\Sigma}\cdot{\bf d}{\bf
d}^{\ast}-\Delta_{0}^{2}\boldsymbol{\Sigma}\cdot{\bf d}^{\ast}{\bf d},$ (43)

Here $a_{0}=\tilde{\omega}-\xi$, $b_{0}=\tilde{\omega}+\xi$,
$X_{0}=b_{0}^{2}-\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}$ and the
denominator
$\widetilde{\mathbf{D}}=(\xi^{2}+\tilde{Q}_{+}^{2})(\xi^{2}+\tilde{Q}_{-}^{2})$,
where $\tilde{Q}_{\pm}$ is,

$\displaystyle\tilde{Q}_{\pm}^{2}$ $\displaystyle=$
$\displaystyle\Delta_{0}^{2}|{\bf
d}|^{2}-\tilde{\omega}^{2}-\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}\pm\sqrt{\Delta_{0}^{4}{\bf
q}\cdot{\bf
q}+4\tilde{\omega}^{2}\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}+4\Delta_{0}^{2}\tilde{\omega}{\bf
q}\cdot\boldsymbol{\Sigma}-4\Delta_{0}^{2}(\boldsymbol{\Sigma}\cdot{\bf
d})(\boldsymbol{\Sigma}\cdot{\bf d}^{\ast})}$ (44)

First we evaluate the trace term,

$\displaystyle\frac{1}{2}\mathrm{Tr}\left[\check{\mathcal{A}}({\bf
k},\omega)\check{\mathbf{j}}_{Q_{i}}\check{\mathcal{A}}({\bf
k},\omega)\check{\mathbf{j}}_{Q_{i}}\right]$ $\displaystyle=$ $\displaystyle
v_{Fi}^{2}\frac{1}{2}\mathrm{Tr}[T_{1}],$ $\displaystyle T_{1}$
$\displaystyle=$
$\displaystyle\hat{\mathcal{A}}_{11}\hat{\mathcal{A}}_{11}-\hat{\mathcal{A}}_{12}\hat{\mathcal{A}}_{21}-\hat{\mathcal{A}}_{21}\hat{\mathcal{A}}_{12}+\hat{\mathcal{A}}_{22}\hat{\mathcal{A}}_{22}.$
(45)

Note, the Green’s function have complex gap function, but the poles of the
Green’s function depend on the sign of $\tilde{\omega}$, which is renormalized
due to disorder. We can express $\check{\mathbf{G}}^{\prime\prime}$ as
$(\check{\mathbf{G}}-\check{\mathbf{G}}^{*})/2i$, and only a combination of
$\check{\mathbf{G}}$ and $\check{\mathbf{G}}^{*}$ will have poles on the
opposite side of the real axis to make a non-zero contribution to the integral
over $\xi$. For $\xi$ integration, it is useful to express the numerator and
the denominator in terms of polynomials. The numerator can be written as,

$\displaystyle g_{4}(\xi)$ $\displaystyle=$ $\displaystyle
b_{0}\xi^{6}+b_{1}\xi^{4}+b_{2}\xi^{2}+b_{3},$ (46) $\displaystyle b_{0}$
$\displaystyle=$ $\displaystyle 1,$ (47) $\displaystyle b_{1}$
$\displaystyle=$ $\displaystyle(|\tilde{\omega}|^{2}-\Delta_{0}^{2}|{\bf
d}|^{2}+|\boldsymbol{\Sigma}|^{2})+\mathrm{Re}\left[\tilde{Q}_{+}^{2}+\tilde{Q}^{2}_{-}\right]$
(48) $\displaystyle b_{2}$ $\displaystyle=$
$\displaystyle\frac{1}{4}|\tilde{Q}_{+}^{2}+\tilde{Q}^{2}_{-}|^{2}+(|\tilde{\omega}|^{2}+|\boldsymbol{\Sigma}|^{2}-\Delta_{0}^{2}|{\bf
d}|^{2})\mathrm{Re}\left[\tilde{Q}_{+}^{2}+\tilde{Q}^{2}_{-}\right]+3\Delta_{0}^{4}{\bf
q}\cdot{\bf
q}+4|\boldsymbol{\Sigma}|^{2}(|\tilde{\omega}|^{2}+\mathrm{Re}[\tilde{\omega}^{2}])$
(49) $\displaystyle-4\Delta_{0}^{2}\left(|\boldsymbol{\Sigma}\cdot{\bf
d}|^{2}+|\boldsymbol{\Sigma}\cdot{\bf
d}^{\ast}|^{2}+\mathrm{Re}[(\boldsymbol{\Sigma}\cdot{\bf
d})(\boldsymbol{\Sigma}\cdot{\bf d}^{\ast})]\right)$
$\displaystyle+4|\tilde{\omega}|^{2}\mathrm{Re}[\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}]+4\Delta_{0}^{2}\mathrm{Re}[\tilde{\omega}{\bf
q}\cdot\boldsymbol{\Sigma}^{\ast}+2\tilde{\omega}{\bf
q}\cdot\boldsymbol{\Sigma}]$ $\displaystyle b_{3}$ $\displaystyle=$
$\displaystyle(|\tilde{\omega}|^{2}-\Delta_{0}^{2}|{\bf
d}|^{2})\left[|\alpha_{+}|^{2}+\Delta_{0}^{4}{\bf q}\cdot{\bf
q}\right]-2\Delta_{0}^{4}{\bf q}\cdot{\bf
q}\mathrm{Re}[\alpha_{+}]+\mathcal{Y}(\boldsymbol{\Sigma})$ (50)
$\displaystyle{\mathcal{Y}}(\boldsymbol{\Sigma})$ $\displaystyle=$
$\displaystyle\left\\{|(\tilde{\omega}^{2}+\Delta_{0}^{2}|{\bf
d}|^{2}-\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma})|^{2}-\Delta_{0}^{4}{\bf
q}\cdot{\bf
q}\right\\}|\boldsymbol{\Sigma}|^{2}-2\Delta_{0}^{2}\mathrm{Re}[(\alpha_{-}-\alpha_{-}^{\ast})\tilde{\omega}{\bf
q}\cdot\boldsymbol{\Sigma}^{\ast}]$ (51)
$\displaystyle-4\Delta_{0}^{2}(|\tilde{\omega}|^{2}+|\boldsymbol{\Sigma}|^{2}-\Delta_{0}^{2}|{\bf
d}|^{2})\mathrm{Re}[\tilde{\omega}{\bf
q}\cdot\boldsymbol{\Sigma}^{\ast}]+2\Delta_{0}^{4}\mathrm{Re}\left[|{\bf
q}\cdot\boldsymbol{\Sigma}|^{2}-({\bf
q}\cdot\boldsymbol{\Sigma})^{2}+\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}{\bf
q}\cdot{\bf q}\right]$
$\displaystyle+2\Delta_{0}^{2}\mathrm{Re}\left[(|\tilde{\omega}|^{2}-\tilde{\omega}^{2})-(|\boldsymbol{\Sigma}|^{2}-\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma})\right](|\boldsymbol{\Sigma}\cdot{\bf
d}|^{2}+|\boldsymbol{\Sigma}\cdot{\bf d}^{\ast}|^{2})-4\Delta_{0}^{2}|{\bf
d}|^{2}|\tilde{\omega}|^{2}|\boldsymbol{\Sigma}|^{2}$
$\displaystyle+4\Delta_{0}^{2}\mathrm{Re}\left[\alpha_{+}(\boldsymbol{\Sigma}^{\ast}\cdot{\bf
d}^{\ast})(\boldsymbol{\Sigma}^{\ast}\cdot{\bf
d})\right]-4\Delta_{0}^{2}\mathrm{Re}[(\tilde{\omega}^{2}-\Delta_{0}^{2}|{\bf
d}|^{2})\boldsymbol{\Sigma}^{\ast}\cdot\boldsymbol{\Sigma}^{\ast}],$

where $\alpha_{\pm}=\tilde{\omega}^{2}-\Delta_{0}^{2}|{\bf
d}|^{2}\pm\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}$, and the denominator
can be expressed as,

$\displaystyle|\mathbf{D}(\xi)|^{2}$ $\displaystyle=$ $\displaystyle
h_{4}(\xi)h_{4}(-\xi),$ (52) $\displaystyle h_{4}(\xi)$ $\displaystyle=$
$\displaystyle(\xi-iQ_{+})(\xi-iQ_{+}^{*})(\xi-iQ_{-})(\xi-iQ_{-}^{*})$ (53)
$\displaystyle=$ $\displaystyle
a_{0}\xi^{4}+a_{1}\xi^{3}+a_{2}\xi^{2}+a_{3}\xi+a_{4},$ (54) $\displaystyle
a_{0}$ $\displaystyle=$ $\displaystyle 1,$ (55) $\displaystyle a_{1}$
$\displaystyle=$
$\displaystyle-i(Q_{+}+Q_{+}^{*}+Q_{-}+Q_{-}^{*})=-2i\mathrm{Re}(Q_{+}+Q_{-})=ic_{1},$
(56) $\displaystyle a_{2}$ $\displaystyle=$
$\displaystyle-|Q_{+}|^{2}-|Q_{-}|^{2}-4\mathrm{Re}(Q_{+})\mathrm{Re}(Q_{-})=-c_{2}$
(57) $\displaystyle a_{3}$ $\displaystyle=$ $\displaystyle
2i|Q_{+}|^{2}\mathrm{Re}(Q_{-})+2i|Q_{-}|^{2}\mathrm{Re}(Q_{+})=-ic_{3},$ (58)
$\displaystyle a_{4}$ $\displaystyle=$ $\displaystyle
Q_{+}Q_{+}^{*}Q_{-}Q_{-}^{*}=|Q_{+}|^{2}|Q_{-}|^{2}=c_{4}.$ (59)

Here $h_{4}(\xi)$ is defined to ensure poles on the upper half of the plane
for the retarded Green’s functions. Using(Gradshteyn and Ryzhik, 2007, Eq.
3.112.5)

$\displaystyle\int_{-\infty}^{\infty}dx\frac{g_{4}(x)}{h_{4}(x)h_{4}(-x)}$
$\displaystyle=$ $\displaystyle
i\pi\frac{b_{0}(-a_{1}a_{4}+a_{2}a_{3})-a_{0}a_{3}b_{1}+a_{0}a_{1}b_{2}+a_{0}b_{3}(a_{0}a_{3}-a_{1}a_{2})/a_{4}}{a_{0}(a_{0}a_{3}^{2}+a_{1}^{2}a_{4}-a_{1}a_{2}a_{3})},$
$\displaystyle=$ $\displaystyle
i\pi\frac{(-a_{1}a_{4}+a_{2}a_{3})-a_{3}b_{1}+a_{1}b_{2}+b_{3}(a_{3}-a_{1}a_{2})/a_{4}}{(a_{3}^{2}+a_{1}^{2}a_{4}-a_{1}a_{2}a_{3})},$
$\displaystyle=$
$\displaystyle\pi\frac{(-c_{1}c_{4}+c_{2}c_{3})+c_{3}b_{1}+c_{1}b_{2}+b_{3}(-c_{3}+c_{1}c_{2})/c_{4}}{(c_{3}^{2}+c_{1}^{2}c_{4}-c_{1}c_{2}c_{3})}$
(61)

Finally, we have the thermal conductivity,

$\frac{\kappa_{ii}}{T}=\int_{-\infty}^{\infty}d\omega\frac{\omega^{2}}{T^{2}}\left(-\frac{dn_{F}(\omega)}{d\omega}\right)\left\langle
N_{0}v_{Fi}^{2}\frac{(-c_{1}c_{4}+c_{2}c_{3})+c_{3}b_{1}+c_{1}b_{2}+b_{3}(-c_{3}+c_{1}c_{2})/c_{4}}{(c_{3}^{2}+c_{1}^{2}c_{4}-c_{1}c_{2}c_{3})}\right\rangle_{FS}$
(62)

##  Appendix B Gap equation and transition temperature

We consider a simple separable pairing potential

$\displaystyle\mathbf{V}_{\alpha\beta;\alpha^{\prime}\beta^{\prime}}({\bf
k},{\bf k}^{\prime})=-V(\mathbf{d}({\bf
k})\cdot\boldsymbol{\sigma}i\sigma_{y})_{\alpha\beta}(\mathbf{d}({\bf
k}^{\prime})\cdot\boldsymbol{\sigma}i\sigma_{y})^{\dagger}_{\alpha^{\prime}\beta^{\prime}}$
(63)

here $\alpha,~{}\beta,~{}\alpha^{\prime},~{}\beta^{\prime}$ denote spin and
$V(>0)$ is the strength of pairing. This pairing potential gives a single
transition temperature for the mixed two component order parameters. In
principle, it is possible to construct an effective pairing potential to get
more than one transition, however, such scenario is excluded in the current
discussion. The gap equation can be written as,

$\displaystyle\Delta_{0}$ $\displaystyle=$
$\displaystyle\pi\lambda_{pair}T\sum_{\omega_{n}}\left[\frac{1}{2}\left(\frac{1}{\tilde{\mathtt{Q}}_{+}}+\frac{1}{\tilde{\mathtt{Q}}_{-}}\right)|{\bf
d}|^{2}-\frac{\Delta_{0}^{2}{\bf q}\cdot{\bf q}-2\boldsymbol{\Sigma}\cdot{\bf
d}\boldsymbol{\Sigma}\cdot{\bf d}^{\ast}+2i\tilde{\omega}_{n}{\bf
q}\cdot\boldsymbol{\Sigma}}{\tilde{\mathtt{Q}}_{+}\tilde{\mathtt{Q}}_{-}\left(\tilde{\mathtt{Q}}_{+}+\tilde{\mathtt{Q}}_{-}\right)}\right]\Delta_{0}$
(64)

where $\lambda_{pair}=N_{0}V$ and $\tilde{\mathtt{Q}}_{\pm}$ is,

$\displaystyle\tilde{\mathtt{Q}}_{\pm}^{2}$ $\displaystyle=$
$\displaystyle\Delta_{0}^{2}|{\bf
d}|^{2}+\tilde{\omega}^{2}_{n}-\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}$
(65) $\displaystyle\pm\sqrt{\Delta_{0}^{4}{\bf q}\cdot{\bf
q}-4\tilde{\omega}^{2}_{n}\boldsymbol{\Sigma}\cdot\boldsymbol{\Sigma}+4i\Delta_{0}^{2}\tilde{\omega}_{n}{\bf
q}\cdot\boldsymbol{\Sigma}-4\Delta_{0}^{2}(\boldsymbol{\Sigma}\cdot{\bf
d})(\boldsymbol{\Sigma}\cdot{\bf d}^{\ast})}$

We solve the gap equation using Matsubara Green’s function. In Matsubara
Frequencies, $\boldsymbol{\Sigma}\propto i\Delta_{0}^{2}$, therefore it will
drop out near $T_{c}$ and the linearized gap equation becomes,

$\Delta_{0}=\lambda(2\pi
T)\sum_{\omega_{n}>0}^{\Omega_{c}}{\Delta_{0}}\left\langle|\mathbf{d}|^{2}\left(\frac{1}{\omega_{n}+\Gamma_{N}}\right)\right\rangle_{FS},$
(66)

where $\Gamma_{N}$ is the normal state scattering rate. Using
$\theta\equiv\langle|\mathbf{d}|^{2}\rangle_{FS}$, we can express $T_{c0}$,

$T_{c0}=1.13\Omega_{c}\exp\left(-\frac{1}{\lambda\theta}\right),$ (67)

and the critical value of $\Gamma_{N}$ to completely kill the
superconductivity as,

$\displaystyle\Gamma_{cric}$ $\displaystyle=$
$\displaystyle\frac{1}{1.13}T_{c0}=0.885T_{c0},$ (68)

which is same as the critical value of impurity scattering rate that kills
$d$-wave superconductivity. The general equation to determine $T_{c}$ is,

$\displaystyle\ln\left(\frac{T_{c}}{T_{c0}}\right)$ $\displaystyle=$
$\displaystyle\Psi\left[\frac{1}{2}\right]-\Psi\left[\frac{1}{2}+\frac{\Gamma_{N}}{2\pi
T_{c}}\right].$ (69)

Here $\Psi$ denotes the digamma function and we have assumed
$T_{c0},\Gamma_{N}\ll\Omega_{c}$. In the weak disorder limit ($\Gamma_{N}\ll
T_{c0}$),

$T_{c}\approx
T_{c0}\left[1-\frac{\pi}{4}\frac{\Gamma_{N}}{T_{c0}}+\frac{7\zeta(3)}{4\pi^{2}}\left(\frac{\Gamma_{N}}{T_{c0}}\right)^{2}\right]$
(70)

##  Appendix C Some additional results

###  C.1 Density of states: Nonunitary states on cylindrical Fermi surface

Fig. A1 shows the low energy density of state (DOS) for the $A_{1u}+irB_{1u}$
state and the $B_{2u}+irB_{3u}$ state for various values of the mixing
parameter $r$. The former is an antiferromagnetic (AF) nonunitary state, while
the latter is a ferromagnetic nonunitary state. For both these states, the
possibility of $\hat{z}$-axis nodes is prevented by the topology of the Fermi
surface, which is assumed to be open along the $\hat{z}$-axis.

Figure A1:  The DOS for the $A_{1u}+irB_{1u}$ state shown in panel (a) with
corresponding $\Delta_{-}$ structure’s in the $xz$ and $yz$ planes in panels
(b) and (c), respectively. The DOS for the $B_{2u}+irB_{3u}$ state is shown in
panel(d) alongside the $\Delta_{-}$ in the $xy$-plane in panel(e) and in the
$yz$-plane in panel(f).  Figure A2:  Zero temperature limit $\kappa
T_{c}/\kappa_{N}T$ for the $A_{1u}+irB_{3u}$ in panels (a) to (d) and for the
$B_{1u}+irB_{2u}$ state in panels (e) to (h) for various values of mixing
parameter $r$. The $s$-wave scattering phase shift is
$\theta_{s}=\tan^{-1}(5)$ for all the panels.

###  C.2 Thermal conductivity for the $A_{1u}+irB_{3u}$ and $B_{1u}+irB_{2u}$
states

Fig. A2 shows the zero temperature limit thermal conductivity as a function of
relative $T_{c}$ reduction $\delta t_{c}(\equiv 1-T_{c}/T_{c0})$ for the
$A_{1u}+irB_{3u}$ in panels (a) to (d) and for the $B_{1u}+irB_{2u}$ states in
panels (e) to (h). The $s$-wave scattering phase shift is $\tan^{-1}(5)$. The
temperature evolution for these states in shown in Fig. A3.

Figure A3:  Thermal conductivity normalized to its normal state value at
$T_{c}$ for the $A_{1u}+irB_{3u}$ shown as a function of temperature
normalized to $T_{c}$ for a scatterers for $\delta t_{c}=0.02$ from panels (a)
to (d) and for $\delta t_{c}=0.15$ in panels (e) to (h) for various values of
mixing parameter $r$. The normalized thermal conductivity for the
$B_{1u}+irB_{2u}$ shown as a function of temperature normalized to $T_{c}$ for
$\delta t_{c}=0.02$ from panels (i) to (l) and for $\delta t_{c}=0.15$ in
panels (m) to (p). The $s$-wave scattering phase shift is $\tan^{-1}(5)$.

Figure A4: Thermal conductivity normalized to its normal state value at
$T_{c}$ for the $A_{1u}+irB_{2u}$ shown as a function of temperature
normalized to $T_{c}$ for a scatterers for $\delta t_{c}=0.02$ from panels (a)
to (d) and for $\delta t_{c}=0.15$ in panels (e) to (h) for various values of
mixing parameter $r$. The normalized thermal conductivity for the
$B_{1u}+irB_{3u}$ shown as a function of temperature normalized to $T_{c}$ for
$\delta t_{c}=0.02$ from panels (i) to (l) and for $\delta t_{c}=0.15$ in
panels (m) to (p). The $s$-wave scattering phase shift is $\tan^{-1}(2)$.

###  C.3 Effect of weak sactterers on thermal conductivity

Fig. A4 shows the temperature dependence of thermal conductivity normalized to
its normal state value for the $A_{1u}+irB_{2u}$ state and for the
$B_{1u}+irB_{3u}$ state. The $s$-wave scattering phase shift is
$\tan^{-1}(2)$. For each state, thermal conductivity is calculated for $\delta
t_{c}=0.02$ and $\delta t_{c}=0.15$.

Figure A5: Panel (a) shows the DOS per spin for the $A_{1u}+irB_{1u}$ state
for different values of the mixing parameter $r$. The inset shows the energy
dependence of DOS at lower energies. The dashed line illustrates the DOS for
unitary $B_{3u}$ state has been shown for comparison. Panel (b) and (c) show
the variation of $\Delta_{-}$ for different values of $r$ in the $xy$ and $yz$
planes, respectively. Panel (d) shows the DOS per spin for the
$B_{1u}+irB_{3u}$ state, and the gap structures in the $xy$ and $xz$ planes
are shown in panels (e) and (f).

##  Appendix D Spherical Fermi surface

Figure A6:  Zero temperature limit $\kappa T_{c}/\kappa_{N}T$ for the
$A_{1u}+irB_{1u}$ state in panels (a) to (c) and for the $B_{2u}+irB_{3u}$
state in panels (d) to (f) for various values of mixing parameter $r$. The
Fermi surface is assumed to be a three-dimensional sphere. The number of nodes
$N_{n}$ is indicated for each case.

Figure A7:  Temperature dependence of $\kappa T_{c}/\kappa_{N}T$ for the
$A_{1u}+irB_{1u}$ state in panels (a) to (c) for $\delta t_{c}=0.02$ and for
for $\delta t_{c}=0.15$ in panels (d) to (f) for various values of mixing
parameter $r$. The Fermi surface is assumed to be a three-dimensional sphere
and $\tan\theta=5$. Temperature dependence of $\kappa T_{c}/\kappa_{N}T$ for
the $B_{2u}+irB_{3u}$ state in panels (g) to (i) for $\delta t_{c}=0.02$ and
for $\delta t_{c}=0.15$ in panels (j) to (l) for various values of mixing
parameter $r$.

Now we consider the two component states on a spherical Fermi surface. Fig.
A5(a) shows the DOS for the $A_{1u}+irB_{1u}$ state, which is an AF and non-
chiral state. The variation of $\Delta_{-}$ for this state is shown in the
Fig. A5 (b) and (c). In case of a spherical Fermi surface, all AF states are
non-chiral and all the FM nonunitary states are chiral. For $r<1$ , the
$A_{1u}+irB_{1u}$ state remains gapped with deep minima along the $\hat{x}$
and $\hat{y}$ directions in the $k_{z}=0$ plane, and for $r=1$, four nodes
appear along the $\hat{x}$ and $\hat{y}$ directions and show quadratic density
of states. For $r>1$, nodes split into eight along the $\pm\hat{z}$
directions, where the polar angle of nodal position becomes
$\tan\theta=\pm\sqrt{r^{2}-1}$. For large values of $r$, nodes move towards
the $\hat{z}$ axis. The DOS remains quadratic at low energies. The DOS remains
the same for other AF nonunitary states, but the location of minima and nodes
changes.

Next, we show the DOS of the FM nonunitary state $B_{1u}+irB_{3u}$ in Fig.
A5(d) with $\Delta_{-}$ for this state in panel (e) and (f). For $r<1$, there
are four point nodes in the $xz$ plane at $\cot\theta=\pm r/\sqrt{1-r^{2}}$
and for $r>1$, there are four nodes in the $xy$ plane at
$\tan\theta=\pm\sqrt{1-r^{2}}$. The DOS is similar to a state with linear line
nodes. For $r=1$, there is a pair of quadratic nodes along the $\hat{x}$-axis
with linear DOS. The DOS for other FM nonunitary states remains qualitatively
same.

Now, we look at the residual thermal conductivity for the two component states
on a spherical Fermi surface. Fig. A6 shows the normalized thermal
conductivity in the zero temperature limit for the $A_{1u}+irB_{1u}$ state in
panels (a) to (c). For the gapped state ($r<1$), this state shows complete
isotropic behavior and the $\hat{z}$-axis remains very closer to the
$xy$-plane thermal response. For $r=1$, this has nodes in the $xy$-plane,
therefore, as a function of disorder, $\hat{z}$-axis is weaker, but as the
disorder level increases, the $\hat{z}$-axis response dominate, because there
are minima along the $xz$ and $yz$ planes, and $xy$-plane also contains a
maxima. For $r>1$, there are quadruple pairs of point nodes in the $xz$ and
$yz$planes, therefore, $\hat{z}$-axis always dominates. For the
$A_{1u}+irB_{2u}$ and the $A_{1u}+irB_{3u}$ states,
$\hat{z}\leftrightarrow\hat{y}$ and $\hat{z}\leftrightarrow\hat{x}$,
respectively, in Fig. A6(a) to (c).

Next, we consider the $B_{2u}+irB_{3u}$ state, which is a chiral state on a
spherical Fermi surface. Fig. A6 (d) to (f) show the zero temperature limit
normalized thermal conductivity. For $r<1$, the nodes are located in the
$yz$-plane, therefore, the $\hat{x}$ direction remains smallest. For small
values of $r$, nodes are closer to the $\hat{y}$-axis, hence this direction
show the largest thermal conductivity. As $r$ increases and reaches closer to
unity, the thermal conductivity along the $\hat{z}$-axis starts to dominate.
For $r=1$, there are quadratic point nodes along the $\hat{z}$-axis, the
thermal conductivity is largest along the nodal directions. For $r>1$, nodes
are located in the $xz$-plane, and the thermal conductivity becomes smallest
along $\hat{y}$ direction. For values closer to unity, $\hat{z}$-direction
shows highest normalized thermal conductivity, but as the $r$ increases, the
nodes move towards $\hat{x}$-axis, which shows highest normalized thermal
conductivity.

Finally, we look at the temperature dependence of the normalized thermal
conductivity for the $A_{1u}+irB_{1u}$, which is shown in Fig. A7(a) to (f).
For the $A_{1u}+irB_{1u}$, the anisotropy as a function of temperature follows
the zero temperature trend. in Fig. A7(g) to (l) show the normalized thermal
conductivity for the $B_{2u}+irB_{3u}$ state, the zero temperature anisotropy
continues at low temperature, but at higher temperature the anisotropy
reverses as seen in the case of cylindrical Fermi surface due to the upper
quasiparticle energy band.

## References

  * Ran  _et al._ (2019a) S. Ran, C. Eckberg, Q.-P. Ding, Y. Furukawa, T. Metz, S. R. Saha, I.-L. Liu, M. Zic, H. Kim, J. Paglione, and N. P. Butch, Science 365, 684 (2019a).
  * Aoki  _et al._ (2019a) D. Aoki, K. Ishida, and J. Flouquet, Journal of the Physical Society of Japan 88, 022001 (2019a).
  * Aoki  _et al._ (2022a) D. Aoki, J.-P. Brison, J. Flouquet, K. Ishida, G. Knebel, Y. Tokunaga, and Y. Yanase, Journal of Physics: Condensed Matter 34, 243002 (2022a).
  * Fujibayashi  _et al._ (2022) H. Fujibayashi, G. Nakamine, K. Kinjo, S. Kitagawa, K. Ishida, Y. Tokunaga, H. Sakai, S. Kambe, A. Nakamura, Y. Shimizu, Y. Homma, D. Li, F. Honda, and D. Aoki, Journal of the Physical Society of Japan 91, 043705 (2022).
  * Matsumura  _et al._ (2023) H. Matsumura, H. Fujibayashi, K. Kinjo, S. Kitagawa, K. Ishida, Y. Tokunaga, H. Sakai, S. Kambe, A. Nakamura, Y. Shimizu, Y. Homma, D. Li, F. Honda, and D. Aoki, Journal of the Physical Society of Japan 92, 063701 (2023).
  * Aoki  _et al._ (2019b) D. Aoki, A. Nakamura, F. Honda, D. Li, Y. Homma, Y. Shimizu, Y. J. Sato, G. Knebel, J.-P. Brison, A. Pourret, D. Braithwaite, G. Lapertot, Q. Niu, M. Vališka, H. Harima, and J. Flouquet, Journal of the Physical Society of Japan 88, 043702 (2019b).
  * Metz  _et al._ (2019) T. Metz, S. Bae, S. Ran, I.-L. Liu, Y. S. Eo, W. T. Fuhrman, D. F. Agterberg, S. M. Anlage, N. P. Butch, and J. Paglione, Phys. Rev. B 100, 220504(R) (2019).
  * Hayes  _et al._ (2024) I. M. Hayes, T. E. Metz, C. E. Frank, S. R. Saha, N. P. Butch, V. Mishra, P. J. Hirschfeld, and J. Paglione,  , arXiv:2402.19353 (2024), 10.48550/arXiv.2402.19353, arXiv:2402.19353 .
  * Volovik and Gorkov (1985) G. E. Volovik and L. P. Gorkov, Zh. Eksp. Teor. Fiz. 88, 1412 (1985).
  * Ueda and Rice (1985) K. Ueda and T. M. Rice, Phys. Rev. B 31, 7114 (1985).
  * Blount (1985) E. I. Blount, Phys. Rev. B 32, 2935 (1985).
  * Ran  _et al._ (2019b) S. Ran, I.-L. Liu, Y. S. Eo, D. J. Campbell, P. M. Neves, W. T. Fuhrman, S. R. Saha, C. Eckberg, H. Kim, D. Graf, F. Balakirev, J. Singleton, J. Paglione, and N. P. Butch, Nature Physics 15, 1250 (2019b).
  * Wang  _et al._ (2020) Z. Wang, J. O. Rodriguez, L. Jiao, S. Howard, M. Graham, G. D. Gu, T. L. Hughes, D. K. Morr, and V. Madhavan, Science 367, 104 (2020).
  * Kitaev (2001) A. Y. Kitaev, Physics-Uspekhi 44, 131 (2001).
  * Wilczek (2009) F. Wilczek, Nature Physics 5, 614 (2009).
  * Ivanov (2001) D. A. Ivanov, Phys. Rev. Lett. 86, 268 (2001).
  * Chung  _et al._ (2011) S. B. Chung, H.-J. Zhang, X.-L. Qi, and S.-C. Zhang, Phys. Rev. B 84, 060510 (2011).
  * Hayes  _et al._ (2021) I. M. Hayes, D. S. Wei, T. Metz, J. Zhang, Y. S. Eo, S. Ran, S. R. Saha, J. Collini, N. P. Butch, D. F. Agterberg, A. Kapitulnik, and J. Paglione, Science 373, 797 (2021).
  * Rosa  _et al._ (2022) P. F. S. Rosa, A. Weiland, S. S. Fender, B. L. Scott, F. Ronning, J. D. Thompson, E. D. Bauer, and S. M. Thomas, Communications Materials 3, 33 (2022).
  * Sakai  _et al._ (2022) H. Sakai, P. Opletal, Y. Tokiwa, E. Yamamoto, Y. Tokunaga, S. Kambe, and Y. Haga, Phys. Rev. Mater. 6, 073401 (2022).
  * Ajeesh  _et al._ (2023) M. O. Ajeesh, M. Bordelon, C. Girod, S. Mishra, F. Ronning, E. D. Bauer, B. Maiorov, J. D. Thompson, P. F. S. Rosa, and S. M. Thomas, Phys. Rev. X 13, 041019 (2023).
  * Azari  _et al._ (2023) N. Azari, M. Yakovlev, N. Rye, S. R. Dunsiger, S. Sundar, M. M. Bordelon, S. M. Thomas, J. D. Thompson, P. F. S. Rosa, and J. E. Sonier, Phys. Rev. Lett. 131, 226504 (2023).
  * Theuss  _et al._ (2024) F. Theuss, A. Shragai, G. Grissonnanche, I. M. Hayes, S. R. Saha, Y. S. Eo, A. Suarez, T. Shishidou, N. P. Butch, J. Paglione, and B. J. Ramshaw, Nature Physics (2024), 10.1038/s41567-024-02493-1.
  * Ambegaokar and Griffin (1965) V. Ambegaokar and A. Griffin, Phys. Rev. 137, A1151 (1965).
  * Kadanoff and Martin (1963) L. P. Kadanoff and P. C. Martin, Annals of Physics 24, 419 (1963).
  * Schmitt-Rink  _et al._ (1986) S. Schmitt-Rink, K. Miyake, and C. M. Varma, Phys. Rev. Lett. 57, 2575 (1986).
  * Hirschfeld  _et al._ (1986) P. Hirschfeld, D. Vollhardt, and P. Wölfle, Solid State Communications 59, 111 (1986).
  * Hirschfeld  _et al._ (1988) P. J. Hirschfeld, P. Wölfle, and D. Einzel, Phys. Rev. B 37, 83 (1988).
  * Ishihara  _et al._ (2023) K. Ishihara, M. Roppongi, M. Kobayashi, K. Imamura, Y. Mizukami, H. Sakai, P. Opletal, Y. Tokiwa, Y. Haga, K. Hashimoto, and T. Shibauchi, Nature Communications 14, 2966 (2023).
  * Marchenko (1986) V. I. Marchenko, Zh. Eksp. Teor. Fiz. 93, 583 (1986).
  * Sigrist and Rice (1987) M. Sigrist and T. M. Rice, Z. Phys. B 68, 9 (1987).
  * Annett (1990) J. F. Annett, Advances in Physics 39, 83 (1990).
  * Fledderjohann and Hirschfeld (1995) A. Fledderjohann and P. Hirschfeld, Solid State Communications 94, 163 (1995).
  * Suetsugu  _et al._ (2024) S. Suetsugu, M. Shimomura, M. Kamimura, T. Asaba, H. Asaeda, Y. Kosuge, Y. Sekino, S. Ikemori, Y. Kasahara, Y. Kohsaka, M. Lee, Y. Yanase, H. Sakai, P. Opletal, Y. Tokiwa, Y. Haga, and Y. Matsuda, Science Advances 10, eadk3772 (2024).
  * Aoki  _et al._ (2022b) D. Aoki, H. Sakai, P. Opletal, Y. Tokiwa, J. Ishizuka, Y. Yanase, H. Harima, A. Nakamura, D. Li, Y. Homma, Y. Shimizu, G. Knebel, J. Flouquet, and Y. Haga, J. Phys. Soc. Jpn 91, 083704 (2022b).
  * Lee  _et al._ (2023) S. Lee, A. J. Woods, P. F. S. Rosa, S. M. Thomas, E. D. Bauer, S.-Z. Lin, and R. Movshovich, “Anisotropic field-induced changes in the superconducting order parameter of ute2,” (2023), arXiv:2310.04938 [cond-mat.supr-con] .
  * Klemens (1955) P. G. Klemens, Proc. Phys. Soc. London, Sect. A 68, 1113 (1955).
  * Bardeen  _et al._ (1959) J. Bardeen, G. Rickayzen, and L. Tewordt, Phys. Rev. 113, 982 (1959).
  * Uher (1990) C. Uher, Journal of Superconductivity 3, 337 (1990).
  * Gradshteyn and Ryzhik (2007) I. S. Gradshteyn and I. M. Ryzhik,  _Table of Integrals, Series, and Products_ (Academic Press, Cambridge, 2007).