# Investigation of Quasi-particle Relaxation in Strongly Disordered
Superconductor Resonators

Jie Hu1 ID , Jean-Marc Matin1 ID , Paul Nicaise1 ID , Faouzi Boussaha1 ID ,
Christine Chaumont1, Michel Piat2, Pham Viet Dung2 and Piercarlo Bonifacio3 ID
1 GEPI, Observatoire de Paris, Université PSL, CNRS, 75014 Paris, France  2
Université de Paris, CNRS, Astroparticule et Cosmologie, F-75013 Paris, France
3 GEPI, Observatoire de Paris, Université PSL, CNRS, 92195 Meudon, France
<EMAIL_ADDRESS>
###### Abstract

In this paper, we investigate the quasi-particle (QP) relaxation of strongly
disordered superconducting resonators under optical illumination at different
bath temperatures with the Rothwarf and Taylor equations and the gap-
broadening theory described by the Usadal equation. The analysis is validated
with various single-photon responses of Titanium Nitride (TiN) microwave
kinetic inductance detectors (MKIDs) under pulsed 405 nm laser illumination.
The QP relaxation in TiN is dominated by QPs with energy below the energy gap
smeared by the disorder, and its duration is still inversely proportional to
the QP density. The QP lifetime versus temperature can be fitted. The
relaxation of the resonator can be further modeled with QP diffusion. The
fitted QP diffusion coefficient of TiN is significantly smaller than expected.
Our result also shows a significant increase in QP generation efficiency as
the bath temperature increases.

  * December, 2023

##  1 Introduction

Strongly disordered superconductors have been widely used for superconducting
resonators for photon detection in astrophysics[1, 2, 3, 4] and qubit readout
in quantum circuits[5]. They are particularly favored for microwave kinetic
inductance detectors[6] (MKIDs) in the optical and near-infrared bands because
of their high quality and high resistivity for better optical efficiency[2, 7,
8].

Compared to aluminum superconducting resonators[9], MKIDs made of strongly
disordered superconductors show a different QP relaxation process under
optical illumination. First, the QP relaxation in these MKIDs shows a very
fast QP relaxation after photon absorption as the bath temperature
($T_{bath}$) decreases[2, 7, 8, 10]. Second, their resonance frequency shift
versus $T_{bath}$ is proportional to the optical power ($P_{opt}$) under
steady illumination[3, 11], while the resonance frequency of MKIDs made of
superconductor described by BCS theory[12] is proportional to[9]
$\sqrt{P_{opt}}$. Third, the response in the optical band is much smaller than
expected[8, 13, 14]. Fourth, the density of states (DoS) of the QPs below the
energy gap ($\Delta$) is not zero due to the energy gap being smeared by the
disorder[15]. Thus, many QPs with energy less than $\Delta$ would be in the
disordered superconductors. However, the response of these QPs to the optical
illumination remains unknown.

In this paper, we investigate the QP relaxation time in strongly disordered
superconductors with the single-photon response of TiN MKIDs at 405 nm. We
show that the QP relaxation in TiN under optical illumination at low
temperatures still follows the Rothwarf and Taylor (RT) equations[16] despite
the observed anomalous response. We have also obtained an analytical solution
of RT equations for superconductors with relatively slow relaxation at low
temperatures with various phonon-trapping factors. We further considered QP
diffusion in the superconductor and fitted the diffusion coefficient
($D_{qp}$) and the pair-breaking coefficient ($\eta$), showing there is a
significant increase in $\eta$ as the bath temperature increases, which can be
the reason for the increasing responsivity of TiN under illumination observed
elsewhere[13]. The fitted $D_{qp}$ is significantly smaller than expected,
which can be why the QPs relax quickly after the photon absorption.

MKIDs in the optical and near-infrared bands are usually high-quality lumped
superconducting resonators[7, 10, 17, 18, 19, 20] that are made of an
interdigitated capacitor (IDC) and a meander. Photons ($\hbar\omega>2\Delta$)
absorbed in the meander break Cooper pairs and generate QPs, recombining into
Cooper pairs and emitting phonons within the QP lifetime. The extra QPs
increase the kinetic inductance and lower the resonance frequency. Such a
response can be read out in the phase and amplitude of the resonator by a
probing tone at the resonance frequency.

##  2 Superconducting Gap Broadening

The gap broadening of disorder superconductors is formulated by the Usadal
equations as[15, 21, 22]

$\displaystyle
iE\sin\theta+\Delta\cos\theta-\alpha_{d}\sin\theta\cos\theta=0,$ (1)

where $\sin\theta$ and $\cos\theta$ are the quasi-classical disorder-averaged
Green’s functions. $E$ is the energy relative to the Fermi level. The QP DoS
can be obtained as $\rho_{qp}(E)=N_{0}\Re{\cos\theta}$. The case
$\alpha_{d}=0$ corresponds to the BCS DoS with a peak at $E=\Delta$. When
$\alpha_{d}>0$, the superconducting gap will be broadened, and there would be
nonzero DoS for $E<\Delta$. It means there are QPs with $E<\Delta$. The
complex superconductivity $\sigma=\sigma_{1}-i\sigma_{2}$ can be obtained via
Nam’s theory[23] as

$\begin{split}\frac{\sigma_{1}}{\sigma_{n}}=&\frac{1}{\hbar\omega}\int_{E_{g}-\hbar\omega}^{-E_{g}}g_{1}(E,E^{\prime})[1-2f(E^{\prime})]\differential
E\\\
&+\int_{E_{g}}^{\infty}g_{1}(E,E^{\prime})[f(E)-f(E^{\prime})]\differential
E\\\ \end{split}$ (2)
$\displaystyle\begin{split}\frac{\sigma_{2}}{\sigma_{n}}&=\frac{1}{\hbar\omega}\int_{E_{g}-\hbar\omega}^{\infty}g_{2}(E,E^{\prime})[1-2f(E^{\prime})]\differential
E\\\ &+\int_{E_{g}}^{\infty}g_{2}(E^{\prime},E)[1-2f(E)]\differential
E,\end{split}$ (3)

with $g_{1}$ and $g_{2}$ as

$\begin{split}g_{1}=&\Re{\cos\theta(E)}\Re{\cos\theta(E^{\prime}])}\\\
&+\Re{\cos\theta(E)}\Re{\cos\theta(E^{\prime})}\end{split}$
$\begin{split}g_{2}=&\Im{\cos\theta(E)}\Re{\cos\theta(E^{\prime})}\\\
&+\Im{i\sin\theta(E)}\Re{i\sin\theta(E^{\prime})}\end{split}$

$E^{\prime}=E+\hbar\omega$ and $f(E)$ the Fermi-Dirac distribution, and
$E_{g}$ is the effective energy gap where the DoS starts to be positive and
$\sigma_{n}$ is the normal-state conductivity, equations 1.3a or 1.3b in [23].

Thus, phonons with energy above $2E_{g}$ can break the cooper pairs in TiN but
with a lower probability.

##  3 Rothwarf and Taylor Equations

RT equations are a general phenomenological description of the dynamics
between QPs and phonons for various superconductors[24, 25, 26, 27, 28] as

$\displaystyle\begin{split}\frac{dn}{dt}&=I+\beta P-Rn^{2}\\\
\frac{dP}{dt}&=-\frac{\beta}{2}P+\frac{Rn^{2}}{2}-\gamma(P-P_{T}^{0})\end{split}$
(4)

where $n=n(t)$ and $P=P(t)$ are the concentrations of QPs and phonons in the
superconductor, $\beta=2/\tau_{B}$ is the pair-breaking rate caused by the
phonon absorption, $\tau_{B}$ is the time a phonon breaks a cooper pair, and
$R$ is the QP recombination rate with the creation of a phonon. $I$ is the
background QP generation rate per unit volume. $P_{T}^{0}$ is the
concentration of phonons in thermal equilibrium at temperature $T$ with $I=0$,
$\gamma=1/\tau_{es}$ is their decay rate $\tau_{es}$ is the time phonons takes
to escape to the substrate. $P_{T}^{0}=Rn_{T}^{2}/\beta$ and $n_{T}$ is the
thermally excited QP density as

$\displaystyle n_{T}(T)=4N_{0}\int_{0}^{\infty}f(E,T)\rho(E,T)dE,$ (5)

where $f(E,T)$ is the Fermi-Dirac distribution and $\rho(E,T)$ is the DoS of
the QPs. $N_{0}$ is the single spin density on the Fermi level. Factor 4
accounts for two possible electronic spins and the integration of DoS from
$-\infty$ to $\infty$. Here all the phonons with $\hbar\omega>2E_{g}$ are
taken into account for the broadened superconducting gap. And we still
consider that the number of the QP generated by a photon with energy $E_{ph}$
as

$\delta n=\frac{\eta E_{ph}}{\Delta V},$ (6)

