# Magnetoconductance and photoresponse properties of disordered NbTiN films

M. Sidorova  Humboldt-Universität zu Berlin, Department of Physics, Newtonstr.
15, 12489 Berlin, Germany  German Aerospace Center (DLR), Institute of Optical
Sensor Systems, Rutherfordstr. 2, 12489 Berlin, Germany  A.D. Semenov  German
Aerospace Center (DLR), Institute of Optical Sensor Systems, Rutherfordstr. 2,
12489 Berlin, Germany  H.-W. Hübers  Humboldt-Universität zu Berlin,
Department of Physics, Newtonstr. 15, 12489 Berlin, Germany  German Aerospace
Center (DLR), Institute of Optical Sensor Systems, Rutherfordstr. 2, 12489
Berlin, Germany  S. Gyger  Department of Applied Physics, KTH Royal Institute
of Technology, SE-106 91 Stockholm, Sweden  S. Steinhauer  Department of
Applied Physics, KTH Royal Institute of Technology, SE-106 91 Stockholm,
Sweden  X. Zhang  State Key Laboratory of Functional Materials for
Informatics, Shanghai Institute of Microsystem and Information Technology,
Chinese Academy of Sciences (CAS), Shanghai 200050, China  CAS Center for
Excellence in Superconducting Electronics, Shanghai 200050, China  A.
Schilling  Physics Institute, University of Zürich, Winterthurerstrasse 190,
8057 Zürich, Switzerland

###### Abstract

We report on the study of phonon properties and electron-phonon coupling in
thin NbTiN films, which are intensively exploited in superconducting devices.
Studied NbTiN films with thicknesses less than 10 nm are disordered with
respect to electron transport, the Ioffe-Regel parameter of
$k_{F}l_{e}=2.5-3.0$ ($k_{F}$ is the Fermi wave vector and $l_{e}$ is the
electron mean free path), and inelastic electron-phonon interaction, the
product $q_{T}l_{e}\ll 1$ ($q_{T}$ is the wave vector of a thermal phonon). By
means of magnetoconductance and photoresponse techniques, we derived the
inelastic electron-phonon scattering rate $1/\tau_{e-ph}$ and determined sound
velocities and phonon heat capacities. In the temperature range from 12 to 20
K, the scattering rate varies with temperature as $1/\tau_{e-ph}\propto
T^{3.45\pm 0.05}$; its value extrapolated to 10 K amounts to approximately 16
ps. Making a comparative analysis of our films and other films used in
superconducting devices, such as polycrystalline granular NbN and amorphous
WSi, we found a systematic reduction of the sound velocity in all these films
by about 50% as compared to the corresponding bulk crystalline materials. A
corresponding increase in the phonon heat capacities in all these films is,
however, less than the Debye model predicts. We attribute these findings to
reduced film dimensionality and film morphology.

##  I Introduction

Thin disordered NbTiN films are intensively exploited in various
superconducting devices. They have become a prominent material for the
fabrication of scientific and commercial superconducting nanowire single-
photon detectors (SNSPDs) [1] and already a decade ago were considered as an
alternative to NbN films in hot-electron bolometers (HEBs) [2]. Both electron
transport properties and phonon-related properties of these films are of high
importance because they determine the performance metrics of practical
devices. For instance, the strength of electron-phonon coupling and the sound
velocity along with the heat capacity of phonons jointly control the relative
strength of electron heating and the cooling rate of electrons and, therefore,
determine the intermediate-frequency (IF) bandwidth in HEBs [3, 4] or impact
the timing jitter in SNSPDs [5]. Although in the literature, it was repeatedly
indicated that NbN-based HEBs demonstrate larger IF bandwidths than NbTiN-
based HEBs (see e.g. the review [2] and references therein), the explicit
reason for this experimental fact has never been clarified. A lack of both
theoretical and empirical description of properties of NbTiN and other films
used in superconducting devices is often caused by their complexity, i.e. high
level of disorder and reduced dimensionality with respect to different
physical phenomena. NbTiN is a compound of two well studied nitrides, NbN and
TiN, with very different electron-phonon (e-ph) scattering times
($\tau_{e-ph}=1.8$ ns for TiN [6] and $\sim 10$ ps for NbN [7] at 10 K). From
this perspective, it is interesting to understand whether the $\tau_{e-ph}$
can be finely tuned by changing the Nb/Ti fraction $x$ in NbxTi1-xN. The
matter is even more complicated since the superconducting properties of NbN
films alone such as transition temperature, resistivity, and critical current
density are strongly dependent on the film stoichiometry [8]. It has been
shown that SNSPDs based on amorphous superconductors, e.g. silicides such as
WSi [9] or MoSi [10], can show favorable properties for the detection of low-
energy infrared photons due to their lower superconducting transition
temperature compared to Nb-based materials. The explicit role, if any, of the
amorphous state in this improvement however remains unclear.

In this study, we derive various properties of NbTiN films at low
temperatures, in particular, the e-ph scattering time, the sound velocity, and
the phonon heat capacity, and carry out a comparative analysis with other
films, NbN and WSi, intensively used in various superconducting devices. We
show that, in all films, the sound velocity is reduced and the phonon heat
capacity is increased as compared to the corresponding bulk crystalline
materials. We suppose that this effect is exclusively controlled by
modification of acoustic phonons due to reduced film dimensionality and film
morphology.

##  II Experiment and results

NbTiN films with thicknesses 6 and 9 nm (samples $\\#$L135 and $\\#$L134,
respectively, in Table 1) were deposited on 270 nm-thick SiO2 layers thermally
grown on Si substrates. Details of the deposition process have been reported
elsewhere [11]. According to studies of films prepared under the same
sputtering conditions [12], the compositions of our NbxTi1-xN films had a
Nb/Ti fraction of approximately $x=0.6$ and polycrystalline granular structure
with diameters of most grains in the range from 4 to 5 nm.

###  II.1 Magneto-transport measurements