with $\eta$ the pair-breaking efficiency and $V$ the volume of the
superconductor, as the QP DoS shows a maximum at $E=\Delta$.

By neglecting the spatial fluctuation of the QP density in TiN[29], the
steady-state ($dn/dt=0$ and $dP/dt=0$) QP density $n_{qp}$ and $P_{T}$ can be
obtained as[30] $n_{qp}=\sqrt{\Gamma
I/R+(n_{T})^{2}},\text{and}~{}P_{T}=I/2\gamma+P_{T}^{0},$ where
$\Gamma=(1+\beta/2\gamma)$ is the so-called phonon trapping factor[31].

Eq. (4) can be approximated as[26, 27]

$\displaystyle n(t)=n_{1}(t)+n_{2}(t)-n_{s}$ (7)

$n_{1}$ describes the dynamic for the superconductor to reach a quasi-
stationary state, in which the QPs density in the superconductor reaches the
maximum $n_{s}$ and starts to recombine back to Cooper pairs. $n_{2}$
describes the relatively slow procedure[27] ($\differential^{2}n/\differential
t^{2}\ll\differential n/\differential t$) in which the superconductor recovers
from the perturbation. At low temperatures ($T\ll T_{c}$, $T_{c}$ is the
critical temperature of the superconductor), $n_{1}$, $n_{2}$ and $n_{s}$ are
solved as follows[27]

$\displaystyle n_{1}(t)$
$\displaystyle=\frac{\beta}{R}\quantity[-\frac{1}{4}-\frac{1}{2\xi}+\frac{1}{\xi}\frac{1}{1-Kexp(-\beta
t/\xi)}]$ (8) $\displaystyle n_{2}(t)$
$\displaystyle=n_{qp}\quantity(\frac{1+\kappa e^{-t/\tau_{qp}}}{1-\kappa
e^{-t/\tau_{qp}}})$ (9) $\displaystyle n_{s}$
$\displaystyle=\frac{\beta}{4R}(\frac{2}{\xi}-1)$ (10)

with the QP relaxation time $\tau_{qp}$ as

$\displaystyle\tau_{qp}=\Gamma/(2n_{qp}R),$ (11)

and

$\displaystyle\frac{1}{\xi}$
$\displaystyle=\sqrt{\frac{1}{4}+\frac{2R}{\beta}\quantity(n_{0}+2P_{0})}$
(12) $\displaystyle K$
$\displaystyle=\frac{(4Rn_{0}/\beta+1)-2\xi^{-1}}{(4Rn_{0}/\beta+1)+2\xi^{-1}}$
(13)

$\kappa=(n_{s}-n_{qp})/(n_{s}+n_{qp})<1$ represents the magnitude of
perturbation of the QP density in the superconductor.

$n_{0}=n_{qp}+\delta n$ and $\delta n=\eta E_{ph}/V\Delta$ is the change of
the initial QP density with $E_{ph}$ the photon energy, $\eta$ the pair-
breaking efficiency and $V$ the volume of the superconductor. $P_{0}$ is the
initial phonon density as

$\displaystyle P_{0}=P_{T}+\frac{\delta P}{1+2\gamma/\beta}$ (14)

where $\delta P$ is the initial change of the phonon density. We introduce a
factor ($1+2\gamma/\beta$) to account for the phonon loss when QP density
reaches the maximum. $\delta P>0$ will lead to an increase in QP density after
photon absorption[26]. Eq. (9) shares the same expression by Wang et al[32]
(Eq. (2)), but with a different meaning of $\kappa$. Eq. (9)-(14) are valid
for different phonon trapping factors and are compared with the numerical
solution of Eq. (4), which are shown in Fig. 1. Here the case with $\delta
P=0$, $\delta P=(1-\eta)E_{ph}/(2V\Delta$) with $\eta=0.57$ are compared,
which is a typical value for Cooper-pair breaking detectors [33], primarily
constrained by the fact that phonons will not further break Cooper pairs when
their energy is less than the binding energy of the Cooper pairs.. The
situation $\gamma/\beta=0$, which corresponds to the original condition[27]
for Eq. (7) to (11) is also included for comparison. This condition is usually
not met for superconductors at extremely low temperatures as $\tau_{B}$ is
relatively long. The QP density increases after photon absorption for $\delta
P>0$, as observed in $\text{MgB}_{2}$[26].

Figure 1: Analytical solution of Eq. (4) at low temperatures. ($--$),
($-\cdot$) and ($\cdot\cdot\cdot\cdot\cdot$) corresponds to $\delta P=0$,
$\delta P=(1-\eta)E/2V\Delta$ and $P_{0}=P_{T}+\delta P$ after initial photon
absorption. $P_{0}=P_{T}+\delta P$ corresponds to $\gamma/\beta=0$, which is
referred as strong bottleneck condition[27]. This condition is usually not met
for superconductors at extremely low temperatures.

Here, we set $\delta P=0$ for simplicity mainly for two reasons. First, after
the superconductor absorbs a photon in the optical band, the created phonons
with energy more significant than the Debye energy will escape to the
substrate immediately[18]. Second, the response time of MKIDs
$\tau_{res}=Q/\pi f_{r}$ is on the order of $1~{}\mu\text{s}$, where $Q$ is
the quality of the resonator, $f_{r}$ is the resonance frequency. Thus, it is
difficult to determine the value of $\delta P$ from pulse relaxation. A
separate measurement of a TiN film with femtosecond spectroscopy[26] is
necessary to determine the value of $\delta P$, the recombination rate $R$,
and the pair breaking rate $\beta$.

##  4 MKIDs Transient Response

As the current distribution in the meander is almost uniform, the first-order
time-dependent phase response of the MKID $\delta\phi$ relative to the center
of its resonance circle on the IQ plane is the convolution[34, 35, 36] between
the detector response and the QP’s temporal and spatial relaxation $n(x,t)$,
which is

$\displaystyle\delta\phi$
$\displaystyle=\int_{0}^{t}\int_{V}\frac{\phi_{0}}{\tau_{res}}e^{-(t-\tau)/\tau_{res}}\delta
n(v,\tau)\differential v\differential\tau,$ (15) $\displaystyle\phi_{0}$
$\displaystyle=\frac{2\alpha
Q\hbar\omega}{\pi\Delta_{0}V}\cdot\frac{\differential\sigma_{2}}{\differential
n_{qp}},$ (16)

where $\omega=2\pi f_{r}$ is the angular frequency, $\alpha$ is the kinetic
inductance fraction, $\Delta_{0}$ is the energy gap of the superconductor at
absolute zero, and $\sigma_{2}$ is the imaginary part of the complex
conductivity $\sigma=\sigma_{1}-j\sigma_{2}$ of the superconductor.
$\differential\sigma_{2}/\differential n_{qp}$ is estimated numerically with
Eq. (3) for TiN.

Table 1: Main Parameters for TiN MKIDs

| $T_{c}$ | d | $L_{k}$ |  $\rho_{n}$* | $Q_{i}$ | $\alpha_{d}$ | Q | $f_{r}$ | $\alpha$
---|---|---|---|---|---|---|---|---|---
| (K) | (nm) | (pH/$\square$) | ($\mu\Omega\cdot$cm) | (k) | ($\Delta_{0}$) | (k) | (GHz) |
$TiN_{1K}$ | 0.84 | 60 | 84 | 330 | 9.0 | 0.12 | 6.3 | 3.182 | $0.98$
$TiN_{2K}$ | 2.1 | 15 | 61 | 140 | 18 | 0.092 | 13.1 | 2.625 | 0.85
$TiN_{4K}$ | 4.0 | 25 | 30 | 210 | 49 | 0.132 | 8.7 | 4.009 | 0.60

  * •

*$\rho_{n}$ is measured at before the superconducting transition.

For sufficiently large $t$ (usually on the order of $10-20~{}\mu s$, depending
on the QP lifetime) after photon absorption, the QPs diffusion can be
neglected. By inserting Eq. (7) into Eq. (15), the temporal phase response of
the resonator solved analytically is expressed as

$\displaystyle\delta\phi=$ $\displaystyle
2n_{qp}\phi_{0}V\sum_{m=1}^{\infty}\frac{\tau_{qp}\kappa^{m}}{\tau_{qp}-m\tau_{res}}\quantity(e^{-mt/\tau_{qp}}-e^{-t/\tau_{res}})$
(17)

We neglected the contribution of $n_{1}(t)-n_{s}$ as the time scale is on the
order of nanoseconds, much shorter than $\tau_{res}$. For
$\tau_{res}\ll\tau_{qp}$, Eq. (17) can be further simplified as

$\displaystyle\delta\phi\approx 2n_{qp}\phi_{0}V\quantity(\frac{\kappa
e^{-t/\tau_{qp}}}{1-\kappa
e^{-t/\tau_{qp}}}-\frac{\kappa}{1-\kappa}e^{-t/\tau_{res}}).$ (18)

The first term in Eq. (18) is the same as the result obtained by Fyhrie et al.
[37] with an assumption that QP lifetime saturates at low temperatures[38].
Eq. (18) has been used to fit the relaxation for QP relaxation for different
MKIDs[37, 10]. The procedure to obtain Eq. (15) to (18) can be found in A

$\tau_{qp}$ fitted only by exponential decay will decrease when $\kappa$
increases. For $\kappa\rightarrow 1$, the higher-order terms
$\kappa^{m}e^{-mt/\tau_{qp}}$ in Eq. (17) become dominant because they have a
weight of $\kappa^{m}$ in the series and lead to faster relaxation of the
resonator with a relaxation time of $\tau_{qp}/m$. Thus, the QP relaxation
shows a second relaxation at low temperatures when the QP density in the
superconductor is low. The resonator follows an exponential decay when
$\kappa$ decreases, as $\kappa\rightarrow 0$, only the first term in the
series in Eq. (17) is essential. In this case, the QP relaxation in MKIDs
approaches the exponential decay, and the relaxation time fitted from
exponential decay will be close $\tau_{qp}$ as the temperature rises.

##  5 MKIDs Characterization

We investigate the QP relaxation in TiN MKIDs with different critical
temperatures by fitting the single photon relaxation at $\lambda=405$ nm. The
detailed parameters for TiN MKIDs are listed in Table 1, where $d$ is the TiN
film thickness. The TiN films are deposited on sapphire substrates by
magnetron sputtering, which are then patterned into resonators using
photolithography and reactive ion etching. The films’ internal quality factor
$Q_{i}$ is relatively low as they are fabricated on sapphire[13, 14, 39],
which is probably due to the lattice mismatch on their interface. The MKIDs
design of $TiN_{1K}$ can be found in the B The design of $TiN_{2K}$ and
$TiN_{4K}$ can be found in our previous publication[13].

Figure 2: The 1-photon and the 0-photon responses averaged from the FWHM of
the corresponding peak in the pulse statistics of the measured MKID for
$P_{opt}=5$ dBm. The photon absorption in the meander is obtained by removing
the 0-photon event calculated by Eq. (19). The inset shows the pulse
statistics illuminated with two different optical powers.

The MKIDs are characterized in a pulse tube-pre-cooled adiabatic
demagnetization refrigerator (ADR)[40]. A niobium cylinder shields the stray
magnetic field and sheets of metglas 2714a around MKIDs. A pulsed 405 nm laser
illuminates the MKIDs through an optical fiber. The MKIDs are read by a
standard homodyne mixing scheme with a readout power of about 2 dB below
bifurcation. The pulse response of the MKID is sampled by an oscilloscope at
100 MHz and processed by a Wiener optimal filter to generate the photon-
counting statistics with the average pulse phase response as the template.
More details on the measurement setup can be found in C.

The measured height pulse statistics with two different optical illumination
powers are shown in the inset of Fig. 2 for $TiN_{1K}$, in which the first two
peaks correspond to 0-photon and 1-photon. We obtain 0-photon and 1-photon
responses by averaging the hundreds of events in the full-width-at-half-
maximum (FWHM) in the corresponding peaks. The 0-photon event corresponds to a
non-zero background phase response that is due to electrical crosstalk from
the surrounding pixels, phonons from the surrounding pixels through the
substrate, as well as the photon absorbed in the Nb ground plane and the Nb
feedline, as its rise time is around 3.5 $\mu\text{s}$ and is much longer than
$\tau_{res}\approx 0.6~{}\mu\text{s}$. The pulse statistics shift towards the
right, which indicates a more significant background response. The 0-photon
can be reduced using a microlens array[2] above the meander or continuous
illumination[18]. We used a pulsed illumination as it better captures the
pulse response’s rising edge and the multi-photon events.

Due to the high resistivity of the TiN film, the diffusion length of the QPs
$l=\sqrt{D_{qp}\tau_{qp}}$ is much shorter than the total length of the
meander, as the QP diffusion coefficient[41] is
$D_{qp}=D_{N}\sqrt{2k_{B}T/\pi\Delta}$, with $D_{N}=\sigma_{n}/e^{2}N_{0}$,
$\sigma_{n}$ the normal state conductivity. Therefore, the 0-photon and the
photon absorption in the meander can be considered independent. The 1-photon
response ($\delta I_{1},~{}\delta Q_{1}$) can be decomposed as the summation
of the 0-photon response ($\delta I_{0}$, $\delta Q_{0}$) and the photon
absorption ($\delta I_{ph}$, $\delta Q_{ph}$) in the meander in the in-phase
and quadrature plane as

$\displaystyle\delta I_{1}=\delta I_{0}+\delta I_{ph},~{}~{}~{}\delta Q=\delta
Q_{0}+\delta Q_{ph}.$ (19)

After removing $\delta I_{0}$ and $\delta Q_{0}$, the phase response of the
photon absorption can be calculated, as shown in Fig. 2, which has a shorter
relaxation time than the 0-photon response and a very fast relaxation after
the photon absorption. And they are the same for the two optical powers.

Figure 3: (A): Measured frequency shift versus $T_{bath}/T_{c}$ for different
MKIDs as shown in Table 1. The dashed line is the fitted frequency shift with
the disorder. The dashed dot line ($-\cdot$) is calculated with the standard
Mattis Bardeen theory[42]. The inset shows the broadened DoS at $T=0.1T_{c}$.
The dashed dot line shows the DoS for the BCS superconductor. (B): Fitted QP
lifetime with Eq. (17) from the single-photon response of different MKIDs
versus $T_{bath}/T_{c}$. The dashed line is fitted with
$\tau_{qp}=\Gamma/2n_{qp}R$ versus temperature. The right axis shows the ratio
of the QPs above $\Delta$ in $n_{qp}$, calculated as
$4N_{0}\int_{\Delta}^{\infty}f(E,T)\rho(E,T)dE/n_{T}(T_{qp})$. $n_{T}(T_{qp})$
is calculated with (5).

We fit the resonance frequency shift $\delta f_{r}/f_{r}$ and the QP
relaxation time $\tau_{qp}$ versus $T_{bath}$ together with the
superconducting gap broadening with the disorder to obtain $\Gamma/R$,
$\alpha_{d}$, the measure of the disorder of the superconductor and a
parameter to calculate the DoS of the QPs in TiN, $I$, the background QP
generation rate in Eq.(̃4) and $T_{qp}$, the temperature of the QPs. The
resonance frequency shift based on the Mattis-Bardeen theory[42] is added for
comparison. The resonance frequency versus the $T_{bath}$ can be fitted as[43]
$\delta f_{r}/f_{r}=\alpha\delta\sigma_{2}/2\sigma_{2}$ for thin film, as
shown in Fig. 3-(A). The inset shows the broadened DoS at $T=T_{c}$. The QP
lifetime shown in Fig. 3-(B) is fitted with phase relaxation starting from the
half of the maximum of the pulse response[10] with Eq. (17). Eq. (11) is
fitted to the measured values of the QP lifetime versus $T_{bath}$.

We assume a Kapiza boundary[44] on the interface between the TiN and the
sapphire as

$\displaystyle P_{d}=\Sigma A_{s}\quantity(T_{qp}^{4}-T_{bath}^{4}),$ (20)

where $\Sigma$ is a material constant, $A_{s}$ is the area of the meander and
$P_{d}$ is the dissipated microwave power in the resonator due to the internal
loss[45]. $T_{bath}$ refers to the temperature of the sapphire substrate and
is considered the same as that of the thermometer. Here, we neglect the
temperature difference between QPs and phonons in the superconductor as we do
not observe significant QP lifetime saturation[46, 47]. $T_{qp}$ is expected
to be higher than phonon temperature in the superconductor when $T_{bath}$ is
low due to the thermal decoupling between the electrons and phonons. The
obtained $T_{qp}$ is about $120$ mK, $273$ mK, and $400$ mK for $TiN_{1K}$,
$TiN_{2K}$ and $TiN_{4K}$ respectively at $T_{bath}=50$ mK. Although the
phonon temperature is overestimated, the equivalent temperature of the QPs
remains reasonable[48, 46]. It can be seen that Fig. 3-(B) that the QPs with
energy less than $\Delta$ dominate, and we can also expect that when
$T_{bath}$ is low compared to $T_{c}$, the QP dynamic of TiN deviates from the
BCS theory as the ratio of QPs with energy less than $\Delta$ increases. The
detailed fitting procedure can be found in the D.

##  6 Quasi Particle Diffusion

We further analyze the QP diffusion by introducing the diffusion term in the
RT equations as[49]

$\displaystyle\begin{split}\frac{dn}{dt}&=D_{qp}\nabla^{2}n+I+\beta
P-Rn^{2}\\\
\frac{dP}{dt}&=D_{ph}\nabla^{2}P-\frac{\beta}{2}P+\frac{Rn^{2}}{2}-\gamma(P-P_{T}^{0})\end{split}$
(21)