Table 1: Transport and superconducting parameters of studied here NbTiN films ($\\#$L134, $\\#$L135) and the parameters (taken from literature) of three other (NbN, WSi) superconducting films. $RRR$ is the ratio between resistances at 300 K and at 25 K. The values of $k_{F}$ and $l_{e}$ are computed from the carrier density $n_{e}$ determined from Hall measurements at 25 K. Material |  Sample |  $d$ |  RRR |  $T_{C0}$ |  $R_{SN}$ |  $D$ |  $N(0)$ |  $B_{C2}(0)$ |  $\xi(0)$ |  $n_{e}$ |  $k_{F}$ |  $l_{e}$ |  $k_{F}l_{e}$
---|---|---|---|---|---|---|---|---|---|---|---|---|---
|  |  (nm) |  |  (K) |  ($\Omega$/sq) |  (cm2/s) |  (J-1m-3) |  (T) |  (nm) |  (m-3) |  (nm-1) |  (nm) |
NbTiN |  L135 |  6.0 |  0.80 |  8.41 |  710.6 |  0.458 |  2.00$\times 10^{47}$ |  13.9 |  4.9 |  4.97$\times 10^{28}$ |  11.4 |  0.22 |  2.5
NbTiN |  L134 |  9.0 |  0.83 |  9.51 |  381.4 |  0.472 |  2.40$\times 10^{47}$ |  15.2 |  4.6 |  5.39$\times 10^{28}$ |  11.7 |  0.26 |  3.0
NbN[7] |  2259 |  5.0 |  0.79 |  10.74 |  529.5 |  0.474 |  3.11$\times 10^{47}$ |  17.4 |  4.4 |  - |  - |  - |  1.2111For these films, the Ioffe-Regel parameter was computed as $k_{F}l_{e}=3Dm_{e}/\hbar$ wth the free electron mass $m_{e}$.
NbN[7] |  A853 |  6.4 |  0.71 |  8.35 |  954.0 |  0.339 |  1.89$\times 10^{47}$ |  19.2 |  4.2 |  - |  - |  - |  0.9111For these films, the Ioffe-Regel parameter was computed as $k_{F}l_{e}=3Dm_{e}/\hbar$ wth the free electron mass $m_{e}$.
WSi[13] |  4 |  5.0 |  0.96 |  4.08 |  357.0 |  0.700 |  3.10$\times 10^{47}$ |  6.7 |  7.0 |  - |  - |  - |  1.8111For these films, the Ioffe-Regel parameter was computed as $k_{F}l_{e}=3Dm_{e}/\hbar$ wth the free electron mass $m_{e}$.

Resistance and Hall measurements were carried out by a standard four-probe
technique in a Physical Property Measurement System (manufactured by Quantum
Design) in magnetic fields applied perpendicular to the film surface. The van
der Pauw method was used to determine the square resistance $R_{S}$ and
eliminate the effect of the planar geometry for two-dimensional (2d)
specimens. Fig. 1 shows the zero-field $R_{S}$ as a function of temperature.
When temperature decreases from 300 K down to about 25 K, $R_{S}$ increases
and reaches the maximum value of 705.9 and 383.8 $\Omega$/sq for the film
$\\#$L135 and $\\#$L134, respectively, that is most likely due to Anderson
localization. In the vicinity of superconducting transition, where the
inequality $\ln(T/T_{C0})\ll 1$ holds, we fitted the $R_{S}(T)$ data with the
theory of fluctuation conductivity of Aslamazov and Larkin [14] and Maki and
Thompson [15, 16]. For 2d films, it is given by

$R_{S}(T)=\left[\dfrac{1}{R_{SN}}+A_{2d}\dfrac{1}{16}\dfrac{e^{2}}{\hbar\ln(T/T_{C0})}\right]^{-1}.$
(1)

Here, $\hbar$ is the reduced Planck constant and $e$ is the elementary charge.
The best-fit values of the normal-state square resistance $R_{SN}$ and the BCS
mean-field transition temperature $T_{C0}$ are listed in Table 1; the best-fit
values of $A_{2d}$ are 2.21 and 2.16 for the film $\\#$L135 and $\\#$L134,
respectively.

Figure 1: Square resistance vs temperature at zero magnetic field (semi-
logarithmic scale). Points: experimental data; curves: the best fits obtained
with Eq. 1 at temperatures, where $\ln(T/T_{C0})\ll 1$, and extrapolated to
higher temperatures.

An external magnetic field suppresses the superconducting transition as it is
shown in Fig. 2(a) for a one representative sample $\\#$L135. From these
measurements, we determined the temperature-dependent upper critical field,
$B_{C2}(T)$, shown in Fig. 2(b), as a preset field at the midpoint-transition
temperature $T_{C}$, i.e. at the temperature where $R_{S}=R_{SN}/2$. From the
best-linear fits of the experimental $B_{C2}(T)$ data at temperatures close to
$T_{C}$ ($T/T_{C}$ between 0.05 and 0.25), we extracted the slopes
$\text{d}B_{C2}/\text{d}T$. We further used them to compute the electron
diffusion coefficient as $D=4k_{B}/(\pi e)[\text{d}B_{C2}/\text{d}T]^{-1}$,
the zero-temperature upper critical field as
$B_{C2}(0)=-0.69T_{C0}\text{d}B_{C2}/\text{d}T$ for the dirty regime, and the
zero-temperature Ginzburg-Landau coherence length $\xi(0)=\sqrt{\Phi_{0}/(2\pi
B_{C2}(0))}$, where $k_{B}$ is the Boltzmann constant and $\Phi_{0}=h/(2e)$ is
the magnetic flux. Total electron density of states at the Fermi energy $N(0)$
was determined via the Einstein relation $N(0)=1/(De^{2}R_{SN}d)$ where $d$ is
the film thickness. All these parameters are listed in Table 1.

Figure 2: (a) Square resistance for the film $\\#$L135 vs temperature at a set
of fixed magnetic fields in steps of 0.5 T. The straight horizontal line drawn
at $R_{S}=R_{SN}/2$ defines the midpoint transition temperature for each
field. (b) Upper critical magnetic field as function of temperature. Points:
experimental data; lines: the best-linear fits, which are used to evaluate $D$
and $B_{C2}(0)$.

Figure 3: Hall resistance vs magnetic field at 25 K.

Magnetic field dependence of the Hall resistance $R_{xy}(B)$, i.e. the
resistance appearing transverse to the current flow, is shown in Fig. 3. It
was measured at 25 K, well above the superconducting transition, where
$R_{xy}$ changes linearly with $B$. From the slope of the $R_{xy}(B)$
dependence, we found the Hall coefficient as $R_{H}=R_{xy}d/B$ that provides
us the carrier density $n_{e}$ according to $R_{H}=1/(en_{e})$. Further,
knowing the $n_{e}$, we determined the Fermi wavevector as
$k_{F}=(3\pi^{2}n_{e})^{1/3}$, the elastic electron mean free path as
$l_{e}=\hbar k_{F}/(n_{e}e^{2}R_{SN}d)$, and, finally, the Ioffe-Regel
parameter $k_{F}l_{e}$. The latter indicates a slightly different level of
disorder in our films. All these parameters are listed in Table 1.

###  II.2 Magnetoconductance

We determined the e-ph scattering time for our NbTiN films by means of the
magnetoconductance technique described in [7]. We measured $R_{S}$, varying
the magnetic field from 0 to 9 T at a set of fixed temperatures from 9 to 20
K. From these data, we computed the dimensionless magnetoconductance as
$\delta
G(B,T)=2\pi^{2}\hbar/e^{2}\left(R_{S}^{-1}(B,T)-R_{S}^{-1}(0,T)\right)$. The
result is shown in Fig. 4 with points for the film $\\#$L135. The change in
the conductance induced by the magnetic field originates from quantum
interference effects and is described by the theory of quantum corrections to
conductivity (see Appendix A). We considered four quantum corrections given by
Eqs. (3-6), which account for the electron weak-localization effect (WL),
Aslamazov-Larkin (AL) and Maki-Thompson (MT) superconducting fluctuations, and
electronic density of states (DOS) fluctuations. At experimental fields below
$B_{C2}(0)$, the correction due to renormalization of the single-particle
diffusion coefficient is small [17, 18] and we have neglected it.

Figure 4: Dimensionless magnetoconductance vs magnetic field for a set of
fixed temperatures indicated in the legend for the film $\\#$L135. Points:
experimental data; black curves: the best fits with a sum of Eqs. (3-6).

We fitted the experimental $\delta G(B,T)$ data with a sum of Eqs. (3-6) using
three independent fitting parameters, namely the temperature-dependent
electron dephasing time $\tau_{\phi}(T)$, the temperature-independent spin-
orbit interaction time $\tau_{s.o.}$, and a constant $C^{*}$ allowing the
latter to take values between 1.0 to 3.0 according to [19, 20, 21, 22]. The
best fits, shown with black curves in Fig. 4 for the representative film
$\\#$L135, were obtained with $C^{*}\approx 1.0$, $\tau_{s.o.}=20\pm 5$ ps,
and a maximum $\tau_{\phi}(T)$ of about 4 ps. For both films, $\tau_{s.o.}$ is
of the order of the maximum $\tau_{\phi}$ that corresponds to the weak spin-
orbit interaction in our NbTiN films. The best-fit values of
$\tau_{\phi}^{-1}(T)$ are plotted in Fig. 5(a).

Figure 5: (a) Electron dephasing rate vs temperature. Points: the best-fit
values of $\tau_{\phi}^{-1}(T)$ obtained by fitting the experimental $\delta
G(B,T)$ data; solid curves: the best fits with Eq. (2). (b and c) Electron-
phonon scattering time vs temperature for the NbTiN film $\\#$L135 (b) and
$\\#$L134 (c) in a double-logarithmic scale. Points: values of
$\tau_{e-ph}(T)$ obtained from the experimental $\tau_{\phi}^{-1}$ data as
$\tau_{e-ph}=(\tau_{\phi}^{-1}-\tau_{e-e}^{-1}-\tau_{e-fl}^{-1})^{-1}$; solid
lines: the best fits with the SM model, Eqs. (7b,7c).

The experimental dependence $\tau_{\phi}^{-1}(T)$ is described with a sum of
rates affiliated with different phase-breaking scattering events in which
electrons are involved. In the absence of extrinsic phase-breaking sources,
these are electron-electron (e-e) [23], electron-phonon (e-ph) [24], and
electron-fluctuation (e-fl) [25] scattering events. Consequently,

$\tau_{\phi}^{-1}=\tau_{e-e}^{-1}+\tau_{e-ph}^{-1}+\tau_{e-fl}^{-1},$ (2a)
where
$\displaystyle\tau_{e-e}^{-1}=\frac{k_{B}T}{\hbar}\frac{1}{2C_{1}}\ln(C_{1}),$
(2b)
$\displaystyle\tau_{e-ph}^{-1}=\alpha_{e-ph}^{-1}\left(T/T_{C0}\right)^{n},$
(2c)
$\displaystyle\tau_{e-fl}^{-1}=\frac{k_{B}T}{\hbar}\frac{1}{2C_{1}}\frac{2\ln(2)}{\ln(T/T_{C0})+C_{2}}.$
(2d)

Here, $C_{1}=\pi\hbar/(R_{SN}e^{2})$ and
$C_{2}=4\ln(2)/[\sqrt{\ln(C_{1})^{2}+128C_{1}/\pi}-\ln(C_{1})]$. The
expression for the e-e dephasing rate, Eq. (2b), accounts for dephasing only
due to Nyquist noise, which dominates in our experimental temperature range
$T\ll\hbar/(k_{B}\tau_{e})\sim 10^{4}$ K ($\tau_{e}=l_{e}^{2}/(3D)$ equals 0.3
and 0.5 fs for $\\#$L135 and $\\#$L135, respectively). It is worth noting that
the expression for the e-ph dephasing rate, Eq. (2c), is valid in a relatively
small temperature range, where $n$ is constant, in general, $n$ varies with
the temperature. We fitted the experimental $\tau_{\phi}^{-1}(T)$ data shown
in Fig. 5(a) with Eq. (2) using two independent fitting parameters
$\alpha_{e-ph}$ and $n$. At $T_{C0}$, $\tau_{e-ph}$ equals $\alpha_{e-ph}$.
The best-fit values of $\alpha_{e-ph}$ and $n$ are listed in Table 2, together
with $\tau_{e-ph}$ values extrapolated to 10 K.

Table 2: Some of the electron and phonon-related parameters of the studied here NbTiN films ($\\#$L134 and $\\#$L135) and films from two other materials NbN and WSi, for comparison. The values of $u_{t}$ and $\rho$ are obtained from fits with the SM model. The values of $\tau_{esc}$ and $c_{e}/c_{ph}$ are obtained from fits with the 2T model. Material |  Sample |  $\alpha_{e-ph}$111At $T=T_{C0}$, $\tau_{e-ph}$ equals $\alpha_{e-ph}$. |  $n$222The exponent in the temperature dependence of $\tau_{e-ph}^{-1}\propto T^{n}$ |  $\tau_{e-ph}$(10 K)333The dephasing time due to e-ph scattering (identical with the single-particle e-ph scattering time [26]) extrapolated to 10 K. |  $\tau_{EP}(T_{C0})$444The e-ph energy relaxation time (proportional to the $\tau_{e-ph}$ with a proportionality coeff. computed with Eq. (11) in [27]). |  $\tau_{esc}$ |  $c_{e}/c_{ph}(T_{C0})$
---|---|---|---|---|---|---|---
|  |  (ps) |  |  (ps) |  (ps) |  (ps) |
NbTiN |  L135 |  31.2 $\pm$ 3.5 |  3.5 $\pm$ 0.1 |  16.9 |  4.9 |  52.6 $\pm$ 1.5 |  0.39 $\pm$ 0.04
NbTiN |  L134 |  17.7 $\pm$ 1.1 |  3.4 $\pm$ 0.1 |  15.0 |  2.9 |  79.6 $\pm$ 4.3 |  0.13 $\pm$ 0.01
NbN[7] |  2559 |  9.3 |  3.5 |  11.9 |  1.4 |  25.9 |  0.83 $\pm$ 0.18
NbN[7] |  A853 |  21.7 |  3.2 |  12.4 |  4.2 |  39.0 |  0.25 $\pm$ 0.03
WSi[13] |  4 |  66 |  3.0 |  4.5 |  14.9 |  - |  1.40 $\pm$ 0.30

Fig. 5(b,c) shows $\tau_{e-ph}$ obtained from the dephasing rates as
$\tau_{e-ph}=(\tau_{\phi}^{-1}-\tau_{e-e}^{-1}-\tau_{e-fl}^{-1})^{-1}$ in the
temperature range where $\tau_{e-ph}^{-1}$ dominates other dephasing rates. In
the dirty limit with respect to e-ph scattering characterized by the product
$q_{T}l_{e}\ll 1$, where $q_{T}=k_{B}T/\hbar u$ is the phonon wave vector and
$u$ is the sound velocity, e-ph scattering is described by a model developed
by Sergeev and Mitin in [24] (the SM model, Appendix B). The model provides
the single-particle e-ph scattering time, which is identical with
$\tau_{e-ph}$ [26] extracted by the magnetoconductance method. We fitted
experimental $\tau_{e-ph}$ data wit the SM model, Eq. (7), using two
independent fitting parameters: the sound velocity of transverse phonons
$u_{t}$ and the mass density $\rho$. We fixed other parameters such as
$k_{F}$, $l_{e}$, and $N(0)$, at their values given in Table. 1 and adopted
for $m_{e}$ the value of the free electron mass. We implemented the lattice
parameter $a_{0}=0.43$ nm [28, 29, 30], the sound velocity of longitudinal
phonons $u_{l}=2u_{t}$ that is approximately valid for a large variety of
materials, and $k\sim 1.0$ (Appendix B). The best-fit values of $u_{t}$ and
$\rho$ are listed in Table 3.

###  II.3 Photoresponse in the time domain

We studied the energy relaxation of nonequilibrium electrons in NbTiN
microbridges by means of the photoresponse technique. The microbridges were
fabricated from the two NbTiN films, $\\#$L135 and $\\#$L134, shaped to a
width of $\sim$90 $\mu$m and a length of $\sim$1 $\mu$m (along the current
flow) with tapered contacts toward the electrodes, which are 21 $\mu$m apart,
to reduce the current-crowding effect. The microbridges were installed in a
continuous flow cryostat, kept in the resistive state at an ambient
temperature slightly larger than $T_{C0}$, and exposed to light pulses via a
quartz window of the cryostat. Light was focused onto the microbridges into a
spot with a diameter larger than $\sim 100~{}\mu$m. The light pulses were
generated by a Ti:Sapphire laser at a wavelength of about 800 nm with a
repetition rate of 80 MHz and a subpicosecond duration. Biasing the
microbridges with a small direct current, we recorded their photoresponse
(amplified voltage transients) to light pulses in the time domain with a
sampling oscilloscope (Fig. 6). The overall bandwidth of our readout was 0.1 -
5 GHz.

Figure 6: Normalized voltage transient vs time for microbridges fabriacated
from NbTiN film $\\#$L135 (a) and $\\#$L134 (b). Solid curves: experimental
data; dashed curves: the best fits with the 2-T model, Eqs. (9-13) in [7].

We described the voltage transients, shown in Fig. 6, with the two-temperature
(2T) model using the formalism given by Eqs. (9-13) in [7]. This formalism
takes into account the effect of a finite readout bandwidth and the signal
ringing caused by impedance mismatch. The later is seen in Fig. 6(a). The 2T
model describes the evolution of electron and phonon effective temperatures
raised by excitation via two coupled time-dependent equations with three
independent parameters: the ratio between electron and phonon heat capacities
$c_{e}/c_{ph}$, the phonon escape time $\tau_{esc}$, and the e-ph energy
relaxation time $\tau_{EP}$. We note here that $\tau_{EP}$ differs from
$\tau_{e-ph}$ (both are in Table 2) by a coefficient
$\tau_{EP}/\tau_{e-ph}<1$, given by Eq. (11) in [27], which depends on the
exponent $n$, e.g., is about 0.6 for $n=2.0$ and 0.1 for $n=4.0$. Hence,
describing the photoresponse with the 2T model, one gets indirect calorimetric
information on heat capacities.

We fitted the experimental transients using $c_{e}/c_{ph}$ and $\tau_{esc}$ as
fitting parameters, along with the fixed values of $\tau_{EP}$. The results
are shown in Fig. 6. The best-fit values of $c_{e}/c_{ph}$ and $\tau_{esc}$,
together with computed values of $\tau_{EP}$, are listed in Table 2.

##  III Discussion

For the studied thin NbTiN films, we found that in the temperature range from
12 to 20 K the e-ph scattering time varies as $\tau_{e-ph}\propto T^{-n}$ with
$n$ between 3.4 and 3.5 and that at 10 K the extrapolated value,
$\tau_{e-ph}$(10 K), falls between 15 and 17 ps. NbTiN is a compound of two
well studied nitrides, NbN and TiN, which exhibit rather different scattering
times. The $\tau_{e-ph}\propto T^{-3.0}$ was found for thin TiN films in the
temperature range from 1.7 to 4.2 K [6] resulting in $\tau_{e-ph}$(10 K) $\sim
1.8$ ns (actually, the authors derived the e-ph energy relaxation time, i.e.,
$\tau_{EP}$ in our notation, which, along with the exponent $n$, we used to
estimate $\tau_{e-ph}$). For NbN films with thicknesses from 5 to 10 nm, $n$
spreading between 3.2 and 3.8 was found for the temperature range from 14 to
30 K along with the extrapolated value $\tau_{e-ph}$(10 K) in the range from
12 to 18 ps. [7]. Consequently, with respect to inelastic e-ph interaction,
our NbxTi1-xN films with the Nb/Ti fraction $x\approx 0.6$ are closer to NbN
than to TiN. This certainly indicates that the strength of e-ph coupling in
NbxTi1-xN films cannot be controlled by a fine change of $x$ in the range from
0.6 to 1.0.

Let us now discuss the $\tau_{e-ph}$ provided by the SM model (Appendix B).
Under the condition $q_{T}l_{e}\ll 1$ and the dominance of vibrating
scattering centers ($k\sim 1$) fulfilled for our films,
$\tau_{e-ph(l)}\gg\tau_{e-ph(t)}$ and, consequently, the total e-ph scattering
rate is controlled by the interaction of electrons with transverse phonons
$\tau_{e-ph}\equiv\tau_{e-ph(t)}$. Therefore, we limit the discussion to the
properties of transverse phonons. Interacting with a phonon, an electron
exchanges momentum $\sim q_{T}$ within an interaction region $\sim
1/q_{T}\propto u_{t}$, so that the smaller the phonon velocity, the smaller
the interaction region, and the shorter the interaction time. In fact,
$\tau_{e-ph(t)}\propto\rho u_{t}^{3}$ (Eq. 7c). With $\rho^{*}=7.5$ g/cm3 and
$u_{t}^{*}=4.7$ nm/ps for bulk crystalline NbTiN computed from first-
principles [29] (Table 3), the SM model provides two orders of magnitude
larger $\tau_{e-ph}$ than the experimental value obtained from
magnetoconductance (Table 2). The best-fit values of $\rho$ and $u_{t}$, which
provide experimetal values of $\tau_{e-ph}$ in the framework of the SM model
are listed in Table 3. The same procedure we applied to fit the $\tau_{e-ph}$
data obtained for thin WSi films by magnetoconductance measurements in [13]
(see Appendix C). The best-fit values of $u_{t}$ and $\rho$ for WSi are also
listed in Table 3.

Table 3: Phonon-related properties. $c_{ph}^{(*)}$ were computed as for 3d Debye phonons using corresponding sound velocities $u_{t}^{(*)}$. $c_{ph}^{exp}$ were obtained from $c_{e}/c_{ph}$ ratios found as the best fits of the 2T model to the photoresponse data taking $c_{e}$ predicted by the Drude model. Material |  Sample |  $u_{t}^{*}$111Computed from first principles for NbTiN [29], NbN [29], and W3Si [31] with crystal structures. |  $u_{l}^{*}$111Computed from first principles for NbTiN [29], NbN [29], and W3Si [31] with crystal structures. |  $\rho^{*}$111Computed from first principles for NbTiN [29], NbN [29], and W3Si [31] with crystal structures. |  $c_{ph}^{*}$ |  $u_{t}$222Best-fit values of the parameters in the SM model. |  $\rho$222Best-fit values of the parameters in the SM model. |  $c_{ph}$ |  $c_{ph}^{exp}$
---|---|---|---|---|---|---|---|---|---
|  |  (nm/ps) |  (nm/ps) |  (g/cm3) |  (J/Km3) |  (nm/ps) |  (g/cm3) |  (J/Km3) |  (J/Km3)
NbTiN |  L135 |  5.3-4.7 |  9.0-8.4 |  6.7-7.5 |  400-500 |  1.8 |  3.0 |  8800 |  2700
NbTiN |  L134 |  |  |  |  500-700 |  1.9 |  3.6 |  10800 |  11000
NbN[7] |  2259 |  4.4 |  8.0 |  8.2 |  1300 |  2.4 |  7.8 |  7800 |  2500
NbN[7] |  A853 |  |  |  |  600 |  2.2 |  5.2 |  4700 |  4000
WSi[13] |  4 |  3.0 |  5.4 |  15.8 |  200 |  2.1333Obtained in Appendix C. |  8.2333Obtained in Appendix C. |  630 |  560

Let us now estimate the phonon escape time using the acoustic mismatch model
[32]. First, we compute the transmission $\bar{\eta}\approx 0.2$ of the
film/substrate interface for phonons using the best-fit values of $u_{t}$ and
$\rho$ (Table 3 for our NbTiN films; for SiO2 substrate, we took them from
Table I in [32]). Further, we find $\tau_{esc}=4d/(\bar{\eta}\bar{u})\approx
60.6$ ps for the film $\\#$L135 and 81.9 ps for the film $\\#$L134, where
$\bar{u}\approx u_{t}$ is the weighted sound velocity. Hence, in NbTiN films
the phonon escape time scales with the film thickness as $\tau_{esc}$[ps]
$\approx 9\,d$[nm]. The computed values of $\tau_{esc}$ are in good agreement
with those found as best fits to the photoresponse data in the framework of
the 2T model (Table 2). This allows us to conclude that the SM model provides
reasonable values of $u_{t}$ and $\rho$.

Further, from the best-fit values of $c_{e}/c_{ph}$ ratios (Table 2) and
$c_{e}=\pi^{2}k_{B}^{2}N(0)T/3$ predicted by the Drude model, we find the
phonon heat capacities $c_{ph}^{exp}$, which are listed in Table 3. Here we
assumed that the Drude model predicts a correct value of $c_{e}$, relying on
the fact that our NbTiN films are three-dimensional (3d) with respect to
electron transport ($d\gg l_{e}$).

It is clearly seen from Table 3 that the best-fit values of $u_{t}$ for thin
NbTiN, NbN, and WSi films are systematically reduced as compared to those of
corresponding bulk crystalline materials (marked with the asterisk in Table 3)
by $(u_{t}^{*}-u_{t})/u_{t}^{*}\sim 40-60\%$. Moreover, for all these films,
both experimental (calorimetric) $c_{ph}^{exp}$ and $c_{ph}$ computed in the
framework of the Debye model for 3d phonons with the reduced sound velocities
as $c_{ph}=2/5\,\pi^{2}k_{B}(k_{B}T/\hbar)^{3}u_{av}^{-3}$, where
$u_{av}=[2/3u_{t}^{-3}+1/3u_{l}^{-3}]^{-1/3}\approx u_{t}$, are systematically
larger than the heat capacities $c_{ph}^{*}$ computed within the same
approximation with the sound velocities of corresponding bulk crystalline
materials.

There are three effects on the sound velocity and phonon heat capacity to be
discussed: (i) depletion of long-wavelength phonon states,(ii) phonon
softening, and (iii) film morphology (amorphous or granular structure).

At low temperatures, the wavelength of a thermal phonon $\lambda_{ph}\approx
2\pi\hbar u/(k_{B}T)$ becomes comparable or larger than the film thickness
that modifies the 3d Debye phonon spectrum. Phonon spectra in a metallic film
with a variable thickness on a semi-infinite substrate were computed and
compared with experimental data [33]. For phonons excited perpendicularly to
the film/substrate interface, by lowering the film thickness at low
temperatures, the authors found a strong depletion of long-wavelength (low-
energy) phonon states beyond the cut-off wavelength $\lambda_{ph}=2d$.
Qualitatively the depletion should not affect the sound velocity but only
cause a reduction of the phonon heat capacity; its impact should grow with the
increase in the ratio $\lambda_{ph}/d$. A similar but even stronger impact of
the depletion on the phonon heat capacity should occur in granular films where
the depletion is almost isotropic. It affects phonon states with all
directions of the wave vector. Even though for all the films in Table 3 the
ratio $\lambda_{ph}/d>1$, the depletion alone does not fully explain our
experimental results.

In a thin mono-crystalline film, phonon softening, i.e. the decrease of the
effective Debye temperature with respect to the bulk value, occurs due to
weakening of ion bounds at film surfaces [34]. Grain boundaries and defects
increase the relative amount of weak bounds and enhance phonon softening. [35,
36]. Clearly, the weakening of bounds leads to a reduction in the mean sound
velocity and a corresponding increase in the phonon heat capacity. The effect
is most pronounced in ultra-thin granular films due to a large surface-to-
volume ratio. We have to note here that there exists an alternative
explanation of this effect relying on the proximity effect between the film
and the oxide layers on the film surfaces [37].

As it is seen from the discussion above, the effect of the film morphology on
the sound velocity and phonon heat capacity is hard to separate from the
effect of the reduced dimensionality. The structure of both NbTiN films [12]
and NbN films [37, 38, 39] is polycrystalline, granular with the mean grain
size of the order of the film thickness or larger. Polycrystalline NbTiN and
NbN films are composites of crystalline grains with amorphous boundaries [40]
while WSi films are rather amorphous (verified by x-ray diffraction in 60 nm-
thick WSi films [41]) without pronounced grains. However, it has been recently
shown that 6 nm-thick WSi films exhibit a preferred orientation rather than an
amorphous state [42]. For longitudinal phonons, the amorphous phase is not
expected to reduce noticeably the sound velocity, while for transverse phonons
it does because of substantial weakening of shear modulus. A large reduction
of transverse sound velocities in amorphous materials, relative to the
crystalline materials, was observed, for instance, in [43]. This effect can be
considered as amorphous phonon softening. It is known [44] that for amorphous
materials it results in an increase in the phonon heat capacity just by the
factor ($c_{ph}\propto u_{t}^{-3}$) the sound velocity is reduced in an
amorphous state as compared to crystalline material.

Comparing calorimetric experimental values of $c_{ph}^{exp}$ with the $c_{ph}$
values computed with reduced sound velocities $u_{t}$ of the SM model (Table
3), we found that for NbN ($\\#$A853), WSi, and NbTiN ($\\#$L134) films
$c_{ph}^{exp}\approx c_{ph}$ that might be an indirect indication of their
amorphous state. Although the statement may sound conceivable for the WSi film
and for the NbN, $\\#$A853 film (see the discussion on this particular film in
[7]), it contradicts to the known polycrystalline structure of the NbTiN,
$\\#$L134 film. For other granular films (NbTiN, $\\#$L135 and NbN,
$\\#$A-2259), $c_{ph}^{exp}$ is noticeably smaller than $c_{ph}$. The most
plausible explanation of this difference is the depletion of long-wavelength
phonon modes in granules. We presume that like for NbN films grown on non-
heated substrates [37], the grain size in NbTiN films increases with the
increase of the film thickness. Hence, in the thicker NbTiN film ($\\#$L134,
$d$=9 nm), the phonon depletion is less pronounced than in the thinner film
($\\#$L135, $d$=6 nm) and the phonon heat capacity in the former is defined by
the phonon softening due to film surfaces. We suppose that in the ultra-thin
polycrystalline films analyzed here, the impact of the film morphology on the
phonon heat capacity depends on the grain size and may noticeably modify the
impact of film dimensionality. Another important observation is that the mass
densities of all thin films reported here are noticeably less than those of
corresponding crystalline materials. This is the direct consequence of the
film morphology since in an amorphous or polycrystalline film host atoms are
packed at a smaller density than in the corresponding crystalline material.

Finally, let us compare the properties of our NbTiN films and NbN films
studied in [7] (samples $\\#$2259 and $\\#$A853) that are relevant to SNSPD
and HEB practical devices. With respect to superconducting and transport
properties, together with e-ph scattering times, both materials are very
similar. The main difference is in their phonon escape times, which, in
NbTiN/SiO2 (this study) $\tau_{esc}$[ps] $\approx 9\,d$[nm] is about twice as
large as in NbN/SiO2 ($\tau_{esc}$[ps] $\approx 5\,d$[nm]) [7]. Besides the
timing jitter, the difference in $\tau_{esc}$ is not expected to somewhat
noticeably affect the performance of an SNSPD. However, this difference
severely affects the performance of an HEB. Specifically, its intermediate-
frequency (IF) bandwidth is controlled by the cooling rate of electrons, in
which the phonon escape rate plays a limiting role [3, 4]. A larger
$\tau_{esc}$ results in a smaller IF bandwidth. This explains why, despite all
the efforts, NbTiN-based HEBs repeatedly demonstrate smaller IF bandwidth than
NbN-based HEBs ([2]).

##  IV Conclusion

We have studied dephasing and energy relaxation of electrons in thin
superconducting NbTiN films on Si/SiO2 substrates with magnetoconductance and
photoresponse techniques, respectively, and compared our results with those
reported in the literature for NbN and WSi films. Our main findings are the
following:

(a) The studied NbTiN films are strongly disordered with respect to electron-
phonon scattering, the product $q_{T}l_{e}(T_{C})\ll 1$, and electron
transport, the Ioffe-Regel parameter $k_{F}l_{e}=2.5-3.0$. In the temperature
range from 12 to 20 K, their inelastic e-ph scattering rate varies as
$\tau_{e-ph}^{-1}\propto T^{n}$ with $n=3.4-3.5$ and amounts (extrapolated) to
15.0-16.9 ps at 10 K. The phonon escape time varies with the film thickness as
$\tau_{esc}$[ps] $\approx 9\,d$[nm].

(b) In all studied and considered films, we have found a systematic reduction
of the sound velocity by $\Delta u/u=40-60$ % as compared to the sound
velocities computed from the first principles for corresponding bulk
crystalline materials.

(c) For all films the experimental calorimetric heat capacities of phonons are
much larger than the heat capacities computed for corresponding crystalline
materials with the 3d Debye model. For amorphous films and relatively thick
polycrystalline films, calorimetric heat capacities are almost equal to those
computed with the reduced sound velocities, while for granular films they are
a few times less. The latter is most plausibly due to the depletion of low
energy phonon modes in granules.

## Acknowledgements:

Authors greatly acknowledge the support of V. Zwiller in the sample
preparation and stimulating discussions with K. Illin on the magnetic
properties and morphology of superconducting granular films.

##  Appendix A Quantum corrections

The theory of quantum corrections to the classical Drude conductance is
applicable for materials with $k_{F}l_{e}>1$. With respect to the
characteristic scales of the theory, i.e. the thermal coherence length
$L_{T}=\sqrt{2\pi\hbar D/(k_{B}T)}$ and the electron dephasing length
$L_{\phi}$, our NbTiN films belong to quasi-2d systems ($d<L_{T},L_{\phi}$).
Therefore, the analytical expressions we use to compute theoretical
(dimensionless) corrections correspond to the 2d limit. For weak spin-orbit
interaction and in the absence of magnetic scattering, the WL correction [45,
19] is given by

$\delta
G^{WL}(B,T)=\frac{3}{2}Y\left(\frac{B}{\frac{4}{3}B_{s.o.}+B_{\phi}}\right)-\frac{1}{2}Y\left(\frac{B}{B_{\phi}}\right).$
(3)

Here, the function $Y(x)=\ln(1/x)+\psi(1/2+1/x)$ is defined via the digamma
function $\psi(x)$. The characteristic fields are defined as
$B_{j}=\hbar/(4eD\tau_{j})$, where the indices ’s.o.’ and ’$\phi$’ stand for
spin-orbit and dephasing scattering, respectively. The AL correction [14] is
given by

$\delta
G^{AL}(B,T)=\frac{\pi^{2}}{2\ln(T/T_{C0})}\left\\{\frac{B_{C}}{B}\left(1-2\frac{B_{C}}{B}\\\
\left[\psi\left(1+\frac{B_{C}}{B}\right)-\psi\left(\frac{1}{2}+\frac{B_{C}}{B}\right)\right]\right)-\frac{1}{4}\right\\}.$
(4)

Here, $B_{C}=C^{*}\hbar/(4eD\tau_{GL})$, where the Ginzburg-Landau time
$\tau_{GL}=(\pi\hbar)/(8k_{B}T\ln(T/T_{C0}))$ represents the lifetime of
Cooper pairs. The MT correction [15, 16, 46, 47] is given by

$\delta
G^{MT}(B,T)=-\beta_{LSA}(T)\left[Y\left(\frac{B}{B_{\phi}}\right)-Y\left(\frac{B}{B_{C}}\right)\right],$
(5)

where the parameter $\beta_{LSA}(T)=2\pi
k_{B}T\hbar^{-1}(1/\tau_{GL}-1/\tau_{\phi})^{-1}$. Finally, the DOS correction
[17, 48] is given by

$\delta
G^{DOS}(B,T)=\frac{28\varsigma(3)}{\pi^{2}}Y\left(\frac{B}{B_{C}}\right),$ (6)

where $\varsigma(3)=1.202$ is the Riemann zeta function. The total theoretical
magnetoconductance, i.e. the sum of the four terms in Eqs. (3-6), is used to
describe the experimental $\delta G(B,T)$ data.

##  Appendix B Inelastic e-ph (single-particle) scattering time

The degree of disorder with respect to e-ph scattering is characterized by the
product $q_{T}l_{e}$. For our NbTiN films, in the experimental temperature
range, $q_{T}l_{e}\ll 1$ that corresponds to the strong disordered regime. The
SM model [24] describes e-ph scattering in disordered metals and provides the
inelastic (single-particle) scattering rate of electrons at the Fermi level
via interaction with phonons of different polarizations:

$\tau_{e-ph}^{-1}=\tau_{e-ph(l)}^{-1}+\tau_{e-ph(l)}^{-1},$ (7a) where
interaction with longitudinal phonons is accounted in
$\tau_{e-ph(l)}^{-1}=\frac{7\pi\zeta(3)}{2\hbar}\frac{\beta_{l}(k_{B}T)^{3}}{(p_{F}u_{l})^{2}}F_{l}(q_{T(l)}l_{e}),$
(7b) and with transverse phonons (two polarizations are taken into account) in
$\tau_{e-ph(t)}^{-1}=3\pi^{2}\frac{\beta_{t}(k_{B}T)^{2}}{(p_{F}u_{t})(p_{F}l_{e})}kF_{t}(q_{T(t)}l_{e}).$
(7c)

The indices ’$l$’ and ’$t$’ denote values associated with longitudinal and
transverse phonon modes. Here, $\beta_{l(t)}=(2E_{F}/3)^{2}(N(0)/(2\rho
u_{l(t)}^{2}))$ is the dimensionless coupling constant,
$E_{F}=(p_{F})^{2}/(2m_{e})$ is the Fermi energy, $p_{F}=\hbar k_{F}$ is the
Fermi momentum, and $m_{e}$ is the electron mass. In Eq. 7b, the integral
$F_{l}(z)=\frac{2}{7\zeta(3)}\int_{0}^{A_{l}}dx\,\Phi_{l}(xz)[N(x)+n(x)]x^{2}$,
where $N(x)$ and $n(x)$ are Bose and Fermi distribution functions, and
$\Phi_{l}(x)=\frac{2}{\pi}\left(\frac{x\arctan(x)}{x-\arctan(x)}-\frac{3}{x}k\right)$
is the Pippard function. In Eq. 7c, the integral
$F_{t}(z)=\frac{4}{\pi^{2}}\int_{0}^{A_{t}}dx\,\Phi_{t}(xz)[N(x)+n(x)]x$,
where $\Phi_{t}(x)=1+k[3x-3(x^{2}+1)\arctan(x)]/[2x^{3}]$. The upper limit of
the integrals $F(z)$ is $A_{l(t)}=(6\pi^{2})^{1/3}(a_{0}\,q_{T,l(t)})^{-1}$.
The parameter $1\geq k\geq 0$ reveals the property of electron scattering
centers ($k=1$ corresponds to scattering centers vibrating in the same way as
the host lattice, e.g. light impurities, $k=0$ to the static scattering
centers, e.g. heavy impurities and rigid boundaries).

The total inelastic e-ph scattering rate given by Eq. (7) is used to describe
the experimental $\tau_{e-ph}$ data shown in Fig. 5(a,b).

##  Appendix C WSi films

We describe the $\tau_{e-ph}(T)$ data obtained by the magnetoconductance
method for WSi films in [13] using the SM model (Appendix B). Following the
notation here, the authors of [13] found
$\alpha_{e-ph}\equiv\tau_{e-ph}(T_{C0})=66$ ps and $n=3$ between 5 and 20 K
which we added to Table 2. Fitting the SM model to these data, we used two
independent parameters $u_{t}$ and $\rho$; other parameters as $D$ and $N(0)$
were fixed at their values provided in [13] (also listed in Table. 1). We used
$a_{0}=0.46$ nm [49], $k_{F}=N(0)\pi^{2}\hbar^{2}/m_{e}$ computed wit the free
electron mass $m_{e}$, $l_{e}\sim 0.2$ nm, and $k\sim 1.0$. The best-fit
values of $u_{t}$ and $\rho$ are listed in Table 3.

## References

  * Esmaeil Zadeh  _et al._ [2021] I. Esmaeil Zadeh, J. Chang, J. W. Los, S. Gyger, A. W. Elshaari, S. Steinhauer, S. N. Dorenbos, and V. Zwiller, Superconducting nanowire single-photon detectors: A perspective on evolution, state-of-the-art, future developments, and applications, Applied Physics Letters 118, 190502 (2021).
  * Shurakov  _et al._ [2015] A. Shurakov, Y. Lobanov, and G. Goltsman, Superconducting hot-electron bolometer: from the discovery of hot-electron phenomena to practical applications, Superconductor Science and Technology 29, 023001 (2015).
  * Semenov  _et al._ [2002] A. D. Semenov, G. N. Gol’tsman, and R. Sobolewski, Hot-electron effect in superconductors and its applications for radiation sensors, Superconductor Science and Technology 15, R1 (2002).
  * Klapwijk and Semenov [2017] T. M. Klapwijk and A. Semenov, Engineering physics of superconducting hot-electron bolometer mixers, IEEE Transactions on Terahertz Science and Technology 7, 627 (2017).
  * Vodolazov [2019] D. Y. Vodolazov, Minimal timing jitter in superconducting nanowire single-photon detectors, Physical Review Applied 11, 014016 (2019).
  * Kardakova  _et al._ [2013] A. Kardakova, M. Finkel, D. Morozov, V. Kovalyuk, P. An, C. Dunscombe, M. Tarkhov, P. Mauskopf, T. Klapwijk, and G. Goltsman, The electron-phonon relaxation time in thin superconducting titanium nitride films, Applied Physics Letters 103, 252602 (2013).
  * Sidorova  _et al._ [2020] M. Sidorova, A. Semenov, H.-W. Hübers, K. Ilin, M. Siegel, I. Charaev, M. Moshkova, N. Kaurova, G. N. Goltsman, X. Zhang, _et al._ , Electron energy relaxation in disordered superconducting NbN films, Physical Review B 102, 054501 (2020).
  * Semenov [2021] A. D. Semenov, Superconducting nanostrip single-photon detectors some fundamental aspects in detection mechanism, technology and performance, Superconductor Science and Technology 34, 054002 (2021).
  * Le Jeannic  _et al._ [2016] H. Le Jeannic, V. B. Verma, A. Cavaillès, F. Marsili, M. D. Shaw, K. Huang, O. Morin, S. W. Nam, and J. Laurat, High-efficiency WSi superconducting nanowire single-photon detectors for quantum state engineering in the near infrared, Optics letters 41, 5341 (2016).
  * Verma  _et al._ [2015] V. B. Verma, B. Korzh, F. Bussieres, R. D. Horansky, S. D. Dyer, A. E. Lita, I. Vayshenker, F. Marsili, M. D. Shaw, H. Zbinden, _et al._ , High-efficiency superconducting nanowire single-photon detectors fabricated from MoSi thin-films, Optics express 23, 33792 (2015).
  * Steinhauer  _et al._ [2020] S. Steinhauer, L. Yang, S. Gyger, T. Lettner, C. Errando-Herranz, K. D. Jöns, M. A. Baghban, K. Gallo, J. Zichi, and V. Zwiller, NbTiN thin films for superconducting photon detectors on photonic and two-dimensional materials, Applied Physics Letters 116, 171101 (2020).
  * Zichi  _et al._ [2019] J. Zichi, J. Chang, S. Steinhauer, K. Von Fieandt, J. W. Los, G. Visser, N. Kalhor, T. Lettner, A. W. Elshaari, I. E. Zadeh, _et al._ , Optimizing the stoichiometry of ultrathin NbTiN films for high-performance superconducting nanowire single-photon detectors, Optics express 27, 26579 (2019).
  * Zhang  _et al._ [2016] X. Zhang, A. Engel, Q. Wang, A. Schilling, A. Semenov, M. Sidorova, H.-W. Hübers, I. Charaev, K. Ilin, and M. Siegel, Characteristics of superconducting tungsten silicide WxSi1-x for single photon detection, Physical Review B 94, 174509 (2016).
  * Aslamasov and Larkin [1968] L. Aslamasov and A. Larkin, The influence of fluctuation pairing of electrons on the conductivity of normal metal, Physics Letters A 26, 238 (1968).
  * Maki [1968] K. Maki, Critical fluctuation of the order parameter in a superconductor. I, Progress of Theoretical Physics 40, 193 (1968).
  * Thompson [1970] R. S. Thompson, Microwave, flux flow, and fluctuation resistance of dirty type-II superconductors, Physical Review B 1, 327 (1970).
  * Glatz  _et al._ [2011] A. Glatz, A. Varlamov, and V. Vinokur, Quantum fluctuations and dynamic clustering of fluctuating cooper pairs, EPL (Europhysics Letters) 94, 47005 (2011).
  * Minkov  _et al._ [2004] G. Minkov, A. Germanenko, and I. Gornyi, Magnetoresistance and dephasing in a two-dimensional electron gas at intermediate conductances, Physical Review B 70, 245423 (2004).
  * Bergmann [1984] G. Bergmann, Quantum corrections to the resistance in two-dimensional disordered superconductors above $T_{C}$: Al, Sn, and amorphous Bi0.9Tl0.1 films, Physical Review B 29, 6114 (1984).
  * Redi [1977] M. H. Redi, Two-dimensional fluctuation-induced conductivity above the critical temperature, Physical Review B 16, 2027 (1977).
  * Abrahams  _et al._ [1971] E. Abrahams, R. Prange, and M. Stephen, Effect of a magnetic field on fluctuations above TC, Physica 55, 230 (1971).
  * Tinkham [1980] M. Tinkham, Introduction to superconductivity, reprinted ed (1980).
  * Altshuler  _et al._ [1982] B. L. Altshuler, A. Aronov, and D. Khmelnitsky, Effects of electron-electron collisions with small energy transfers on quantum localisation, Journal of Physics C: Solid State Physics 15, 7367 (1982).
  * Sergeev and Mitin [2000] A. Sergeev and V. Mitin, Electron-phonon interaction in disordered conductors: Static and vibrating scattering potentials, Physical Review B 61, 6041 (2000).
  * Brenig  _et al._ [1985] W. Brenig, M.-c. Chang, E. Abrahams, and P. Wölfle, Inelastic scattering time above the superconductivity transition in two dimensions: Dependence on disorder and magnetic field, Physical Review B 31, 7001 (1985).
  * Rammer and Schmid [1986] J. Rammer and A. Schmid, Destruction of phase coherence by electron-phonon interactions in disordered conductors, Physical Review B 34, 1352 (1986).
  * Il’in  _et al._ [1998] K. Il’in, N. Ptitsina, A. Sergeev, G. Gol’tsman, E. Gershenzon, B. Karasik, E. Pechen, and S. Krasnosvobodtsev, Interrelation of resistivity and inelastic electron-phonon scattering rate in impure NbC films, Physical Review B 57, 15623 (1998).
  * Makise  _et al._ [2010] K. Makise, H. Terai, M. Takeda, Y. Uzawa, and Z. Wang, Characterization of NbTiN thin films deposited on various substrates, IEEE transactions on applied superconductivity 21, 139 (2010).
  * Arockiasamy  _et al._ [2016] M. Arockiasamy, M. Sundareswari, and M. Rajagopalan, Ductility behaviour of cubic titanium niobium nitride ternary alloy: a first-principles study, Indian Journal of Physics 90, 149 (2016).
  * Hazra  _et al._ [2018] D. Hazra, N. Tsavdaris, A. Mukhtarova, M. Jacquemin, F. Blanchet, R. Albert, S. Jebari, A. Grimm, A. Konar, E. Blanquet, _et al._ , Superconducting properties of NbTiN thin films deposited by high-temperature chemical vapor deposition, Physical Review B 97, 144518 (2018).
  * Pan and Guan [2017] Y. Pan and W. Guan, Probing the balance between ductility and strength: transition metal silicides, Physical Chemistry Chemical Physics 19, 19427 (2017).
  * Kaplan [1979] S. B. Kaplan, Acoustic matching of superconducting films to substrates, Journal of Low Temperature Physics 37, 343 (1979).
  * Frick  _et al._ [1975] W. Frick, D. Waldmann, and W. Eisenmenger, Phonon emission spectra of thin metallic films, Applied physics 8, 163 (1975).
  * Lang and Jiang [2005] X. Lang and Q. Jiang, Finite size effect on critical transition temperature of superconductive nanosolids, Solid state communications 134, 797 (2005).
  * Ma  _et al._ [2010] W. Ma, X. Zhang, and K. Takahashi, Electrical properties and reduced debye temperature of polycrystalline thin gold films, Journal of Physics D: Applied Physics 43, 465301 (2010).
  * Cheng  _et al._ [2015] Z. Cheng, L. Liu, S. Xu, M. Lu, and X. Wang, Temperature dependence of electrical and thermal conduction in single silver nanowire, Scientific reports 5, 1 (2015).
  * Engel  _et al._ [2006] A. Engel, A. Semenov, H.-W. Hübers, K. Il’in, and M. Siegel, Electric noise and local photon-induced nonequilibrium states in a current-carrying nanostructured superconductor, New Frontiers in Superconductivity Research 6, 153 (2006).
  * Ilin  _et al._ [2008] K. Ilin, R. Schneider, D. Gerthsen, A. Engel, H. Bartolf, A. Schilling, A. Semenov, H.-W. Hübers, B. Freitag, and M. Siegel, Ultra-thin NbN films on Si: crystalline and superconducting properties, in  _Journal of Physics: Conference Series_ , Vol. 97 (IOP Publishing, 2008) p. 012045.
  * Lin  _et al._ [2013] S.-Z. Lin, O. Ayala-Valenzuela, R. D. McDonald, L. N. Bulaevskii, T. G. Holesinger, F. Ronning, N. R. Weisse-Bernstein, T. L. Williamson, A. H. Mueller, M. A. Hoffbauer, _et al._ , Characterization of the thin-film NbN superconductor for single-photon detection by transport measurements, Physical Review B 87, 184507 (2013).
  * Zhang  _et al._ [2009] H. Zhang, D. J. Srolovitz, J. F. Douglas, and J. A. Warren, Grain boundaries exhibit the dynamics of glass-forming liquids, Proceedings of the National Academy of Sciences 106, 7735 (2009).
  * Baek  _et al._ [2011] B. Baek, A. E. Lita, V. Verma, and S. W. Nam, Superconducting a-WxSi1-x nanowire single-photon detector with saturated internal quantum efficiency from visible to 1850 nm, Applied Physics Letters 98, 251105 (2011).
  * Jin  _et al._ [2019] J. Jin, F. Fu, X. Jia, L. Kang, Z. Wang, X. Tu, L. Zhang, B. B. Jin, J. Chen, W. Xu, _et al._ , Preparation and characterization of ultrathin WSi films for superconducting nanowire single-photon detectors, IEEE Transactions on Applied Superconductivity 29, 1 (2019).
  * Golding  _et al._ [1972] B. Golding, B. Bagley, and F. Hsu, Soft transverse phonons in a metallic glass, Physical Review Letters 29, 68 (1972).
  * Bergmann [1976] G. Bergmann, Amorphous metals and their superconductivity, Physics Reports 27, 159 (1976).
  * Hikami  _et al._ [1980] S. Hikami, A. I. Larkin, and Y. Nagaoka, Spin-orbit interaction and magnetoresistance in the two dimensional random system, Progress of Theoretical Physics 63, 707 (1980).
  * Larkin [1980] A. Larkin, Reluctance of two-dimensional systems, JETP lett 31, 219 (1980).
  * dos Santos and Abrahams [1985] J. L. dos Santos and E. Abrahams, Superconducting fluctuation conductivity in a magnetic field in two dimensions, Physical Review B 31, 172 (1985).
  * Larkin and Varlamov [2003] A. Larkin and A. Varlamov, Fluctuation phenomena in superconductors, in  _The Physics of Superconductors_ (Springer, 2003) pp. 95–231.
  * Mattheiss [1992] L. Mattheiss, Calculated structural properties of CrSi2, MoSi2, and WSi2, Physical Review B 45, 3252 (1992).


# Linear State Estimation via 5G C-RAN Cellular Networks using Gaussian Belief
Propagation

Mirsad Cosovic, Dejan Vukobratovic, Vladimir Stankovic M. Cosovic is with
Schneider Electric DMS NS, Novi Sad, Serbia (e-mail: mirsad.cosovic@schneider-
electric-dms.com). D. Vukobratovic is with Department of Power, Electronic and
Communications Engineering, University of Novi Sad, Novi Sad, Serbia (e-mail:
dejanv@uns.ac.rs). V. Stankovic is with Department of Electronic and
Electrical Engineering, University of Strathclyde, Glasgow, UK
(e-mail:vladimir.stankovic@eee.strath.ac.uk).

###### Abstract

Machine-type communications and large-scale information processing
architectures are among key (r)evolutionary enhancements of emerging fifth-
generation (5G) mobile cellular networks. Massive data acquisition and
processing will make 5G network an ideal platform for large-scale system
monitoring and control with applications in future smart transportation,
connected industry, power grids, etc. In this work, we investigate a
capability of such a 5G network architecture to provide the state estimate of
an underlying linear system from the input obtained via large-scale deployment
of measurement devices. Assuming that the measurements are communicated via
densely deployed cloud radio access network (C-RAN), we formulate and solve
the problem of estimating the system state from the set of signals collected
at C-RAN base stations. Our solution, based on the Gaussian Belief-Propagation
(GBP) framework, allows for large-scale and distributed deployment within the
emerging 5G information processing architectures. The presented numerical
study demonstrates the accuracy, convergence behavior and scalability of the
proposed GBP-based solution to the large-scale state estimation problem.

##  I Introduction

With transition towards fifth generation (5G), mobile cellular networks are
evolving into ubiquitous systems for data acquisition and information
processing suitable for monitoring and control of large-scale systems. At the
forefront of this evolution is the transformation of radio access network
(RAN) to support massive-scale machine-type communications (MTC) [1] and
transformation of core network (CN) to support large-scale centralized or
distributed information processing through Cloud-RAN (C-RAN) and Fog-RAN
(F-RAN) architecture [2], [3]. MTC services in 5G will offer both massive-
scale data acquisition from various machine-type devices through massive MTC
(mMTC) service, but also, provide ultra reliable and low-latency communication
(URLLC) service for mission-critical applications [4]. Complemented with
ultra-dense RAN deployment and flexible and virtualized signal processing
architecture, novel 5G network services that are particularly suitable for
large-scale system monitoring and control of various _smart infrastructures_
are emerging [5].

In this work, we focus on a generic state estimation problem placed in the
context of a future 5G-inspired C-RAN-based cellular network. We consider an
underlying large-scale physical system characterized by the state vector
$\mathbf{s}$ that contains values of $N$ system state variables. The state
variables are observed through the set of $M$ measurements $\mathbf{x}$ of
physical quantities collected at the measurement devices spread across the
system. This paper considers linear system model in which measured quantities
are linear functions of the (sub)set of state variables. Further, we assume
measurements are wirelessly communicated across C-RAN-based cellular network.
In C-RAN, large number of spatially distributed remote radio heads (RRH)
constitutes an ultra-dense RAN infrastructure that receives signals from
densely populated MTC devices (e.g., the measurement devices under
consideration) [2]. The signal vector $\mathbf{y}$ collected at RRHs is
forwarded via backhaul links to a central C-RAN location where it is fed into
a collection of base-band units (BBU) for signal detection and recovery. In
the standard C-RAN signal detection problem, the goal is to recover the signal
$\mathbf{x}$ transmitted by the set of MTC devices from the signal
$\mathbf{y}$ received at RRHs and gathered centrally at BBUs [6][7]. However,
in this paper, focusing on widely applicable linear system state estimation,
we extend this goal and investigate the problem of recovering the system state
$\mathbf{s}$ directly from the signal $\mathbf{y}$ collected across the C-RAN.

The problem we observe represents a concatenation of the two well-studied
problems: the linear system state estimation problem (see, e.g., [8], for the
case of power system state estimation) and the problem of uplink signal
detection in C-RAN [6]. For the joint problem, it is straightforward to derive
(and implement at a central location) the standard minimum mean-square error
(MMSE) estimator, however, such a solution comes with prohibitive
$O(N^{3})$-complexity that hinders its application for large-scale systems. By
exploiting inherent sparsity within both of the component problems, an
approximate MMSE solution for each problem can be obtained using the tools
from probabilistic graphical models, as recently investigated for both (power
system) state estimation [9] and uplink signal detection in C-RANs [7]. In
particular, an instance of the Belief-Propagation (BP) algorithm, called
Gaussian BP (GBP) [10], can be applied to produce an exact MMSE estimate with
$O(N)$-complexity, thus scaling the MMSE solution to large-scale system
scenarios.

In this paper, we motivate, formulate and solve the linear system state
estimation problem considered jointly with the signal detection problem in
C-RAN-based cellular networks. We cast the problem of estimating the system
state $\mathbf{s}$ from the received vector $\mathbf{y}$ into an equivalent
maximum a-posteriori (MAP) problem, and place it into the framework of a
popular class of probabilistic graphical model called factor graphs. The state
estimate $\hat{\mathbf{s}}$ is then derived as a solution of the GBP algorithm
applied over a specific bi-layer structure of the factor graph. Throughout the
paper, we use state estimation in power systems with the measurements
collected via 5G-inspired C-RAN network as a running example. Our initial
numerical results demonstrate the viability of the proposed approach, both in
terms of accuracy and convergence.

The paper is organized as follows. In Section II, we present the joint state
estimation and C-RAN uplink communication system model. In Section III, this
model is mapped into a corresponding factor graph, and the state estimate is
obtained via GBP. Section IV provides numerical results of the proposed GBP
state estimator. The paper is concluded in Section V.

##  II System Model

Figure 1: System model example: State estimation in Smart Grid via C-RAN-based
mobile cellular network.

We consider a generic state estimation problem where a set of state variables
$\mathbf{s}$ is estimated from a set of observed noisy linear measurements
$\mathbf{x}$. However, unlike the traditional setup where the measurements in
$\mathbf{x}$ are assumed available at a central node, here we assume they are
transmitted via radio access network (RAN) of a mobile cellular system based
on a cloud-RAN (C-RAN) architecture. The received signal $\mathbf{y}$ is
collected at a large-number of densely deployed remote radio heads (RRHs) and
jointly processed at the C-RAN base-band units (BBUs). The problem we consider
is that of estimating the system state $\mathbf{s}$ from the received signal
$\mathbf{y}$.

Linear system measurements model: We consider a system described via the set
of $N$ state variables
$\mathbf{s}=(s_{1},s_{2},\ldots,s_{N})^{T}\in\mathbb{C}^{N\times 1}$. The
system is observed via a set of $M$ measurements
$\mathbf{x}=(x_{1},x_{2},\ldots,x_{M})^{T}\in\mathbb{C}^{M\times 1}$. Each
measurement is a linear function of the state variables additionally corrupted
by the additive noise, i.e., $x_{i}=\mathbf{a}_{i}\cdot\mathbf{s}+n_{i}$,
$1\leq i\leq M$, where
$\mathbf{a}_{i}=(a_{i,1},a_{i,2},\ldots,a_{i,N})\in\mathbb{C}^{1\times N}$ is
a vector of coefficients, while $n_{i}\in\mathbb{C}$ is a complex random
variable. Overall, the system is represented via noisy linear observation
model $\mathbf{x}=\mathbf{A}\cdot\mathbf{s}+\mathbf{n}$, where the matrix
$\mathbf{A}\in\mathbb{C}^{M\times N}$ contains vectors $\mathbf{a}_{i},1\leq
i\leq M,$ as rows, while
$\mathbf{n}=(n_{1},n_{2},\ldots,n_{M})^{T}\in\mathbb{C}^{M\times 1}$ is a
vector of additive noise samples. We assume noise samples $n_{i}$ are
independent identically distributed (i.i.d.) Gaussian random variables with
variance $\sigma_{n_{i}}^{2}$. For simplicity, we assume the measurement noise
variances are equal, i.e., $\sigma_{n_{i}}^{2}=\sigma_{n}^{2},1\leq i\leq M$.
In other words, $\mathbf{n}$ represents a complex Gaussian random vector with
zero means $\mathbf{\mu}_{\mathbf{n}}=\mathbf{0}$ and the covariance matrix
$\mathbf{\Sigma}_{\mathbf{n}}=\sigma_{n}^{2}\mathbf{I}\in\mathbb{C}^{M\times
M}$ ($\mathbf{I}$ is an identity matrix).

C-RAN uplink communication model: In the standard state estimation models,
measurements are either assumed available, or they are communicated to the
central node, where the state estimation problem is solved. Communication
models typically involve point-to-point communication links between the
measurement devices and the central node, affected by communication
impairments such as delays, packet losses, limited bit rates, etc. In the
cellular networks context, this assumes reservation of uplink resources and
subsequent non-orthogonal transmission, which typically incurs significant
communication delays. Inspired by the recent evolution of massive MTC and
ultra-dense C-RAN architectures in upcoming 5G mobile cellular networks, in
this work, we consider different grant-less and non-orthogonal communication
model, as we detail next. Note that the following C-RAN cellular network model
could provide ultra-low latency for the state estimation application under
consideration. In other words, such an architecture could produce the system
state estimate at the central network node with very low delay after the
measurements are acquired, which is crucial for emerging mission-critical 5G
MTC use cases [11].

In the mobile cellular system under consideration, measurements are collected
by the measurement devices that we refer to as MTC user equipment (MTC-UE). We
consider uplink transmission of $M$ single-antenna MTC UEs towards $L$ single-
antenna RRHs. We assume both MTC UEs and RRHs are randomly and uniformly
distributed across a given geographic area (note that this placement model is
somewhat refined in the numerical results section). The signal
$\mathbf{x}\in\mathbb{C}^{M\times 1}$, representing the set of collected noisy
measurements, is transmitted111Note that, at this point, one can insert
specific linear modulation scheme $\mathbf{x}_{m}=f_{m}(\mathbf{x})$. For
simplicity, we assume $f_{m}(\mathbf{x})=\mathbf{x}$. by $M$ MTC UEs, while
the signal $\mathbf{y}=(y_{1},y_{2},\ldots,y_{L})\in\mathbb{C}^{L\times 1}$ is
received at $L$ RRHs, where $\mathbf{y}=\mathbf{H}\cdot\mathbf{x}+\mathbf{m}$.
The matrix $\mathbf{H}\in\mathbb{C}^{L\times M}$ represents the channel
matrix, where $h_{i,j}$ represents a complex channel coefficient between the
$j$-th MTC UE and the $i$-th RRH, while
$\mathbf{m}=(m_{1},m_{2},\ldots,m_{L})\in\mathbb{C}^{L}$ is a vector of
additive noise samples. As for the measurement process, for the communication
process we also assume noise samples $m_{i}$ are i.i.d. zero-mean Gaussian
random variables with variance $\sigma_{m}^{2}$, i.e., the mean value and the
covariance matrix of $\mathbf{m}$ is given as
$\mathbf{\mu}_{\mathbf{m}}=\mathbf{0}$ and
$\mathbf{\Sigma}_{\mathbf{m}}=\sigma_{m}^{2}\mathbf{I}\in\mathbb{C}^{L\times
L}$, respectively.

Linear system state estimation problem: From the received signal $\mathbf{y}$
collected across RRHs, we are interested in finding an estimate
$\mathbf{\hat{s}}$ of the state vector $\mathbf{s}$. In this paper, we focus
on the centralized C-RAN architecture, where all the BBUs are collocated at
the central C-RAN node. Thus, due to availability of $\mathbf{y}$ at the
central location, we consider centralized algorithms for the state estimation
problem. However, the solution we propose in this paper is based on GBP
framework, thus it is easily adaptable to a distributed scenario that we refer
to as fog-RAN (F-RAN), where BBUs are distributed across different geographic
locations closer to the MTC UEs. We refer the interested reader to our recent
overview of distributed algorithms for solving the state estimation problem in
the context of upcoming 5G cellular networks [5].

A common centralized approach to solve the state estimation problem is to
provide the minimum mean-square error (MMSE) estimate. One can easily obtain
the MMSE estimate $\mathbf{\hat{s}}$ for the underlying linear model in the
form:

$\mathbf{\hat{s}}=\big{[}\mathbf{\Sigma}_{s}^{-1}+(\mathbf{H}\mathbf{A})^{H}\mathbf{\Sigma}^{-1}\cdot(\mathbf{H}\mathbf{A})\big{]}^{-1}(\mathbf{H}\mathbf{A})^{H}\mathbf{\Sigma}^{-1}\mathbf{y},$
(1)

where $(\cdot)^{H}$ is the conjugate-transpose matrix operation,
$\mathbf{\Sigma}=\mathbf{H}\mathbf{\Sigma}_{\mathbf{n}}\mathbf{H}^{H}+\mathbf{\Sigma}_{\mathbf{m}}$,
and where we assume prior distribution of $\mathbf{s}$ is i.i.d. Gaussian with
mean $\mathbf{\mu}_{\mathbf{s}}=\mathbf{0}$ and
$\mathbf{\Sigma}_{\mathbf{s}}=\sigma_{s}^{2}\mathbf{I}\in\mathbb{C}^{N\times
N}$. However, solving (1) scales as $O(N^{3})$ which makes linear MMSE state
estimation inapplicable in large-scale systems which are of interest in this
paper.

In the following, we cast the MMSE estimation into an equivalent MAP state
estimator as follows:

$\mathbf{\hat{s}}=\arg\max_{\mathbf{s}\in\mathbb{C}^{N}}P(\mathbf{s}|\mathbf{y}),$
(2)

where $P(\mathbf{s}|\mathbf{y})$ is the posterior probability of the state
$\mathbf{s}$ after the signal $\mathbf{y}$ is observed at the RRHs. As we will
demonstrate in the sequel, if certain sparsity arguments are applicable in the
system model under consideration, the solution of the MAP problem can be
efficiently calculated using the framework of factor graphs and belief-
propagation (BP) algorithms.

System model example (Smart Grid): Before continuing, it is useful to consider
an example of the above state estimation setup. We consider the state
estimation problem in an electric power system, where the goal is to estimate
the state of the power system $\mathbf{s}$, containing complex voltages of $N$
system buses, via the set of measurements $\mathbf{x}$ obtained using
measurement devices. Measurement devices are geographically distributed across
the power system and we assume they are equipped with wireless cellular
interfaces, i.e., they represent MTC UEs connected to the C-RAN based cellular
network. The signal $\mathbf{y}$ received from the set of RRHs densely
deployed across the cellular network coverage area is processed centrally
within the C-RAN system architecture. Fig. 1 illustrates the smart grid
example that we will further refine in Section IV and use as a running example
throughout this paper.

##  III State Estimation via Gaussian Belief Propagation

In this section, we provide a solution to the combined state estimation and
uplink signal detection problem defined in (2), by applying factor graphs and
GBP framework. In fact, for both constituent scenarios: the conventional state
estimation (in case of power systems) and the uplink signal detection in
C-RAN, the GBP has already been proposed and analyzed (see details in [9] and
[7]). Thus in this work, we propose using GBP in a joint and combined setup of
extracting the system state directly from the observed C-RAN signals.

We note that the various properties of the proposed GBP approach (e.g.,
complexity, convergence, etc.) will strongly depend on the structure of the
underlying factor graph. For example, in terms of complexity (we will come
back to convergence in the next section), for GBP to scale well to large-scale
systems, it is fundamental that both matrices $\mathbf{A}$ and $\mathbf{H}$
defining the two linear problems are sparse, i.e., that for both $\mathbf{A}$
and $\mathbf{H}$, the number of non-zero entries scales as $O(N)$. In many
real-world scenarios, the sparsity typically arises from geographic
constraints and reflect _locality_ that is typically present in both the
measurement and the communication part of the system model. More detailed
account on the sparsity of matrices $\mathbf{A}$ and $\mathbf{H}$ clearly
depends on the specific scenario under consideration, and we relegate these
details to Section IV where we will explicitly deal with the smart grid
example introduced earlier.

Factor Graph System Representation: Assuming that the system state
$\mathbf{s}$ has a Gaussian prior, and given that the measurement and
communication noise is assumed Gaussian, the MAP problem can be rewritten as
follows:

$\displaystyle\mathbf{\hat{s}}$ $\displaystyle=$
$\displaystyle\arg\max_{\mathbf{s}\in\mathbb{C}^{N}}P(\mathbf{s}|\mathbf{y})\propto\arg\max_{\mathbf{s}\in\mathbb{C}^{N}}P(\mathbf{s},\mathbf{y})=$
(3) $\displaystyle=$
$\displaystyle\arg\max_{\mathbf{s}\in\mathbb{C}^{N},\mathbf{x}\in\mathbb{C}^{M}}P(\mathbf{s},\mathbf{x},\mathbf{y})$
(4)

Note that the distribution $P(\mathbf{s},\mathbf{x},\mathbf{y})$ is jointly
Gaussian, and in addition, due to the problem structure where $\mathbf{s}$ and
$\mathbf{y}$ are conditionally independent given $\mathbf{x}$, we obtain:

$\displaystyle\mathbf{\hat{s}}$ $\displaystyle=$
$\displaystyle\arg\max_{\mathbf{s}\in\mathbb{C}^{N},\mathbf{x}\in\mathbb{C}^{M}}P(\mathbf{y}|\mathbf{s},\mathbf{x})P(\mathbf{s},\mathbf{x})=$
(5) $\displaystyle=$
$\displaystyle\arg\max_{\mathbf{s}\in\mathbb{C}^{N},\mathbf{x}\in\mathbb{C}^{M}}P(\mathbf{y}|\mathbf{x})P(\mathbf{x}|\mathbf{s})P(\mathbf{s}).$
(6)

As noted before, in many real-world systems of interest, a measurement $x_{j}$
is a linear function of a small subset of _local_ state variables
$\mathbf{s}_{\mathcal{N}(x_{j})}$, where $\mathcal{N}(x_{j})$ is the index set
of the state variables that affect $x_{j}$, and
$\mathbf{s}_{\mathcal{N}(x_{j})}=\\{s_{i}|i\in\mathcal{N}(x_{j})\\}$. In other
words, the row-vector $\mathbf{a}_{j}$ has non-zero components only on a small
number of positions indexed by the set $\mathcal{N}(x_{j})$, thus making the
matrix $\mathbf{A}$ sparse222More precisely, the number of state variables
that affect certain measurement is limited by a constant, independently of the
size $N$ of the system.. Using this fact and the fact that the measurements
$x_{j}$ are mutually independent, we obtain:

$P(\mathbf{x}|\mathbf{s})=\prod_{j=1}^{M}P(x_{j}|\mathbf{s}_{\mathcal{N}(x_{j})}).$
(7)

On the other hand, in the C-RAN communication part, although in theory the
received signal $y_{i}$ depends on all the transmitted symbols in
$\mathbf{x}$, the channel coefficients between a RRH and a geographically
distant MTC UE can be considered negligible, thus leading to matrix
$\mathbf{H}$ sparsification [7]. Upon distance-based sparsification proposed
in [7], the received symbol $y_{i}$ depends only on a small number of symbols
$\mathbf{x}_{\mathcal{N}(y_{i})}$, where $\mathcal{N}(y_{i})$ is the index set
of symbols transmitted by the set of MTC UEs in geographic proximity of the
$i$-th RRH. Taking the channel sparsification into account, we obtain:

$P(\mathbf{y}|\mathbf{x})=\prod_{i=1}^{L}P(y_{i}|\mathbf{x}_{\mathcal{N}(y_{i})}).$
(8)

Finally, assuming that the state vector is apriori given as a set of i.i.d
Gaussian random variables, we obtain the final factorized form of the initial
MAP problem:

$\mathbf{\hat{s}}=\arg\max_{\mathbf{s}\in\mathbb{C}^{N},\mathbf{x}\in\mathbb{C}^{M}}\prod_{i=1}^{L}P(y_{i}|\mathbf{x}_{\mathcal{N}(y_{i})})\cdot\\\
\prod_{j=1}^{M}P(x_{j}|\mathbf{s}_{\mathcal{N}(x_{j})})\cdot\prod_{k=1}^{N}P(s_{k}).$
(9)

The factor graph representation of the MAP problem follows the factorization
presented in (9) and is illustrated in Fig. 2. Factor graph
$\mathcal{G}=\mathcal{G}(\mathcal{V}\cup\mathcal{F},\mathcal{E})$ is a
bipartite graph consisting of the set of variable nodes $\mathcal{V}$, the set
of factor nodes $\mathcal{F}$, and the set of edges $\mathcal{E}$. In our
setup, the set $\mathcal{V}$ can be further divided as
$\mathcal{V}=\mathcal{S}\cup\mathcal{X}\cup\mathcal{Y},$ where
$\mathcal{S}=\\{s_{1},s_{2},\ldots,s_{N}\\}$ is the set of state nodes,
$\mathcal{X}=\\{x_{1},x_{2},\ldots,s_{M}\\}$ is the set of measurement nodes,
while $\mathcal{Y}=\\{y_{1},y_{2},\ldots,y_{L}\\}$ is the set of received
symbol nodes. The set of factor nodes can be divided as
$\mathcal{F}=\mathcal{F}_{\mathbf{H}}\cup\mathcal{F}_{\mathbf{A}}\cup\mathcal{F}_{\mathbf{y}}\cup\mathcal{F}_{\mathbf{x}}\cup\mathcal{F}_{\mathbf{s}},$
where
$\mathcal{F}_{\mathbf{H}}=\\{f_{\mathbf{h}_{1}},f_{\mathbf{h}_{2}},\ldots,f_{\mathbf{h}_{L}}\\}$
and
$\mathcal{F}_{\mathbf{A}}=\\{f_{\mathbf{a}_{1}},f_{\mathbf{a}_{2}},\ldots,f_{\mathbf{a}_{M}}\\}$
represent factor nodes that capture linear relationships between variable
nodes described by the rows of matrices $\mathbf{H}$ and $\mathbf{A}$,
respectively. In addition,
$\mathcal{F}_{\mathbf{y}}=\\{f_{y_{1}},f_{y_{2}},\ldots,f_{y_{L}}\\}$ and
$\mathcal{F}_{\mathbf{s}}=\\{f_{s_{1}},f_{s_{2}},\ldots,f_{s_{N}}\\}$
represent the factor nodes that provide inputs due to observations of
$\mathbf{y}$ and the prior knowledge about $\mathbf{x}$, respectively, while
$\mathcal{F}_{\mathbf{x}}=\\{f_{x_{1}},f_{x_{2}},\ldots,f_{x_{M}}\\}$ serve as
virtual inputs needed for initialization of measurement nodes. Similarly, the
set of edges $\mathcal{E}$ can be divided as
$\mathcal{E}=\mathcal{E}_{\mathbf{H}}\cup\mathcal{E}_{\mathbf{A}}\cup\mathcal{E}_{\mathbf{y}}\cup\mathcal{E}_{\mathbf{x}}\cup\mathcal{E}_{\mathbf{s}}$,
where
$\mathcal{E}_{\mathbf{H}}\cup\mathcal{E}_{\mathbf{y}}\cup\mathcal{E}_{\mathbf{x}}$
and
$\mathcal{E}_{\mathbf{A}}\cup\mathcal{E}_{\mathbf{s}}\cup\mathcal{E}_{\mathbf{x}}$
can be considered as the set of edges of two bipartite subgraphs
$\mathcal{G}_{\mathbf{H}}=(\mathcal{Y}\cup\mathcal{X}\cup\mathcal{F}_{\mathbf{H}}\cup\mathcal{F}_{\mathbf{y}}\cup\mathcal{F}_{\mathbf{x}},\mathcal{E}_{\mathbf{H}}\cup\mathcal{E}_{\mathbf{y}}\cup\mathcal{E}_{\mathbf{x}})$
and
$\mathcal{G}_{\mathbf{A}}=(\mathcal{X}\cup\mathcal{S}\cup\mathcal{F}_{\mathbf{A}}\cup\mathcal{F}_{\mathbf{s}}\cup\mathcal{F}_{\mathbf{x}},\mathcal{E}_{\mathbf{A}}\cup\mathcal{E}_{\mathbf{s}}\cup\mathcal{E}_{\mathbf{x}})$,
obtained as the subgraphs of $\mathcal{G}$ induced from the set of factor
nodes
$\mathcal{F}_{\mathbf{H}}\cup\mathcal{F}_{\mathbf{y}}\cup\mathcal{F}_{\mathbf{x}}$
and
$\mathcal{F}_{\mathbf{A}}\cup\mathcal{F}_{\mathbf{s}}\cup\mathcal{F}_{\mathbf{x}}$,
respectively333As noted in footnote 1, if the signal $\mathbf{x}$ is modulated
prior to transmission, one can easily add an additional “layer” to the factor
graph in Fig. 2 containing a set of $M$ modulated signal variable nodes
$\mathcal{X}_{m}$ connected via modulation factor nodes $\mathcal{F}_{m}$ with
the corresponding measurement nodes $\mathcal{X}$..

Figure 2: Factor graph representation of the system model.

As noted earlier, the state estimation problem using GBP over factor graph
$\mathcal{G}_{\mathbf{A}}$, and the uplink C-RAN signal detection problem
using GBP over factor graph $\mathcal{G}_{\mathbf{H}}$, have been recently
investigated in detail in [9] and [7], respectively.

Gaussian Belief-Propagation and GBP Messages: To estimate the state variables
$\mathbf{s}$, we apply message-passing GBP algorithm [10]. GBP operates on the
factor graph $\mathcal{G}$ by exchanging messages between factor nodes and
variable nodes in both directions. As a general rule, at any variable or
factor node, an outgoing message on any edge is obtained as a function of
incoming messages from all other edges, using the message calculation rules
presented below. In general, the underlying factor graph describing joint
state estimation and signal detection problem (Fig. 2) will contain cycles,
thus the resulting GBP will be iterative, which means that all nodes will
iteratively repeat message updates on all of the outgoing edges according to a
given message-passing schedule. We provide details on message-passing
schedule, correctness and convergence of GBP on loopy graphs later in this
section.

Let us consider a variable node $v_{i}\in\mathcal{V}$ incident to a factor
node $f_{j}\in\mathcal{F}$. Let $\mathcal{N}(v_{i})$ denote the index set of
factor nodes incident to $v_{i}$, and $\mathcal{N}(f_{j})$ denote the index
set of variable nodes incident to $f_{j}$. We denote messages from $v_{i}$ to
$f_{j}$ and from $f_{j}$ to $v_{i}$ as $\mu_{v_{i}\to
f_{j}}(v_{i})=(m_{v_{i}\to f_{j}},\sigma^{2}_{v_{i}\to f_{j}})$ and
$\mu_{f_{j}\to v_{i}}(v_{i})=(m_{f_{j}\to v_{i}},\sigma^{2}_{f_{j}\to
v_{i}})$, respectively. Note that, in the GBP scenario, all messages exchanged
across the factor graph represent Gaussian distributions defined by the
corresponding mean-variance pairs $(m,\sigma^{2})$. Thus to describe
processing rules in a variable and a factor node, it is sufficient to provide
equations that map input $(m,\sigma^{2})$-pairs into the
output$(m,\sigma^{2})$-pair, as detailed below.

_Message from a variable node to a factor node:_ the equations below are used
to calculate $\mu_{v_{i}\to f_{j}}(v_{i})=(m_{v_{i}\to
f_{j}},\sigma^{2}_{v_{i}\to f_{j}})$:

$\displaystyle m_{v_{i}\to f_{j}}$
$\displaystyle=\Bigg{(}\sum_{k\in\mathcal{N}(v_{i})\setminus
j}\cfrac{m_{f_{k}\to v_{i}}}{\sigma_{f_{k}\to
v_{i}}^{2}}\Bigg{)}\sigma_{v_{i}\to f_{j}}^{2}$ (10a)
$\displaystyle\cfrac{1}{\sigma_{v_{i}\to f_{j}}^{2}}$
$\displaystyle=\sum_{k\in\mathcal{N}(v_{i})\setminus
j}\cfrac{1}{\sigma_{f_{k}\to v_{i}}^{2}}.$ (10b)

_Message from a factor node to a variable node:_ In the setup under
consideration, factor nodes represent linear relations between variable nodes.
Thus, e.g., for a factor node $f_{j}$, we can write the corresponding linear
relationship as:

$\begin{gathered}f_{j}(\mathbf{v}_{\mathcal{N}(f_{j})})=C_{i}v_{i}+\sum_{k\in\mathcal{N}(f_{j})\setminus
i}C_{k}v_{k}.\end{gathered}$ (11)

With this general notation, the equations below provide $\mu_{f_{j}\to
v_{i}}(v_{i})=(m_{f_{j}\to v_{i}},\sigma^{2}_{f_{j}\to v_{i}})$:

$\displaystyle m_{f_{j}\to v_{i}}$
$\displaystyle=\cfrac{1}{C_{i}}\Bigg{(}\sum_{k\in\mathcal{N}(f_{j})\setminus
i}C_{k}m_{v_{k}\to f_{j}}\Bigg{)}$ (12a) $\displaystyle\sigma_{f_{j}\to
v_{i}}^{2}$
$\displaystyle=\cfrac{1}{C_{i}^{2}}\Bigg{(}\sum_{k\in\mathcal{N}(f_{j})\setminus
i}C_{k}^{2}\sigma_{v_{k}\to f_{j}}^{2}\Bigg{)}.$ (12b)

_Calculation of marginals:_ Applying the above rules in variable and factor
nodes of the factor graph results in the sequence of updates of messages
exchanged across the edges of the graph. To complete description of loopy GBP,
we need to define message initialization at the start, and message scheduling
during the course of each iteration, which is done next. After sufficient
number of GBP iterations, the final marginal distributions of the random
variables corresponding to variable nodes is obtained as:

$\displaystyle\hat{m}_{v_{i}}$
$\displaystyle=\Bigg{(}\sum_{k\in\mathcal{N}(v_{i})}\cfrac{m_{f_{k}\to
v_{i}}}{\sigma_{f_{k}\to v_{i}}^{2}}\Bigg{)}\sigma_{v_{i}}^{2}$ (13a)
$\displaystyle\cfrac{1}{\hat{\sigma}_{v_{i}}^{2}}$
$\displaystyle=\sum_{k\in\mathcal{N}(v_{i})}\cfrac{1}{\sigma_{f_{k}\to
v_{i}}^{2}}.$ (13b)

GBP Message-Passing Schedule, Correctness and Convergence: We adopt standard
synchronous GBP schedule in which variable node processing is done in the
first half-iteration, followed by the factor node processing in the second
half-iteration. The iterations are initialized by input messages from
$\mathcal{F}_{\mathbf{y}}$ generated from the received signal $\mathbf{y}$,
and initial messages from $\mathcal{F}_{\mathbf{x}}$ and
$\mathcal{F}_{\mathbf{s}}$ that follow certain prior knowledge (as detailed in
the next section).

GBP performance on linear models defined by loopy factor graphs is fairly well
understood. For example, if the GBP converges, it is known that the GBP
solution will match the solution of the MMSE estimator. The convergence
criteria can also be derived in a straightforward manner, by deriving
recursive fixed point linear transformations that govern mean value and
variance updates through the iterations and investigating spectral radius of
such transformations. Due to space restrictions, we leave the details of the
convergence analysis in our scenario for the future work.

##  IV Numerical Case Study: Smart Grid State Estimation in 5G C-RAN

In this section, we specialize our state estimation setup for a case study in
which we perform power system state estimation by collecting measurements via
5G-inspired C-RAN.

Power system state estimation - DC model: For the sake of simplicity, in the
following, we consider the linear DC model of a power system. The DC model is
an approximate model obtained as a linear approximation of the non-linear AC
model that precisely follows the electrical physical laws of the power system.
In the DC model, the power system containing $N$ buses is described by $N$
state variables $\mathbf{s}=(s_{1},s_{2},\ldots,s_{N})^{T}$, where each state
variable $s_{i}=\theta_{i}$ represents the voltage angle $\theta_{i}$ (in the
DC model, the magnitudes of all voltage phasors are assumed to have unit
values). In the DC model, the measurements include only active power flow
$P_{rk}$ at the branch $(r,k)$ between the bus $r$ and the bus $k$, active
power injection $P_{r}$ into the bus $r$, and the voltage angle $\theta_{r}$.
Collecting $M$ of such arbitrary measurements across the power system, we
obtain the measurement vector $\mathbf{x}=(x_{1},x_{2},\ldots,x_{M})^{T}$,
where each measurement $x_{i}\in\\{P_{rk},P_{r},\theta_{r}\\}$ is a linear
function444More precisely, we have that
$P_{rk}=-b_{rk}(\theta_{r}-\theta_{k})$ and
$P_{r}=-\sum_{k\in\mathcal{N}_{r}}b_{rk}(\theta_{r}-\theta_{k})$, where
$\mathcal{N}_{r}$ is the set of adjacent buses of the bus $r$ and $b_{rk}$ is
susceptance of the branch $(r,k)$. of the (sub)set of state variables
$\mathbf{s}$, additionally corrupted by additive Gaussian noise of fixed
(normalized) noise standard deviation of $\sigma_{n}$ per unit (p.u.). The
noisy measurements $\mathbf{x}$ are then transmitted via C-RAN network as
described below.

We illustrate the methodology using the IEEE test bus case with 30 buses shown
in Fig. 3 ($N=29$, since one of the bus voltage angles is set to the reference
value zero) that we use in the simulations. The example set of $M$
measurements is selected in such a way that the system is observable with the
redundancy $M/N$. For each simulation scenario, we generate 1000 random
(observable) measurement configurations.

Figure 3: The IEEE 30 bus test case divided into disjoint sub-rectangles.

C-RAN cellular network model: The set of $M$ MTC-UEs simultaneously transmit
their measurements to the set of $L$ RRHs during a given allocated time-
frequency slot shared by all MTC-UEs. We assume MTC-UEs and RRHs are placed
uniformly at random following independent Poisson Point Process (PPP) in a
unit-square area, however, with slight refinement of the PPP placement
strategy. Namely, to account for neighboring relations within logical topology
of IEEE 30 bus test case, we first divide a unit-square into $w\times q$
disjoint sub-rectangles as shown in Fig. 3, and then we assign $M$ MTC-UEs to
one of $w\cdot q$ sub-rectangles. We also balance the number of RRHs per sub-
rectangle, thus allocating $\sim{L}/(w\cdot q)$ RRHs per sub-rectangle.
Finally, all RRHs and MTC-UEs allocated to a given sub-rectangle are placed
using the PPP within a given sub-rectangle.

After the placement, we assume $M$ MTC-UEs transmit their signals
$\mathbf{x}$, where each measured signal is normalized to its expected
normalization value555We assume normalization constants are known in advance
at MTC-UEs and C-RAN nodes, either as a prior knowledge or by long-term
averaging.. For the channel coefficients between the MTC UEs and RRHs, we
assume the following model:

$h_{i,j}=\gamma_{i,j}d_{i,j}^{-\alpha},$ (14)

where $\gamma_{i,j}$ is the i.i.d. Rayleigh fading coefficient with zero mean
and unit variance, $d_{i,j}$ is distance between $i$-th MTC-UE and $j$-th RRH
and $\alpha$ is the path loss exponent. We use channel sparsification approach
proposed in [7], with threshold distance set to $d_{0}={\sqrt{w^{2}+q^{2}}}$
(i.e., equal to the diagonal length of each sub-rectangle). The received
signal $\mathbf{y}=(y_{1},y_{2},\ldots,y_{L})$ collected at $L$ RRHs is
additionally corrupted by additive Gaussian noise, whose standard deviation is
selected so as to provide fixed and pre-defined signal-to-noise ratio (SNR)
value. Finally, noisy received signal $\mathbf{y}$ is forwarded via high-
throughput backhaul links to C-RAN BBUs.

GBP-based State Estimation: Using the approach presented in Section III, we
apply GBP across the factor graph illustrated in Fig. 2 to recover the state
estimate $\mathbf{x}$ from the received signal $\mathbf{y}$. More precisely,
for each random measurement configuration, we generate the part of the factor
graph $\mathcal{G}_{\mathbf{A}}$ and, similarly, from known MTC-UE and RRH
random positions, we derive666We note that, in case the small-scale fading is
included in the model, one can assume that the channel state information is
available at the C-RAN. the part of the factor graph
$\mathcal{G}_{\mathbf{H}}$. Upon reception of $\mathbf{y}$, the GBP runs until
it converges. We adopt a synchronous scheduling of GBP messages where messages
are synchronously flooded from the factor nodes to variable nodes and back
within a single GBP iteration. For a linear model, it is well known that if
the GBP converges, it will converge to the minimum mean-square error (MMSE)
estimate of the state $\mathbf{x}$.

Simulation Results: In the first set of experiments, we fix the relative RRH
density ${L}/{M}=1$, $\mathrm{SNR}=10$, and redundancy ${M}/{N}=3$. We
investigate the accuracy of the GBP solution of state estimate as a function
of different values of measurements noise $\sigma_{n}=$ $\\{10^{-1},$
$10^{-2},$ $10^{-3},$ $10^{-4}\\}$.

$10^{-1}$$10^{-2}$$10^{-3}$$10^{-4}$$0.00$$2.00$$4.00$$6.00$$8.00$Measurement
noise $\sigma_{n}(\mathrm{p.u.})$$\mathrm{RMSE}$ Figure 4: The root mean
square error of estimate vector of power system state variables obtained by
C-RAN and without C-RAN model.

Fig. 4 shows the root mean square error
$\mathrm{RMSE}={({1}/{N})||\mathbf{\hat{s}}_{\mathrm{c}}-\mathbf{\hat{s}}_{\bar{\mathrm{c}}}||_{2}}$,
where $\mathbf{\hat{s}}_{\mathrm{c}}$ and
$\mathbf{\hat{s}}_{\bar{\mathrm{c}}}$ are estimate vectors of power system
state variables obtained with and without C-RAN model discussed in this paper,
respectively, for different values of measurement noise $\sigma_{n}$. For the
case without C-RAN model, we assume measurements $\mathbf{x}$ are available at
BBUs as they are, i.e., without additional noise or errors. In practice, this
could be obtained via standard grant-based uplink procedures where each MTC UE
is allocated separate orthogonal resources. However, such a strategy incurs
significant delay as the underlying system scales, due to message exchange
delay, resource allocation delay, as well as ARQ-based error-correction
strategies. Note that the C-RAN model described in this paper admits very low
latency as all MTC UEs transmit their signals immediately and concurrently.
According to the box plot in Fig. 4, the C-RAN approach is able to reach
nearly identical solution as the approach without C-RAN (e.g.,
$\mathrm{RMSE}\to 0$), if the value of measurements noise is sufficiently low.
Note that the typical value (standard deviation) of the measurement noise for
devices located across a power system are in the range between
$10^{-2}\,\mathrm{p.u.}$ and $10^{-3}\,\mathrm{p.u.}$, for legacy measurement
devices, and between $10^{-4}\,\mathrm{p.u.}$ and $10^{-5}\,\mathrm{p.u.}$,
for phasor measurement units. Consequently, the presented approach is suitable
for the state estimation in power systems.

In the next simulation experiment, we investigated the system observability as
a function of the number of RRHs $L$ deployed in the system, for different
values of redundancy $M/N$. We start with $L/M=0.2$ and increase the RRH
density in order to evaluate its effect on the system observability.

$0.2$$0.4$$0.6$$0.8$$1$$1.2$$0$$50$$100$Relative RRH density
${L}/{M}$Unobservable Topologies (%)$M/N$1.52.02.53.03.4 Figure 5: The
fraction of unobservable system topologies for different values of measurement
redundancies $M/N$ versus relative RRH density ${L}/{M}$.

Fig. 5 shows the fraction of instances GBP was not able to converge due to
insufficient rank of the underlying system as a function of the number of base
stations $L$. Note that the fundamental condition for the system to have full
rank is that $L\geq N$. By slightly expanding this condition, we get
$(L/M)\cdot(M/N)\geq 1$. For all the points in Fig. 5 for which this condition
is not satisfied, the system is unobservable. If the condition is satisfied,
then in each simulation run, a random measurement configuration is verified to
provide an observable system, thus the rank insufficiency may only appear as a
consequence of the C-RAN topology and the channel matrix sparsification.
According to Fig. 5, for the parameters used in our simulations, we can see
that GBP generally performs well, however, in the region where
$(L/M)\cdot(M/N)$ is slightly above $1$, rank insufficiency may deteriorate
the performance.

Overall, for systems with large number of MTC UEs $M$, the state can be
estimated with relatively small number of RRHs $L$. In contrast, for the
scenario where a number of MTC UEs is small, the number of RRHs must be
increased for successful reconstruction. In addition, simulation results point
to capability of the proposed scheme to provide successful reconstruction if
the underlying system is observable (i.e., of full rank), while the accuracy
of reconstruction (i.e., the accuracy of the state estimator) will depend on
the parameters such as channel sparsification, SNR, measurement standard
deviations and number of MTC UEs and RRHs. We leave detailed study of these
inter-dependencies for our future work.

##  V Conclusions

Motivated by the development of 5G massive MTC and large-scale distributed 5G
C-RAN architecture, in this paper, we proposed a scalable and efficient linear
state estimation framework. The proposed framework is based on the GBP
algorithm and jointly combines linear state estimation with signal detection
in 5G C-RANs. The advantage of GBP solution is accuracy that matches the MMSE
estimation, low complexity due to lack of scheduling MTC-UE transmissions, low
latency due to simultaneous data transfer, scalability to large-scale systems
(due to the fact that the underlying factor graph is usually sparse), and ease
of parallelization and distributed implementation in future distributed F-RAN
architectures. For the future work, we aim to provide rigorous convergence
analysis of GBP in the presented framework, motivated by similar analysis in
[7] and [9], and provide extensive numerical simulation study.

## References

  * [1] A. Rico-Alvarino, M. Vajapeyam, H. Xu, X. Wang, Y. Blankenship, J. Bergman, T. Tirronen, and E. Yavuz, “An overview of 3GPP enhancements on machine to machine communications,” _IEEE Communications Magazine,_ vol. 54, no. 6, pp. 14-–21, June 2016.
  * [2] A. Checko, H. L. Christiansen, Y. Yan, L. Scolari, G. Kardaras, M. S. Berger, and L. Dittmann, “Cloud RAN for mobile networks – a technology overview,” _IEEE Communications Surveys Tutorials,_ vol. 17, no. 1, pp. 405–-426, 2015.
  * [3] F. Bonomi, R. Milito, J. Zhu, and S. Addepalli, “Fog computing and its role in the internet of things,” in _Proc. ACM Workshop on Mobile cloud computing (MCC),_ pp. 13–-16, 2012.
  * [4] H. Shariatmadari, R. Ratasuk, S. Iraji, A. Laya, T. Taleb, R. Juntti, and A. Ghosh, “Machine-type communications: current status and future perspectives toward 5G systems,” _IEEE Communications Magazine,_ vol. 53, no. 9, pp. 10-–17, September 2015.
  * [5] M. Cosovic, A. Tsitsimelis, D. Vukobratovic, J. Matamoros and C. Anton-Haro, ”5G Mobile Cellular Networks: Enabling Distributed State Estimation for Smart Grids,” in _IEEE Communications Magazine_ , vol. 55, no. 10, pp. 62-69, October 2017.
  * [6] C. Fan, Y. Zhang, X. Yuan,“Advances and challenges towards scalable cloud radio access network,” _IEEE Communications Magazine,_ vol. 54, no. 6, pp. 29–35, June 2016.
  * [7] C. Fan, Y. Zhang, X. Yuan,“Scalable uplink signal detection in C-RANs via Randomized Gaussian Message Passing,” _IEEE Trans. Wireless Communications,_ vol. 16, no. 8, pp. 5187–5200, August 2017.
  * [8] A. Monticelli, “Electric power system state estimation,” _Proceedings of the IEEE,_ vol. 88, no. 2, pp. 262–-282, 2000.
  * [9] M. Cosovic, D. Vukobratovic, “Distributed Gauss-Newton Method for AC State Estimation Using Belief Propagation,” submitted, arXiv: https://arxiv.org/abs/1702.05781v2.
  * [10] H. A. Loeliger, J. Dauwels, J. Hu, S. Korl, L. Ping, and F. R. Kschischang, “The factor graph approach to model-based signal processing,” _Proc. of the IEEE,_ vol. 95, no. 6, pp. 1295–1322, 2007.
  * [11] Lien S.-Y., Hung S.-C., Chen K.-C., Liang Y.-C., “Ultra-low latency ubiquitous connections in heterogeneous cloud radio networks,” _IEEE Wireless Communications Magazine,_ vol. 22, no. 3, pp. 22–31, Jun. 2015.