$D_{qp}$ and $D_{ph}$ are the diffusion coefficients for the QPs and phonons.

It is observed that as long as the QP lifetime $\tau_{qp}$ is kept the same,
the solution of Eq. (4) will be the same for different $\Gamma$. Thus, with
the fitted $\Gamma/R$, $T_{qp}$, $I$ and $\alpha_{d}$, $n_{qp}$ can be
calculated. $R$ can be obtained if we assume a value for $\beta$ and $\eta$.
As we don’t observe high energy resolving power in our MKIDs, $\Gamma$ is set
to be on the order of 1. We further assume $D_{ph}=0$ the phonon diffusion is
negligible. In this case, by neglecting the process of 2D QP diffusion to 1D
after photon absorption, the measured pulse response of MKIDs can be fitted
with the result numerically obtained by solving Eq. (4) and (15) in 1D with
two fitting parameters, the quasi-particle diffusion coefficient $D_{qp}$ and
the pair breaking efficiency $\eta$ with a proper initial condition.

The initial QPs are assumed to be Gaussian distributed with $\sigma_{0}$ on
the order of the thickness of the film. The Gaussian distribution is taken as
the QP density should be continuous and symmetrical after the initial photon
absorption, and it is a good approximation of the Dirac function. $D_{qp}$ is
almost independent of $\sigma_{0}$, while $\eta$ shows a strong dependence as
the initial relaxation rate of QP after the photon absorption is proportional
to $n_{qp}^{2}$. The $\sigma_{0}$ is chosen to have $\eta$ smaller than 0.6.
We consider that $\tau_{qp}$ for $TiN_{4K}$ is too short and that the initial
2D QP diffusion is essential. The 1D diffusion model doesn’t hold for
$TiN_{4K}$. We fix[50] $N_{0}=6.0\times 10^{10}~{}eV^{-1}\mu\text{m}^{-3}$, a
value needed to be determined experimentally for our film by measuring the
diffusion constant in the normal state. The detailed fitting procedure can be
found in E.

Figure 4: (A): Measured and fitted resonator relaxation for $TiN_{1K}$ and
$TiN_{2K}$ at $T_{bath}=50$ mK with 1-D diffusive RT equations. The inset
shows the pulse in the log scale. (B): The fitted diffusion coefficient versus
normalized temperature. The solid line is calculated as
$D_{qp}=D_{N}\sqrt{2k_{B}T/\pi\Delta}$. The dashed line is calculated as
$D_{qp}=\frac{4}{n_{qp}}\int_{\Delta}^{\infty}D_{E}f(E)\rho_{qp}(E)dE$. (C)
The fitted photon to QP conversion efficiency $\eta$ at different
temperatures. The error bar corresponds to one standard deviation $\sigma$ and
is obtained when $\chi^{2}$ increases by about 2.3, corresponding to $\sigma$
for two parameter fitting.

The single photon response fitted with Eq. (4) and (15) is shown in Fig.
(4)-(A) for $TiN_{1K}$ and $TiN_{2K}$ at $T_{bath}=50$ mK, which shows that
the fast QP relaxation term can be well modeled.

The fitted $D_{qp}$ versus $T_{bath}$ is shown in Fig. (4)-(B). It is about
two orders smaller than the theoretical value for two reasons. First, we have
considered the QPs with energy less than $\Delta$ and still assume the QP
density follows the Fermi-Dirac distribution. The number of these QPs
significantly outnumbered the QPs with $E>\Delta$, as is shown in Fig. 3-(B).
The diffusion coefficient of the QPs is energy-dependent, such as for the BCS
superconductor[51],

$\displaystyle D_{E}=D_{N}\sqrt{1-(\Delta/E)^{2}}.$ (22)

Thus, the diffusion coefficient of the low-energy QPs should be significantly
lower than those with energy larger than $\Delta$. Second, there is a spatial
distribution of the $\Delta$ in the TiN[52], which makes the QPs more
difficult to diffuse. The $D_{qp}$ we obtained can be interpreted as the
average $D_{qp}$ over the QPs with different energy. As QP diffusion
dispersion is not known for TiN, Eq. (22) is taken for calculation. Then, the
averaged diffusion coefficient can be calculated as

$\displaystyle
D_{qp}=\frac{4}{n_{qp}}\int_{\Delta}^{\infty}D_{E}f(E)\rho_{qp}(E)dE,$ (23)

which is in qualitative agreement with the fitted QP diffusion coefficient.
Essentially, we treat the QPs with $E<\Delta$ as ”trapped” QPs, which are
consistent with the arguments for the anomalous response in TiN[7, 3].

A direct measurement of the QP diffusion coefficient of TiN at low
temperatures is needed to investigate the diffusion coefficient further. The
low diffusion coefficient can be the reason for the response of the MKIDs made
of disordered superconductors being small[13, 14].

We focus on how $\eta$ evolves with $T_{bath}$ in Fig. 4-(C). $\eta$ tends to
increase when $T_{bath}$ increases especially for $TiN_{2K}$. One of the
possible reasons for the rise of $\eta$ is due to the presence of the two-
level system (TLS) on the interface between the TiN and the sapphire, which
absorbs part of the photon energy[19]. As $T_{bath}$ increases, TLS saturates,
and then $\eta$ increases. Another possible reason is when the temperature is
low, the phonon generated by the recombination of the QPs with energy less
than $\Delta$ has a lower probability of breaking Cooper pairs. When
$T_{bath}$ increases, the percentage of QPs with energy above $\Delta$
increases, increasing $\eta$. Additionally, an increase in $T_{bath}$ boosts
the $D_{qp}$, resulting in more QPs being detected by the resonator. With the
$n^{2}$ term influencing the recombination rate in the RT equation, a higher
$D_{qp}$ translates to slower QP relaxation. In this case, $\eta$ will
increase as the temperature increases, consistent with the previous
publications[3, 13, 14], where an increase of responsivity of the MKIDs has
been observed when $T_{bath}$ rises. The fitted $\eta$ shows large uncertainty
as $\eta$ is mainly related to the maximum of the pulse, which usually shows
large uncertainty.

We show how $\eta$ and $D_{qp}$ affect the response of MKIDs in Fig. 5 and
Fig. 6. It can be seen that $D_{qp}$ significantly affects the pulse response.

The increase in $\eta$ versus $T_{bath}$ could be one of the reasons why the
frequency change of TiN MKIDs versus optical illumination is proportional to
$P_{opt}$ in the continuous illumination in the millimeter wave range. For the
continuous illumination, the change in QPs in the superconductor is $\delta
n_{qp}=\eta P_{opt}\tau_{qp}/(\Delta V)$. With Eq. (11),
$\tau_{qp}=\Gamma/[2(n_{qp}+\delta n_{qp})R]$. Then, $\delta n_{qp}$ can be
solved as

$\displaystyle\delta
n_{qp}=-\frac{n_{qp}}{2}+\frac{1}{2}\sqrt{n_{qp}^{2}+\frac{2\eta(P_{opt})P_{opt}\Gamma}{R\Delta
V}}.$ (24)

Here we put $\eta$ as a function of $P_{opt}$ to indicate its temperature
dependence. If $\eta$ increases with $P_{opt}$, the responsivity of MKIDs
would increase with $P_{opt}$, as is observed in previous publications[3, 11].
We consider it necessary to solve the full kinetic equations[53] to
investigate better the rising edge of the pulse response of MKIDs and better
understand the temperature dependence of $\eta$.

Figure 5: MKIDs response with different $\eta$ for $TiN_{2K}$ at $T_{bath}=50$
mK. The inset shows the maximum of the pulse versus $\eta$. $D_{qp}$ is fixed
with value shown in Fig. 4-(B) Figure 6: MKIDs response with different
$D_{qp}$ for $TiN_{2K}$ at $T_{bath}=50$ mK. The inset shows the maximum of
the pulse versus $D_{qp}$. $\eta$ is fixed with value shown in Fig. 4-(C)

##  7 Conclusion

In conclusion, we have modeled the single-photon relaxation process of TiN
MKIDs on sapphire with the diffusive Rothwarf and Taylor equations and the
gap-broadening theory based on the Usadal equation for disordered
superconductors. Our results qualitatively show the QPs with energy less than
$\Delta$ dominate the QP relaxation. The QP diffusion coefficient leads to
fast relaxation of QPs after the photon absorption observed in different
materials.

## Acknowledgement

The authors would like to thank Damien Prêle and Manuel Gonzalez from APC,
Université Paris Cité, for discussion about measurement, as well as Florent
Reix, Josiane Firminy, and Thibaut Vacelet from Paris Observatory for assembly
and mounting of the devices. The authors would also like to thank Dr. Eduard
Driessen from IRAM and Dr Wei-Tao Lv from the Department of Physics at the
Chinese University of Hong Kong for discussing disordered superconductors.
This work is supported by the European Research Council (ERC) through Grant
835087 (SPIAKID) and the UnivEarthS Labex program.

##  Appendix A MKIDs Responsivity

The response of MKIDs to excess QPs can be expressed as

$\displaystyle\differential\phi=\frac{\differential\phi}{\differential\omega_{r}}\cdot\frac{\differential\omega_{r}}{\differential
L}\cdot\frac{\differential
L}{\differential\sigma_{2}}\cdot\frac{\differential\sigma_{2}}{\differential
n_{qp}}\cdot\frac{\delta N_{qp}}{V}$ (25)

with[54] $\differential\phi/d\omega_{r}$ the change of phase due to the change
of the resonance frequency, and $\differential\omega_{r}/dL$ the resonance
frequency change due to the change of the inductance as

$\displaystyle\frac{\differential\phi}{\differential\omega_{r}}$
$\displaystyle=\frac{4Q}{\omega_{r}}$ (26)
$\displaystyle\frac{\differential\omega_{r}}{\differential L}$
$\displaystyle=-\frac{\alpha\omega_{r}}{2L_{k}}$ (27)

with $L=L_{g}+L_{k}$ the total inductance in the resonator, $\alpha=L_{k}/L$
the fraction of the kinetic inductance and $L_{g}$ is the geometrical
inductance. $\alpha$ is close to 1 for TiN. $\omega_{r}$ is the resonant
frequency. For a thin film superconductor, the kinetic inductance is
$L_{s}=1/(\sigma_{2}\omega d$) and $d$ is the film thickness[43], thus,

$\displaystyle\frac{\differential L}{\differential\sigma_{2}}$
$\displaystyle=-\frac{L_{k}}{\sigma_{2}}$ (28)

With[43] $\sigma_{2}/\sigma_{n}\approx\pi\Delta_{0}/(\hbar\omega)$ and
$\sigma_{n}$ the normal state conductivity, Eq. (25) can be expressed as

$\displaystyle\differential\phi$ $\displaystyle=\frac{2\alpha
Q\hbar\omega}{\pi\Delta_{0}\sigma_{n}}\cdot\frac{\differential\sigma_{2}}{\differential
n_{qp}}\cdot\frac{\delta N_{qp}}{V}$ (29)

Here, $\delta N_{qp}$ is the change of the QP number in the resonator, and it
can be expressed as the convolution with the detector response as

$\displaystyle\delta
N_{qp}=\int_{0}^{t}\int_{V}\frac{e^{-(t-\tau)/\tau_{res}}}{\tau_{res}}\delta
n(v,\tau)dxd\tau$ (30)

where V is the volume of the meander, $\tau_{res}=Q/\pi f_{r}$ is the response
time of the resonator. By substituting Eq. (29) into Eq. (30), Eq. (15) in the
main text can be obtained, as

$\displaystyle\delta\phi$
$\displaystyle=\int_{0}^{t}\int_{S}\frac{\phi_{0}}{\tau_{res}}e^{-(t-\tau)/\tau_{res}}\delta
n(s,\tau)\differential s\differential\tau.$ $\displaystyle\phi_{0}$
$\displaystyle=\frac{2\alpha
Q\hbar\omega}{\pi\Delta_{0}V}\cdot\frac{\differential\sigma_{2}}{\differential
n_{qp}}$

For TiN, we calculate the $\differential\sigma_{2}/\differential n_{qp}$
numerically with Eq. (3) in the main text. For a BCS superconductor,
$\differential\sigma_{2}/\differential n_{qp}$ can be calculated as

$\displaystyle\frac{\differential\sigma_{2}}{\differential n_{qp}}$
$\displaystyle=-\sigma_{n}\frac{\pi S_{2}}{2N_{0}\hbar\omega},$ (31)
$\displaystyle S_{2}$ $\displaystyle=1+\sqrt{\frac{2\Delta_{0}}{\pi
k_{B}T}}e^{-\xi}I_{0}(\xi)$ (32)

with $\xi=\hbar\omega/2k_{B}T$ and $I_{0}$ is the Bessel function of the first
kind.

For sufficiently large $t$ (usually on the order of 10-20 $\mu$s, depending on
the QP lifetime) after photon absorption, the spatial distribution of the
$\delta n(x,t)$ can be neglected. Thus, with Eq. (9), $\delta n(t)$ can be
expressed as

$\displaystyle\begin{split}\delta n(t)&=n_{qp}\left(\frac{1+\kappa
e^{-t/\tau_{qp}}}{1-\kappa e^{-t/\tau_{qp}}}\right)-n_{qp}\\\ &=\frac{2\kappa
n_{qp}e^{-t/\tau_{qp}}}{1-\kappa e^{-t/\tau_{qp}}}\\\
&=2n_{qp}\sum_{m=1}^{\infty}\kappa^{m}e^{-mt/\tau_{qp}}\end{split}$ (33)

By substituting Eq. (33) into Eq. (15), we obtain Eq. (17) as

$\displaystyle\begin{split}\delta\phi&=\frac{\phi_{0}V}{\tau_{res}}\int_{0}^{t}e^{-(t-\tau)/\tau_{res}}2n_{qp}\sum_{m=1}^{\infty}\kappa^{m}e^{-mt/\tau_{qp}}dt\\\
&=\frac{2n_{qp}V\phi_{0}}{\tau_{res}}\sum_{m=1}^{\infty}\kappa^{m}e^{-t/\tau_{res}}\int_{0}^{t}e^{-\frac{m\tau}{\tau_{qp}}+\frac{\tau}{\tau_{res}}}dt\\\
&=2n_{qp}V\phi_{0}\sum_{m=1}^{\infty}\kappa^{m}\frac{\tau_{qp}}{\tau_{qp}-m\tau_{res}}\left(e^{-mt/\tau_{qp}}-e^{-t/\tau_{res}}\right)\end{split}$

For $\tau_{qp}\gg\tau_{res}$, we obtain Eq. (18) as

$\displaystyle\begin{split}\delta\phi&\approx
2n_{qp}V\phi_{0}\sum_{m=1}^{\infty}\kappa^{m}\left(e^{-mt/\tau_{qp}}-e^{-t/\tau_{res}}\right)\\\
&=2n_{qp}V\phi_{0}\quantity(\frac{\kappa e^{-t/\tau_{qp}}}{1-\kappa
e^{-t/\tau_{qp}}}-\frac{\kappa}{1-\kappa}e^{-t/\tau_{res}})\end{split}$ (34)

##  Appendix B TiN1K MKIDs Design

The detailed design of MKIDs with $TiN_{1K}$ is shown in Fig. 7. The measured
TiN MKID is one of the first pixels in a 1000-pixel array. The meander line is
made of sputtered 60 nm TiN on sapphire with a film resistivity of about
$\rho_{n}=330~{}\mu\Omega\cdot\text{cm}$ ($R_{s}\approx
55.0~{}\Omega/\square$) and its size is $40.5\times 38.5~{}\mu\text{m}^{2}$
with a strip width of $2.5~{}\mu$m and a gap of $0.5~{}\mu$m. The total length
of the meander is about $520~{}\mu m$. The IDC is made of niobium, and its
size is about $120\times 90~{}\mu\text{m}^{2}$ with both strip width and
spacing of $1.2~{}\mu m$. The niobium film is 100 nm with a critical
temperature of around $9.2~{}K$ and a kinetic inductance of about
$0.2~{}pH/\square$.

Figure 7: MKIDs design for $TiN_{1K}$. The meander is made of TiN. The rest of
the resonator and the feedline are made of Niobium.

##  Appendix C Measurement Setup

The detailed measurement setup for the MKIDs is shown in Fig. 8. The MKIDs are
characterized in a two-stage pulse tube pre-cooled adiabatic demagnetization
refrigerator (ADR) at Laboratoire Astropaticule & Cosmologie (APC). The stray
magnetic field is shielded by a niobium cylinder that is 1.5 mm thick and
sheets of metglas 2714a around MKIDs. The niobium cylinder is thermalized on
the 1 K stage and encloses the entire 0.1 K stage. The MKIDs are read by a
standard homodyne mixing scheme. The input signal is generated by a signal
generator (SMA100A). It is first attenuated by a programmable attenuator
(RCDAT-8000-60) at room temperature and then attenuated 20 dB, 10 dB, and 20
dB (XMA-2082 series) on 4 K, 1 K, and 100 mK before being fed into MKIDs. The
output signal from MKIDs is first amplified by an LNA (LNF-LNC0.3_14B) on the
4K stage and amplified further by two room temperature amplifiers
(ZVA-183W-S+). The signal is down-converted to DC by an IQ mixer (AD(0.4-6.0))
and then sampled by an oscilloscope (HDO6034). Two double DC blocks (SD3463)
are placed between the 4K and 1K stages to operate the heat switch in the
cryostat. The readout power is estimated to be around -100 dBm at the input of
the detector. The MKIDs array is illuminated by an optical fiber placed 35 mm
above the pixels. The 405 nm laser (LP405C1) is modulated by a 250 Hz pulse
from the pulse generator (TG5012A) of a width of 50 ns. The output power of
the laser is estimated to be a few pW outside the cryostat, attenuated by a
digital step attenuator (DD-100 series from OZ optics). The pulse response of
the MKID is sampled by an oscilloscope (HDO6034) at 100 MHz.

Figure 8: Detailed measurement Setup for the MKIDs

##  Appendix D Resonance frequency shift and QP lifetime fitting

We first fit the quality factor $Q$, resonant frequency $f_{r}$ and the
coupling quality factor $Q_{c}$ from the transmission of the resonator
$S_{21}$ as

$\displaystyle S_{21}=ae^{-2\pi
jf\tau}[1-\frac{Q/Q_{c}e^{j\phi_{0}}}{1+2jQ(f-f_{r})/f_{r}}]$ (35)

where $a$ is the amplitude of the transmission. $\tau$ is the time delay on
the cable. $\phi_{0}$ is the parameter to account for the mismatch of the
resonator to the circuit.

The QP lifetime is fitted from the pulse relaxation from the half of the
maximum using Eq. (17) in the main text, which we rewrite here as

$\displaystyle\delta\phi$
$\displaystyle=2Vn_{qp}\phi_{0}\tau_{qp}\sum_{m=1}^{\infty}\frac{\kappa^{m}}{\tau_{qp}-m\tau_{res}}\quantity(e^{-mt/\tau_{qp}}-e^{-t/\tau_{res}}).$

We show an example of fitting the QP lifetime $\tau_{qp}$ for $TiN_{2K}$ at
$T_{bath}=50~{}$mK using Eq. (17) in Fig. 9.

Figure 9: QP lifetime $\tau_{qp}$ fitting for $TiN_{2K}$ at $T_{bath}=50$ mK.

We assume the QPs and the phonons to be of the same temperature and assume a
Kapiza boundary on the interface between the TiN and the sapphire substrate as

$\displaystyle P_{d}=\Sigma A_{s}\quantity(T_{qp}^{4}-T_{bath}^{4}),$

where $\Sigma A_{s}$ is a material constant, and $P$ is related to the power
absorbed in the resonator[45] as

$\displaystyle
P_{d}=\eta_{r}\frac{P_{r}}{2}\frac{4Q^{2}}{Q_{i}Q_{c}}\frac{Q_{i}}{Q_{i,qp}}$
(36)

where $P_{r}$ is the readout power for the resonator, $\eta_{r}$ is the
fraction of the energy in the resonator that is dissipated, $Q_{i},Q_{c},Q$ is
the internal quality factor, the coupling quality factor, and the loaded
quality factor respectively. $Q=Q_{i}Q_{c}/(Q_{i}+Q_{c})$. We consider $Q_{i}$
of the MKIDs to be dominated by the QPs and $Q_{i,qp}=Q_{i}$.

Thus, $T_{qp}$ can be calculated as

$\displaystyle T_{qp}=(\frac{\eta_{r}}{\sum
A_{s}}\frac{P_{r}}{2}\frac{4Q^{2}}{Q_{i}Q_{c}}+T_{bath}^{4})^{1/4}$ (37)

We set $\eta_{r}/\Sigma A_{s}$ as a fitting parameter and used $Q_{i}$ and
$Q_{c}$ fitted from $S_{21}$ to obtain $T_{qp}$.

We followed the procedure in the thesis by Coumou[55] to calculate the density
of states of the QPs with the disorder $\alpha_{d}$. We rewrite the Usadal
equation Eq. 1 as

$\displaystyle iE\sin\theta+\Delta\cos\theta-\alpha_{d}\sin\theta\cos\theta=0$

To calculate $\Delta(T)$, the Matsubara representation of these equations is
used as

$\displaystyle\Delta\cos\theta_{n}-\omega_{n}\sin\theta_{n}-\alpha\sin\theta_{n}\cos\theta_{n}=0$
(38)

with $\theta_{n}=\theta(\omega_{n})$ and

$\displaystyle\Delta\ln(\frac{T_{c}}{T})=2\pi
k_{B}T\sum_{m=0}^{\infty}(\frac{\Delta}{\omega_{m}}-\sin\theta_{n})$ (39)

Where $\omega_{m}=(2m+1)\pi k_{B}T$ are the Matsubara frequencies and
$m=0,~{}1,~{}2,...$. Eq. (38) and Eq. (39) need to be solved iteratively for
all the frequencies until the convergence is reached. Once the $\Delta(T)$ is
solved, $\sin\theta(E)$ and $\cos\theta(E)$ can be obtained from Eq. 1. With
$\sin\theta(E)$ and $\cos\theta(E)$, the imagery part of complex conductivity
$\sigma=\sigma_{1}-j\sigma_{2}$ can be calculated with Eq. (3) as

$\displaystyle\begin{split}\frac{\sigma_{2}}{\sigma_{n}}&=\int_{E_{g}-\hbar\omega}^{\infty}g_{2}(E,E^{\prime})[1-2f(E^{\prime})]\differential
E\\\ &+\int_{E_{g}}^{\infty}g_{2}(E^{\prime},E)[1-2f(E)]\differential
E\end{split}$

with
$g_{2}=\Im{\cos\theta(E)}\Re{\cos\theta(E^{\prime})}+\Im{i\sin\theta(E)}\Re{i\sin\theta(E^{\prime})}$,
$E=E+\hbar\omega$ and $E_{g}$ is effective energy gap with the DoS starts to
be positive. With the calculated density of state
$\rho(E)=N_{0}Re(\cos(\theta(E)))$, the QP density in the superconductor can
be calculated with Eq. (5) at $T_{qp}$ as

$\displaystyle\begin{split}n_{T}(T_{qp})&=2N_{0}\int_{-\infty}^{\infty}f(E,T_{qp})\rho(E,T_{qp})dE,\\\
n_{qp}&=\sqrt{\frac{\Gamma}{R}I+n_{T}^{2}}\end{split}$ (40)

With Eq. (11), $\tau_{qp}=\Gamma/(2n_{qp}R)$, and $\delta
f/f=\alpha\delta\sigma_{2}/2\sigma_{2}$, $\tau_{qp}$ and the resonance
frequency shift can be fitted together with $F_{1}=\Gamma/R$,
$F_{2}=\eta_{r}/(\sum A_{s})$, $F_{3}=\alpha_{d}$ and $F_{4}=I$ as the fitting
parameters. We find $\gamma I/R$ is negligible.

##  Appendix E Pulse fitting

In this section, we describe how we fit the single photon pulse response of
MKIDs. We rewrite Eq.(21) for easy reference.

$\displaystyle\begin{split}\frac{dn}{dt}&=D_{qp}\nabla^{2}n+I+\beta
P-Rn^{2}\\\ \frac{dP}{dt}&=D_{ph}\nabla^{2}P-\beta
P/2+Rn^{2}/2-\gamma(P-P_{T}^{0})\end{split}$

To solve Eq.(21), we need to know the $D_{qp}$, $\beta$, $R$, $I$ and
$P_{T}^{0}$. From the fitting of $\tau_{qp}$ versus temperature, we have
obtained value $F_{1}=\Gamma/R$, which is

$\displaystyle F_{1}=(1+\frac{\beta}{2\gamma})/R$ (41)

As long as we keep the same $F_{1}$, the solution of Eq.(21) remains the same.
This is reasonable as once $F_{1}$ is fixed, the $\tau_{qp}$ will be fixed
with the same QP density.  We assume that $\tau_{es}=0.35$ ns and
$\tau_{B}=0.25$ ns to keep the $\beta/\gamma$ on the order of 1 for both
$TiN_{1K}$ and $TiN_{2K}$ as we don’t observe $R$ with our detector[10]. In
Fig. 10, we show when the $\beta/\gamma$ is kept the same, the solution is the
same

Figure 10: Comparison with simulation result with same $\beta/\gamma$ but
different $\beta$ and $\gamma$

$\tau_{es}$ strongly depends on the acoustic mismatch between the
superconductor and the substrate[56]. It is almost constant versus
temperature[57]. We list the measured values $\tau_{es}$ of different commonly
used superconductors for easy reference. For a film with a thickness of 60 nm,
it is estimated to be 0.48 ns for niobium[57], around 0.5 ns for NbN depending
on substrates[58], and 0.19 ns for aluminum on silicon[30]. $\tau_{B}$ is also
considered a constant as its change can be negligible in the measured
temperature range[56]. In this case, we obtain $\tau_{0}=248~{}\text{ns}$ and
$28~{}\text{ns}$ for $TiN_{1K}$ and $TiN_{2K}$ respectively, to be on the
order of the measured values[50], which is calculated as[46]

$\displaystyle\tau_{0}=(\frac{2\Delta}{k_{B}T_{c}})^{3}\frac{1}{4\Delta
N_{0}R}$ (42)

Here, we list the detailed parameters for fitting the pulse response in Table
2.

With the $P_{T}^{0}=Rn_{T}^{2}/\beta$ calculated from Eq. (5), $I$ fitted D,
and an assumption that $D_{ph}=0$, we can fit Eq.(21) for $D_{qp}$ and the
pair-breaking efficiency $\eta$ with a proper initial condition.

Table 2: Parameters for pulse fitting for $TiN_{1K}$ and $TiN_{2K}$ @
$T_{bath}=50~{}$mK

Parameters | $TiN_{1K}$ | $TiN_{2K}$ | comment
---|---|---|---
$T_{c}$ | 0.84 K | 2.1 K | The critical temperature
$\alpha_{d}$ | 0.12 | 0.092 | The measure of the disorder of the flim
$d$ | 60 nm | 15 nm | The film thickness
$\tau_{B}$a | 0.25 ns | 0.25 ns | The pair-breaking time
$f_{r}$ | 3.182 GHz | 2.625 GHz | The resonance frequency
$\tau_{0}$ | 248 ns | 28 ns | The electron and phonon interaction time
$\tau_{es}$a | 0.35 ns | 0.35 ns | The phonon escape time
$Q$ | 6.3k | 13.1k | The quality factor of the resonator.
$V$ | 78$~{}\mu\text{m}^{3}$ | 2.55 $\mu\text{m}^{3}$ | The volume of the inductor.
$R$ | $5.90\times 10^{-18}\text{s}^{-1}$ | $2.23\times 10^{-17}\text{s}^{-1}$ | The QP recombination rate
$t_{qp}$ | 121 mK | 273 mK | The equivalent temperature of the QPs
$\sigma_{0}$b | 100 nm | 60 nm | The variance of the initial QP distribution.

  * a

These two values are assumed to obtain the $\Gamma$. Then, by fitting the QP
lifetime versus $T_{bath}$ with Eq. (11), $R$ can be obtained.  $\tau_{es}$
for both films are kept the same for simplicity as we are not able to
determine $\tau_{es}$ in both from our measurement.

  * b

The value is obtained to set $\eta$ to be smaller than 0.6.

###  E.1 Initial Condition

The QPs generated by the absorbed photon in MKIDs is

$\displaystyle\delta N_{qp}=\frac{\eta E_{ph}}{\Delta}$ (43)

where $E_{ph}=3.06~{}\text{eV}$ for a 405 nm photon.

We assume the QPs are Gaussian-distributed with $\sigma_{0}$. Then, the
initial density of the QPs as

$\displaystyle\delta
N_{qp}=wd\int_{-l/2}^{l/2}\frac{1}{\sigma_{0}\sqrt{2\pi}}e^{-\frac{x^{2}}{2\sigma_{0}^{2}}}dx$
(44)

where $w$ is the width of the strip, $d$ is the film thickness and $l$ is the
length of the meander. The Gaussian distribution is taken for simplicity as
the QP density should be continuous and symmetrical after the initial photon
absorption, meanwhile it is a good approximation to the Dirac function.

###  E.2 1D-Diffusive RT equation

We solve the 1D-Diffusive RT equation for $n_{qp}$ and keep only $\eta$ and
$D_{qp}$ as the fitting parameter. We show the fitting $\chi^{2}$ for the
pulse $TiN_{2K}$ at $T_{bath}=50$ mK in Fig. 11.

Figure 11: Fitting $\chi^{2}$ for $TiN_{2K}$ at $T_{bath}=50$ mK

## References

## References

  * [1] Leduc H G, Bumble B, Day P K, Eom B H, Gao J, Golwala S, Mazin B A, McHugh S, Merrill A, Moore D C, Noroozian O, Turner A D and Zmuidzinas J 2010 Applied Physics Letters 97 102509 ISSN 0003-6951 URL http://aip.scitation.org/doi/abs/10.1063/1.3480420
  * [2] Mazin B A, Meeker S R, Strader M J, Szypryt P, Marsden D, van Eyken J C, Duggan G E, Walter A B, Ulbricht G, Johnson M, Bumble B, O’Brien K and Stoughton C 2013 Publications of the Astronomical Society of the Pacific 125 1348 URL https://dx.doi.org/10.1086/674013
  * [3] Bueno J, Coumou P C J J, Zheng G, de Visser P J, Klapwijk T M, Driessen E F C, Doyle S and Baselmans J J A 2014 Applied Physics Letters 105 192601 ISSN 0003-6951 1077-3118
  * [4] Zobrist N, Coiffard G, Bumble B, Swimmer N, Steiger S, Daal M, Collura G, Walter A B, Bockstiegel C, Fruitwala N, Lipartito I and Mazin B A 2019 Applied Physics Letters 115 213503 URL https://doi.org/10.1063/1.5127768
  * [5] Chang J B, Vissers M R, Córcoles A D, Sandberg M, Gao J, Abraham D W, Chow J M, Gambetta J M, Beth Rothwell M, Keefe G A, Steffen M and Pappas D P 2013 Applied Physics Letters 103 012602 ISSN 0003-6951 URL https://doi.org/10.1063/1.4813269
  * [6] Day P K, LeDuc H G, Mazin B A, Vayonakis A and Zmuidzinas J 2003 Nature 425 817–821 ISSN 1476-4687 URL https://doi.org/10.1038/nature02037
  * [7] Gao J, Vissers M R, Sandberg M O, da Silva F C S, Nam S W, Pappas D P, Wisbey D S, Langman E C, Meeker S R, Mazin B A, Leduc H G, Zmuidzinas J and Irwin K D 2012 Applied Physics Letters 101 142602 ISSN 0003-6951 1077-3118
  * [8] Guo W, Liu X, Wang Y, Wei Q, Wei L F, Hubmayr J, Fowler J, Ullom J, Vale L, Vissers M R and Gao J 2017 Applied Physics Letters 110 212601 ISSN 0003-6951 1077-3118
  * [9] de Visser P J, Baselmans J J, Bueno J, Llombart N and Klapwijk T M 2014 Nat Commun 5 3130 ISSN 2041-1723 (Electronic) 2041-1723 (Linking) URL https://www.ncbi.nlm.nih.gov/pubmed/24496036
  * [10] Zobrist N, Clay W H, Coiffard G, Daal M, Swimmer N, Day P and Mazin B A 2022 Physical Review Letters 129 017701
  * [11] Hubmayr J, Beall J, Becker D, Cho H M, Devlin M, Dober B, Groppi C, Hilton G C, Irwin K D, Li D, Mauskopf P, Pappas D P, Van Lanen J, Vissers M R, Wang Y, Wei L F and Gao J 2015 Applied Physics Letters 106 073505 ISSN 0003-6951 1077-3118
  * [12] Bardeen J, Cooper L N and Schrieffer J R 1957 Physical Review 108 1175–1204 URL http://link.aps.org/doi/10.1103/PhysRev.108.1175
  * [13] Boussaha F, Hu J, Nicaise P, Martin J M, Chaumont C, Dung P V, Firminy J, Reix F, Bonifacio P, Piat M and Geoffray H 2023 Applied Physics Letters 122 ISSN 0003-6951 1077-3118
  * [14] Jie H, Paul N, Faouzi B, Jean-Marc M, Christine C, Alexine M, Florent R, Josiane F, Thibaut V, Dung P V, Michel P, Elisabetta C and Piercarlo B 2024 Journal of Low Temperature Physics 214 113–124 ISSN 1573-7357 URL https://doi.org/10.1007/s10909-023-03018-5
  * [15] Driessen E F C, Coumou P C J J, Tromp R R, de Visser P J and Klapwijk T M 2012 Physical Review Letters 109 107003 URL http://link.aps.org/doi/10.1103/PhysRevLett.109.107003
  * [16] Rothwarf A and Taylor B N 1967 Physical Review Letters 19 27–30 ISSN 0031-9007
  * [17] Beldi S, Boussaha F, Hu J, Monfardini A, Traini A, Levy-Bertrand F, Chaumont C, Gonzales M, Firminy J, Reix F, Rosticher M, Mignot S, Piat M and Bonifacio P 2019 Opt Express 27 13319–13328 ISSN 1094-4087 (Electronic) 1094-4087 (Linking) URL https://www.ncbi.nlm.nih.gov/pubmed/31052858
  * [18] de Visser P J, de Rooij S A, Murugesan V, Thoen D J and Baselmans J J 2021 Phys. Rev. Applied 16(3) 034051 URL https://link.aps.org/doi/10.1103/PhysRevApplied.16.034051
  * [19] Hu J, Boussaha F, Martin J M, Nicaise P, Chaumont C, Beldi S, Piat M and Bonifacio P 2021 Applied Physics Letters 119 ISSN 0003-6951 1077-3118
  * [20] Nicaise P, Hu J, Martin J M, Beldi S, Chaumont C, Bonifacio P, Piat M, Geoffray H and Boussaha F 2022 Journal of Low Temperature Physics 209 1242–1248 ISSN 1573-7357 URL https://doi.org/10.1007/s10909-022-02789-7
  * [21] Coumou P C J J, Driessen E F C, Bueno J, Chapelier C and Klapwijk T M 2013 Physical Review B 88 ISSN 1098-0121 1550-235X
  * [22] Lyu W T, Zhi Q, Hu J, Li J and Shi S C 2023 Chinese Physics B URL http://iopscience.iop.org/article/10.1088/1674-1056/ad03dc
  * [23] Nam S B 1967 Physical Review 156 487–493 URL https://link.aps.org/doi/10.1103/PhysRev.156.487
  * [24] Twerenbold D 1986 Phys Rev B Condens Matter 34 7748–7759 ISSN 0163-1829 (Print) 0163-1829 (Linking) URL https://www.ncbi.nlm.nih.gov/pubmed/9939456
  * [25] Lindgren M, Currie M, Williams C, Hsiang T, Fauchet P, Sobolewski R, Moffat S, Hughes R, Preston J and Hegmann F 1999 Applied Physics Letters 74 853–855
  * [26] Demsar J, Averitt R D, Taylor A J, Kabanov V V, Kang W N, Kim H J, Choi E M and Lee S I 2003 Phys Rev Lett 91 267002 ISSN 0031-9007 (Print) 0031-9007 (Linking) URL https://www.ncbi.nlm.nih.gov/pubmed/14754080
  * [27] Kabanov V V, Demsar J and Mihailovic D 2005 Phys Rev Lett 95 147002 ISSN 0031-9007 (Print) 0031-9007 (Linking) URL https://www.ncbi.nlm.nih.gov/pubmed/16241687
  * [28] Beck M, Klammer M, Lang S, Leiderer P, Kabanov V V, Gol’Tsman G and Demsar J 2011 Physical Review Letters 107 177007
  * [29] Sacepe B, Chapelier C, Baturina T I, Vinokur V M, Baklanov M R and Sanquer M 2008 Phys Rev Lett 101 157006 ISSN 0031-9007 (Print) 0031-9007 (Linking) URL https://www.ncbi.nlm.nih.gov/pubmed/18999631
  * [30] de Rooij S 2020 Quasiparticle Dynamics in Optical MKIDs: Single Photon Response and Temperature Dependent Generation-Recombination Noise Master thesis TU Delft
  * [31] Wilson C M and Prober D E 2004 Physical Review B 69 ISSN 1098-0121 1550-235X
  * [32] Wang C, Gao Y Y, Pop I M, Vool U, Axline C, Brecht T, Heeres R W, Frunzio L, Devoret M H, Catelani G, Glazman L I and Schoelkopf R J 2014 Nat Commun 5 5836 ISSN 2041-1723 (Electronic) 2041-1723 (Linking) URL https://www.ncbi.nlm.nih.gov/pubmed/25518969
  * [33] Kozorezov A G, Volkov A F, Wigmore J K, Peacock A, Poelaert A and den Hartog R 2000 Physical Review B 61 11807–11819 pRB URL https://link.aps.org/doi/10.1103/PhysRevB.61.11807
  * [34] Zmuidzinas J 2012 Annual Review of Condensed Matter Physics 3 169–214 ISSN 1947-5454 URL http://dx.doi.org/10.1146/annurev-conmatphys-020911-125022
  * [35] Martinez M, Cardani L, Casali N, Cruciani A, Pettinari G and Vignati M 2019 Physical Review Applied 11 ISSN 2331-7019
  * [36] Kouwenhoven K, Fan D, Biancalani E, de Rooij S A, Karim T, Smith C S, Murugesan V, Thoen D J, Baselmans J J and de Visser P J 2022 arXiv preprint arXiv:2207.05534 URL https://arxiv.org/abs/2207.05534
  * [37] Fyhrie A, Zmuidzinas J, Glenn J, Day P, LeDuc H G and McKenney C 2018 Progress towards ultra sensitive kids for future far-infrared missions: a focus on recombination times vol 10708 ed Zmuidzinas J and Gao J R International Society for Optics and Photonics (SPIE) p 107083A URL https://doi.org/10.1117/12.2312867
  * [38] Barends R, Baselmans J J, Yates S J, Gao J R, Hovenier J N and Klapwijk T M 2008 Phys Rev Lett 100 257002 ISSN 0031-9007 (Print) 0031-9007 (Linking) URL https://www.ncbi.nlm.nih.gov/pubmed/18643694
  * [39] Vissers M R, Gao J, Wisbey D S, Hite D A, Tsuei C C, Corcoles A D, Steffen M and Pappas D P 2010 Applied Physics Letters 97 232509 ISSN 0003-6951 1077-3118
  * [40] Hu J, Salatino M, Traini A, Chaumont C, Boussaha F, Goupil C and Piat M 2020 Journal of Low Temperature Physics 199 355–361 ISSN 1573-7357 URL https://doi.org/10.1007/s10909-019-02313-4
  * [41] Narayanamurti V, Dynes R C, Hu P, Smith H and Brinkman W F 1978 Physical Review B 18 6041–6052 URL https://link.aps.org/doi/10.1103/PhysRevB.18.6041
  * [42] Mattis D C and Bardeen J 1958 Physical Review 111 412–417 URL http://link.aps.org/doi/10.1103/PhysRev.111.412
  * [43] Gao J, Zmuidzinas J, Vayonakis A, Day P, Mazin B and Leduc H 2008 Journal of Low Temperature Physics 151 557–563 ISSN 0022-2291 1573-7357
  * [44] Pollack G L 1969 Reviews of Modern Physics 41 48–81 rMP URL https://link.aps.org/doi/10.1103/RevModPhys.41.48
  * [45] de Visser P J, Goldie D J, Diener P, Withington S, Baselmans J J and Klapwijk T M 2014 Phys Rev Lett 112 047004 ISSN 1079-7114 (Electronic) 0031-9007 (Linking) URL https://www.ncbi.nlm.nih.gov/pubmed/24580483
  * [46] Guruswamy T 2018 Nonequilibrium behaviour and quasiparticle heating in thin film superconducting microwave resonators Ph.d. thesis University of Cambridge
  * [47] Thomas C N, Withington S, Sun Z, Skyrme T and Goldie D J 2020 New Journal of Physics 22 073028 ISSN 1367-2630
  * [48] Budoyo R P, Hertzberg J B, Ballard C J, Voigt K D, Kim Z, Anderson J R, Lobb C J and Wellstood F C 2016 Physical Review B 93 ISSN 2469-9950 2469-9969
  * [49] Rajauria S, Courtois H and Pannetier B 2009 Physical Review B 80 ISSN 1098-0121 1550-235X
  * [50] Kardakova A, Finkel M, Morozov D, Kovalyuk V, An P, Dunscombe C, Tarkhov M, Mauskopf P, Klapwijk T M and Goltsman G 2013 Applied Physics Letters 103 ISSN 0003-6951 1077-3118
  * [51] Martinis J M, Ansmann M and Aumentado J 2009 Phys Rev Lett 103 097002 ISSN 0031-9007 (Print) 0031-9007 (Linking) URL https://www.ncbi.nlm.nih.gov/pubmed/19792820
  * [52] Escoffier W, Chapelier C, Hadacek N and Villégier J C 2004 Phys. Rev. Lett. 93(21) 217005 URL https://link.aps.org/doi/10.1103/PhysRevLett.93.217005
  * [53] Vodolazov D Y 2017 Physical Review Applied 7 ISSN 2331-7019
  * [54] Mazin B A 2005 Microwave kinetic inductance detectors Ph. d thesis California Institute of Technology
  * [55] Coumou P C J J 2015 Electrodynamics of strongly disordered superconductors Ph. d thesis TU Delft, Delft University of Technology
  * [56] Kaplan S B, Chi C C, Langenberg D N, Chang J J, Jafarey S and Scalapino D J 1976 Physical Review B 14 4854–4873 URL https://link.aps.org/doi/10.1103/PhysRevB.14.4854
  * [57] Oktasendra F, Berdiyorov G, Mekki A and Maneval J 2016 IEEE Transactions on Applied Superconductivity 26 1–4 ISSN 1051-8223
  * [58] Sidorova M, Semenov A, Hübers H W, Ilin K, Siegel M, Charaev I, Moshkova M, Kaurova N, Goltsman G N and Zhang X 2020 Physical Review B 102 054501