# A Reinforcement Learning Approach for Wildfire Tracking with UAV Swarms

Carles Diaz-Vilor,1 Angel Lozano,2 and Hamid Jafarkhani1 1C. Diaz-Vilor and H.
Jafarkhani are with the Center for Pervasive Communications and Computing,
Univ. of California, Irvine. Their work was supported in part by NSF Award
CNS-2209695. Emails: {cdiazvil<EMAIL_ADDRESS>2A. Lozano is with Universitat
Pompeu Fabra (UPF), Barcelona. His work was supported by MICIU/AEI
/10.13039/501100011033 under the Maria de Maeztu Units of Excellence Programme
(CEX2021-001195-M), by ICREA, and by the Fractus-UPF Chair on 6G. Email:
<EMAIL_ADDRESS>
###### Abstract

Suitably equipped with cameras and sensors, uncrewed aerial vehicles (UAVs)
can be instrumental for wildfire prediction, tracking, and monitoring,
provided that uninterrupted connectivity can be guaranteed even if some of the
ground access points (APs) are damaged by the fire itself. A cell-free network
structure, with UAVs connecting to a multiplicity of APs, is therefore ideal
in terms of resilience. This work proposes a trajectory optimization framework
for a UAV swarm tracking a wildfire while maintaining cell-free connectivity
with ground APs. Such optimization entails a constant repositioning of the
multiplicity of UAVs as the fire evolves to ensure that the best possible view
is acquired and transmitted reliably, while respecting altitude limits,
avoiding collisions, and proceeding to recharge batteries as needed. Given the
complexity and time-varying nature of this multi-UAV trajectory optimization,
reinforcement learning is leveraged, specifically the twin-delayed deep
deterministic policy gradient algorithm. The approach is shown to be highly
effective for wildfire tracking and coverage and could be likewise applicable
to survey other natural and man-made phenomena, including weather events,
earthquakes, or chemical spills.

###### Index Terms:

UAV, wildfire, tracking, cell-free networks, reinforcement learning, TD3

##  I Introduction

Wildfires are a growing concern worldwide and cause major environmental and
economic damage. It is therefore crucial to predict, track, and monitor them
to facilitate the actions of firefighters and emergency responders [1, 2, 3].
Additionally, the study of wildfires can help to identify areas at risk and to
understand the nature of fire propagation, along with possible causes and
factors that contribute to their spread [4, 5, 6, 7].

An attractive means to monitor, track, and sense events is deploying uncrewed
aerial vehicles (UAVs) given the (_i_) availability of cameras and sensors,
(_ii_) wireless connectivity to ground access points (APs), and (_iii_) easy
control and maneuverability [8, 9, 10, 11, 12]. UAVs can fly at low altitudes
and collect data ranging from high-resolution images to heat signatures,
subsequently conveying those to a ground network of APs [13, 14, 15, 16, 17,
18]. We hasten to emphasize that UAVs surveying a wildfire act as end devices
in the network, a role that is dual to that of flying base stations in the
UAV-assistance paradigm [19, 20, 21].

The transmission of the gathered data is contingent on the connectivity.
Regardless of whether APs are damaged by fire, the flow of information should
be guaranteed. This relates to the network’s resiliency and, in that respect,
it is relevant that wireless systems are evolving from traditional cellular
structures towards cell-free arrangements. These are especially appealing when
a high degree of reliability is needed, given that users can then communicate
with multiple APs [22, 23, 24, 25, 26, 27, 28, 29, 30, 31]. Initial results
have demonstrated the efficacy and benefits of a cell-free architecture for
UAVs [32, 33, 34, 35, 36] and hence this is the structure adopted here.

The design of deployments and trajectories for UAVs in cellular and wireless
sensor networks is a problem of growing interest [37, 38, 39, 40, 41, 42, 43,
44, 45, 46, 47, 48, 49, 50, 51, 52, 53], yet, to the best of our knowledge,
this is the first work that considers a multi-UAV trajectory optimization to
track an event with cell-free connectivity. This gives rise to various
challenges, mainly related to the complexity and time-varying nature of the
problem, which precludes the use of classical optimization techniques.
Compounding this complexity, the lifespan of UAV batteries is limited and
thus, once a UAV is low on energy, a recharge is needed. Authors in [54]
provide an exhaustive review of different techniques to recharge or replace
UAV batteries, along with the required time. Altogether, the problem at hand
breaks down into two stages of _tracking_ and _charging_. Two decidedly
nonconvex optimization problems arise, one for each stage, with the switching
from tracking to charging based on the remaining energy at each UAV.

Reinforcement learning (RL) provides a framework where an agent learns optimal
policies based on interactions with the environment and on feedback in the
form of rewards [55]. For settings with a large number of states and/or
actions, strategies combining RL with neural networks have been devised,
converging to the so-called deep Q-learning (DQL) [56, 57, 58, 59, 60, 61].
These methods have already been considered for UAV trajectory optimizations
[62, 63, 64, 65, 66]. A recent actor-critic algorithm named twin-delayed deep
deterministic policy gradient (TD3) has been shown to perform better than its
predecessors [67] in terms of stability, exploration capabilities, handling of
continuous actions, and sampling efficiency. The field of UAV trajectory
optimization has greatly benefited from this new algorithm [68, 69, 70, 71],
which is applied in this work as well. The main contributions of the paper are
as follows:

  * •

An analytical framework is set forth to describe the tracking of a wildfire by
a swarm of UAVs equipped with cameras. The connectivity between UAVs and APs
is through a cell-free network, with Rician fading, channel estimation, and
minimum mean square error (MMSE) reception explicitly modeled.

  * •

Based on communication and mechanical constraints, the multi-UAV tracking and
charging optimization problems are formulated.

  * •

For each of the two problems, the UAV trajectory and the transmit power
optimizations are confronted using the TD3 approach. While, for the tracking
stage, the objective is to monitor the wildfire perimeter ensuring a correct
reception of the video/images in a timely fashion, during the charging stage
the goal is to reach a charging point with the minimum energy expenditure.

  * •

The impact on the wildfire coverage of parameters such as the number of UAVs,
the allowed flying altitudes, and the UAVs’ energy is established. A similar
analysis is conducted for the charging problem.

The remainder of the manuscript is organized as follows. Sec. II presents the
camera and communication models. In Sec. III, the cell-free connectivity is
introduced while Sec. IV formalizes the two optimization problems. Sec. V
subsequently focuses on the solution of those problems while numerical results
are discussed in Sec. VI. Concluding remarks are provided in Sec. VII.

$\displaystyle\mathcal{B}_{m}^{(n)}=\Big{\\{}\boldsymbol{v}=(v_{x},v_{y}):|x_{m}^{(n)}-v_{x}|\leq
h_{m}^{(n)}\tan(\alpha_{1})\text{ and }|y_{m}^{(n)}-v_{y}|\leq
h_{m}^{(n)}\tan(\alpha_{2})\Big{\\}}$ (1)

##  II System Model

The system under consideration features $M$ UAVs, each equipped with a video
camera. The $m$th UAV is located at
$\boldsymbol{q}_{m}^{(n)}=\big{(}x_{m}^{(n)},y_{m}^{(n)},h_{m}^{(n)}\big{)}$,
where the height $h_{m}^{(n)}$ is the distance from the focal point of its
camera’s lens to the ground, $n$ is the time index, and the duration of each
time slot is $\delta$. The UAVs are served by $L$ APs with locations
$\boldsymbol{q}_{\ell}=({x}_{\ell},{y}_{\ell},h_{\ell})$. Additionally, the
system contains $C$ UAV charging stations located at $\boldsymbol{c}_{c}=({\rm
x}_{c},{\rm y}_{c},{\rm h}_{c})$.

###  II-A Camera Model

Consider a planar environment $\mathcal{F}\in\mathbb{R}^{2}$ and define the
field of view (FoV) as the area projected over $\mathcal{F}$ that a camera
captures. For a UAV with a downward-facing camera, a rectangular FoV
$\mathcal{B}_{m}^{(n)}$ is defined in (1) atop the next page, where
$\alpha_{1}$ and $\alpha_{2}$ represent the two halfview angles associated
with the perpendicular edges of a rectangle (see Fig. 1). The notion of _area
per pixel_ expresses the tradeoff between the quality of the image and the
dimension of the FoV: higher-resolution images correspond to smaller FoVs, and
vice versa. From classical optics, the area per pixel is [72]

$\displaystyle
f\big{(}\boldsymbol{q}_{m}^{(n)},\boldsymbol{v}\big{)}=\left\\{\begin{array}[]{ll}a\,\left(b-h_{m}^{(n)}\right)^{2}&\boldsymbol{v}\in\mathcal{B}_{m}^{(n)}\\\
\infty&\text{otherwise}\\\ \end{array}\right.,$ (4)

where $b$ and $a$ depend on the camera capabilities, representing the focal
length of the lens and the area of a pixel on the lens divided by $b^{2}$,
respectively, while $\boldsymbol{v}$ represents the position on $\mathcal{F}$.
At higher altitudes, $f(\cdot)$ is larger because each pixel covers a bigger
area on the FoV, and the resolution is lower. At lower altitudes, it is the
other way around. Outside the FoV, by definition, there are no pixels and
$f(\cdot)$ is therefore infinite.

The camera’s sensor consists of $I$ equispaced pixels. From sheer geometry,

$\displaystyle
I=\frac{4\big{(}b-h_{m}^{(n)}\big{)}^{2}\tan(\alpha_{1})\tan(\alpha_{2})}{f\big{(}\boldsymbol{q}_{m}^{(n)},\boldsymbol{v}\big{)}}.$
(5)

For a 24-bit RGB color system, the number of bits to represent a picture is
$24I$ and, compressing the image with a ratio of $\rho$, the number of bits to
be transmitted per image is

$\displaystyle B$ $\displaystyle=24I\rho.$ (6)

Figure 1: Rectangular FoV of the $m$th UAV for given $\alpha_{1}$ and
$\alpha_{2}$ over a planar region $\mathcal{F}$.

###  II-B Channel Model

Upon imaging the environment, UAVs establish wireless connections with the APs
to convey the captured images. The air-to-ground channel coefficient between
the $m$th UAV and the $\ell$th AP is denoted by ${g}_{m,\ell}^{(n)}$, drawn
from a Rician distribution such that [73, Sec. 3.4.1]

$\displaystyle
g_{m,\ell}^{(n)}=\sqrt{\frac{\beta_{0}}{\big{(}d^{(n)}_{m,\ell}\big{)}^{\kappa}(K_{m,\ell}^{(n)}+1)}}\left(\sqrt{K_{m,\ell}^{(n)}}e^{j\psi_{m,\ell}^{(n)}}+{a}_{m,\ell}^{(n)}\right),$
(7)

where $\beta_{0}$ and $\kappa$ are, respectively, the pathloss at a reference
distance of $1$ m and the pathloss exponent; in turn, the distance is
$d_{m,\ell}^{(n)}$ and the Rician factor is

$K_{m,\ell}^{(n)}=A_{1}\exp\\!\left(A_{2}\arcsin\\!\left(\frac{h_{\ell}-h_{m}^{(n)}}{d_{m,\ell}^{(n)}}\right)\right)$
(8)

for environment-dependent parameters $A_{1}$ and $A_{2}$ [74]. The phase of
the LoS component, $\psi_{m,\ell}^{(n)}$, is uniformly random to reflect
drifting [29, 30] whereas the small-scale fading is

${a}_{m,\ell}^{(n)}\sim\mathcal{N}_{\mathbb{C}}(0,1).$ (9)

Hence, the channel power gain is given by

$\displaystyle{r}_{m,\ell}^{(n)}=\frac{\beta_{0}}{\big{(}d^{(n)}_{m,\ell}\big{)}^{\kappa}}.$
(10)

###  II-C Channel Acquisition

Unlike most of the UAV-related literature, this paper incorporates imperfect
channel estimates to enhance the model’s realism. Given that the system under
consideration features a small number of UAVs, orthogonal pilots can be
allocated to each of those UAVs while keeping the overhead at bay; pilot
contamination is thus not an issue. After observing the length-$\tau$ pilot
sequence transmitted by the $m$th UAV, the $\ell$th AP can compute the MMSE
channel estimate $\hat{g}_{m,\ell}^{(n)}$, which satisfies

${g}_{m,\ell}^{(n)}=\hat{g}_{m,\ell}^{(n)}+\tilde{g}_{m,\ell}^{(n)}$ (11)

such that $\hat{g}_{m,\ell}^{(n)}$ is zero-mean with [75, 76]

$\displaystyle\gamma_{m,\ell}^{(n)}$
$\displaystyle=\mathbb{E}\big{\\{}|\hat{g}_{m,\ell}^{(n)}|^{2}\big{\\}}$ (12)
$\displaystyle=\frac{\big{(}r_{m,\ell}^{(n)}\big{)}^{2}}{{r}_{m,\ell}^{(n)}+\frac{\sigma^{2}}{p^{\text{t}}\tau}},$
(13)

given $p^{\text{t}}$ as the pilot transmit power and $\sigma^{2}$ as the noise
power at the receiver. Finally, the estimation error
$\tilde{g}_{m,\ell}^{(n)}$ is zero-mean with variance

$c_{m,\ell}^{(n)}={r}_{m,\ell}^{(n)}-\gamma_{m,\ell}^{(n)}.$ (14)

###  II-D Energy Consumption Model

UAVs are governed by certain dynamics relating location, speed, and
acceleration. A first-order Taylor expansion yields [45],[46]

$\displaystyle\boldsymbol{q}_{m}^{(n+1)}=\boldsymbol{q}_{m}^{(n)}+\boldsymbol{v}_{m}^{(n)}\delta+\frac{1}{2}\boldsymbol{a}_{m}^{(n)}\delta^{2},$
(15)

and

$\displaystyle\boldsymbol{v}_{m}^{(n+1)}=\boldsymbol{v}_{m}^{(n)}+\boldsymbol{a}_{m}^{(n)}\delta,$
(16)

where $\boldsymbol{v}_{m}^{(n)}$ and $\boldsymbol{a}_{m}^{(n)}$ denote the
speed and acceleration vectors. For a quad-rotor UAV moving in 3D with a given
climb angle $\tau_{m}^{(n)}$, and denoting the magnitude of the speed vector
by $v_{m}^{(n)}$, i.e., $v_{m}^{(n)}=\|\boldsymbol{v}_{m}^{(n)}\|$, the power
consumption can be modeled as [50, 64]

$\displaystyle P_{m}^{(n)}$
$\displaystyle=P_{0}\Bigg{(}1+\frac{3v_{m}^{(n)2}}{U_{\rm
tip}^{2}}\Bigg{)}+P_{i}\Bigg{(}\sqrt{1+\frac{v_{m}^{(n)4}}{4v_{0}^{4}}}-\frac{v_{m}^{(n)2}}{2v_{0}^{2}}\Bigg{)}^{\\!1/2}$
$\displaystyle\quad+\frac{1}{2}d_{0}\iota
sAv_{m}^{(n)3}+mgv^{(n)}_{m}\sin\tau_{m}^{(n)},$ (17)

where $P_{0}$ and $P_{i}$ represent the blade profile and induced powers,
respectively, $U_{\rm tip}$ is the tip speed of the rotor blade, and $v_{0}$
is the mean rotor induced velocity. In addition, $d_{0}$ and $s$ are the
fuselage drag ratio and rotor solidity, respectively, while $\iota$ denotes
the air density and $A$ the rotor disc area.

As advanced, the UAV’s mission is divided into tracking and charging stages.
The switching occurs at $n=n_{m}$, once the UAV’s energy falls below a
threshold that is a function of the distance between the UAV and the closest
charging point (see Sec. VI). Consequently, the UAV needs to head for a
charging station and the remaining energy, denoted by $E_{m}^{(n_{m})}$, must
satisfy

$\displaystyle
E_{m}^{(n_{m})}-\\!\\!\\!\\!\sum\limits_{n=n_{m}}^{n_{m}+N_{c_{m}}}\\!\\!P_{m}^{(n)}\delta\geq
0,$ (18)

to avoid running out of energy, where $N_{c_{m}}$ is the number of slots it
takes for the $m$th UAV to reach charging station $\boldsymbol{c}_{c}$.

Other constraints include the maximum and minimum UAV altitudes, a maximum
velocity and acceleration, and the avoidance of collisions:

$\displaystyle h_{\rm min}\leq h_{m}^{(n)}$ $\displaystyle\leq h_{\rm max}$
(19) $\displaystyle{v}_{m}^{(n+1)}$ $\displaystyle\leq V_{\rm max}\delta$ (20)
$\displaystyle\|\boldsymbol{v}_{m}^{(n+1)}-\boldsymbol{v}_{m}^{(n)}\|$
$\displaystyle\leq A_{\rm max}\delta$ (21) $\displaystyle{\rm d}_{m,j}^{(n)}$
$\displaystyle\geq D_{\rm safe},$ (22)

where $V_{\rm max}$ and $D_{\rm safe}$ are the maximum UAV velocity and
minimum safety distance among UAVs, respectively, ${\rm d}_{m,j}^{(n)}$ is the
distance between the $m$th and $j$th UAVs, and $A_{\rm max}$ is the maximum
acceleration.

##  III Cell-free Connectivity

On a given time-frequency resource unit, the number of active APs is
$L^{(n)}$, which is a function of time. The uplink channel matrix is then

$\boldsymbol{{G}}^{(n)}=\begin{pmatrix}\boldsymbol{g}_{1}^{(n)},\dots,\boldsymbol{g}_{M}^{(n)}\end{pmatrix},$
(23)

where $\boldsymbol{g}_{m}^{(n)}\in\mathbb{C}^{L^{(n)}\times 1}$ is the channel
between the $m$th UAV and the active APs, satisfying

$\boldsymbol{{G}}^{(n)}=\boldsymbol{{\hat{G}}}^{(n)}+\boldsymbol{{\tilde{G}}}^{(n)},$
(24)

where $\boldsymbol{{\hat{G}}}^{(n)}$ and $\boldsymbol{{\tilde{G}}}^{(n)}$ are
the channel estimate and error matrices, respectively. In general, subsets of
APs can provide service to subsets of UAVs [26, 32, 36], yet this work
considers full connectivity, for the sake of the simplicity of notation,
pooling the observations from the $L^{(n)}$ APs into

$\displaystyle\boldsymbol{y}^{(n)}$
$\displaystyle=\boldsymbol{{G}}^{(n)}\boldsymbol{x}^{(n)}+\boldsymbol{n}^{(n)}$
(25)
$\displaystyle=\underbrace{\boldsymbol{{\hat{G}}}^{(n)}\boldsymbol{x}^{(n)}}_{\text{signal}}+\underbrace{\boldsymbol{{\tilde{G}}}^{(n)}\boldsymbol{x}^{(n)}+\boldsymbol{n}^{(n)}}_{\text{
effective noise: $\boldsymbol{n}_{\rm e}^{(n)}$ }},$ (26)

where

$\boldsymbol{x}^{(n)}=\left(\sqrt{p_{1}^{(n)}}s_{1}^{(n)},\dots,\sqrt{p_{M}^{(n)}}s_{M}^{(n)}\right)^{\\!\mathrm{T}},$
(27)

with unit power symbols $s_{m}^{(n)}$ while $p_{m}^{(n)}$ denotes the transmit
powers and the noise is
$\boldsymbol{n}^{(n)}\sim\mathcal{N}_{\mathbb{C}}(\boldsymbol{0},\sigma^{2}\boldsymbol{{I}})$.
The effective noise $\boldsymbol{n}_{\rm e}^{(n)}$ has covariance
$\boldsymbol{\Sigma}^{(n)}=\boldsymbol{{D}}^{(n)}+\sigma^{2}\boldsymbol{{I}}$
given

$\displaystyle\boldsymbol{{D}}^{(n)}$
$\displaystyle=\mathbb{E}\\!\left\\{\big{(}\boldsymbol{{\tilde{G}}}^{(n)}\boldsymbol{x}^{(n)}\big{)}\big{(}\boldsymbol{{\tilde{G}}}^{(n)}\boldsymbol{x}^{(n)}\big{)}^{*}\right\\}$
(28)
$\displaystyle=\text{diag}\\!\left\\{\sum\limits_{m=1}^{M}c_{m,1}^{(n)}p_{m}^{(n)},\dots,\sum\limits_{m=1}^{M}c_{m,L}^{(n)}p_{m}^{(n)}\right\\}.$
(29)

Upon observing $\boldsymbol{y}^{(n)}$, the combiner that maximizes the signal-
to-interference-plus-noise ratio (SINR) is the MMSE filter, achieving a value
of [26]

$\displaystyle\mathrm{SINR}_{m}^{(n)}$
$\displaystyle=\boldsymbol{\hat{g}}_{m}^{(n)*}\left(\sum\limits_{n\neq
m}\boldsymbol{\hat{g}}_{n}^{(n)}\boldsymbol{\hat{g}}_{n}^{(n)*}p_{n}^{(n)}+\boldsymbol{\Sigma}^{(n)}\\!\right)^{\\!\\!-1}\\!\\!\\!\boldsymbol{\hat{g}}_{m}^{(n)}\,p_{m}^{(n)},$
(30)

giving a spectral efficiency of

$\mathrm{SE}_{m}^{(n)}=\left(1-\frac{\tau}{\tau_{\rm
c}}\right)\mathbb{E}\\{\log_{2}(1+\mathrm{SINR}_{m}^{(n)})\\},$ (31)

where $\tau/\tau_{\rm c}$ accounts for the pilot overhead and $\tau_{c}$
denotes the total number of resource units within a coherence block. Random
matrix theory offers a way to circumvent the numerical evaluations of the
above expectation, providing stable spectral efficiency forms that depend
solely on large-scale parameters.

###  III-A Large-Dimensional Analysis

When evaluating (31) for $M,L^{(n)}\to\infty$, convergence to nonrandom limits
is assured provided that

$\displaystyle\boldsymbol{\Gamma}_{m}^{(n)}$
$\displaystyle=\mathbb{E}\Big{\\{}\boldsymbol{\hat{g}}_{m}^{(n)}\boldsymbol{\hat{g}}_{m}^{(n)*}\Big{\\}}$
(32)
$\displaystyle=\mathrm{diag}\big{\\{}\gamma_{m,\ell}^{(n)}\;\forall\ell\big{\\}}$
(33)

satisfies some technical conditions. Specifically, the inverse of the
resolvent matrix in (30) must exist, which is ensured by the presence of
$\boldsymbol{\Sigma}^{(n)}$, while $\boldsymbol{\Gamma}_{m}^{(n)}$ must have
uniformly bounded spectral norm, meaning that the received power does not
concentrate on a subset of dimensions as the network grows large.

######  Theorem 1.

For $M,L^{(n)}\to\infty$ with MMSE reception,
$\mathrm{SINR}_{m}^{(n)}-\overline{\mathrm{SINR}}_{m}^{(n)}\to 0$ almost
surely (a.s.) with

$\overline{\mathrm{SINR}}_{m}^{(n)}=\sum\limits_{\ell=1}^{L^{(n)}}\frac{\gamma_{m,\ell}^{(n)}p_{m}^{(n)}}{\sum\limits_{i\neq
m}\frac{\gamma_{i,\ell}^{(n)}}{1+e_{i}}p_{i}^{(n)}+\sum\limits_{\forall
i}c_{i,\ell}^{(n)}p_{i}^{(n)}+\sigma^{2}}.$ (34)

The coefficients $e_{j}$ are obtained iteratively with
$e_{j}=\lim_{s\to\infty}e_{j}^{(s)}$, $e_{j}^{(0)}=L^{(n)}$, and

$e_{j}^{(s)}=p_{j}^{(n)}\,\mathrm{tr}\Bigg{[}\boldsymbol{\Gamma}_{j}^{(n)}\bigg{(}\sum\limits_{i\neq
j}^{M}\frac{\boldsymbol{\Gamma}_{i}^{(n)}}{1+e_{i}^{(s-1)}}\,p_{i}^{(n)}+\boldsymbol{\Sigma}^{(n)}\bigg{)}^{\\!-1}\Bigg{]}.$
(35)

###### Proof.

See Appendix B. ∎

As anticipated, (34) is a deterministic quantity that depends only on large-
scale parameters. From the continuous mapping theorem [77], it follows that

$\mathrm{SE}_{m}^{(n)}-\mathrm{\overline{SE}}_{m}^{(n)}\stackrel{{\scriptstyle
a.s.}}{{\to}}0,$ (36)

where

$\mathrm{\overline{SE}}_{m}^{(n)}=\left(1-\frac{\tau}{\tau_{\rm
c}}\right)\log_{2}\\!\left(1+\overline{\mathrm{SINR}}_{m}^{(n)}\right)$ (37)

is the asymptotic spectral efficiency.

##  IV Problem Formulation

Let us now turn to optimizing the UAV trajectories and the transmit powers. In
the tracking stage, the UAVs in the set $\mathcal{S}_{1}^{(n)}$ have
sufficient energy to track while, in the charging stage, the UAVs in the set
$\mathcal{S}_{2}^{(n)}$ attempt to reach one of the $C$ charging points as
expeditiously as possible. Given the different nature of the two stages,
separate optimizations are devised.

###  IV-A Tracking Stage

During the tracking stage, UAVs periodically capture images, to be forwarded
to the APs every $N$ slots. To keep the information fresh on average, the
received image bit content must satisfy

$\displaystyle\delta W\mathrm{\overline{SE}}_{m}^{(n)}\geq\frac{B}{N}$ (38)

where $B$ is given in (6) and $W$ is the transmission bandwidth. Note that
this conservatively assumes that successive transmissions involve completely
new images.

With (38) as a constraint that adds to those on altitude, speed, and
collision-avoidance, the optimization of trajectories and powers can be cast
as the minimization of (a weighted version of) the area per pixel. As
discussed in Sec. II-A, the area per pixel reflects the tradeoff between image
quality and FoV. Before proceeding, though, the definition of area per pixel
must be extended from a single UAV to multiple UAVs, as their FoV may
intersect. If multiple cameras cover a certain point $\boldsymbol{v}$, the
area per pixel is [72]

$\displaystyle\tilde{f}\Big{(}\boldsymbol{q}_{m}^{(n)}:m\in\mathcal{S}_{1}^{(n)},\boldsymbol{v}\Big{)}=\frac{1}{\sum_{i\in\mathcal{N}_{\boldsymbol{v}}}f\big{(}\boldsymbol{q}_{i}^{(n)},\boldsymbol{v}\big{)}^{-1}+\Delta},$
(39)

where $\mathcal{N}_{\boldsymbol{v}}$ is the set of UAVs whose cameras cover
$\boldsymbol{v}$ while $\Delta$ is a regularization factor that prevents the
function from diverging if none of the UAVs covers $\boldsymbol{v}$.
Therefore, the joint optimization of the UAVs in tracking mode emerges as

$\displaystyle\min_{\boldsymbol{q}_{m}^{(n)},p_{m}^{(n)}}$
$\displaystyle\int_{\mathcal{F}}\tilde{f}\Big{(}\boldsymbol{q}_{m}^{(n)}:m\in\mathcal{S}_{1}^{(n)},\boldsymbol{v}\Big{)}\psi(\boldsymbol{v},n)d\boldsymbol{v},$
(40) s.t. $\displaystyle p_{m}^{(n)}\leq p_{\rm max}$
$\displaystyle\eqref{eq:linloc},\eqref{eq:linvel},\eqref{ct:altitude}\text{-}\eqref{ct:collision_av},\eqref{ct:Rate}$

where $\psi:\mathbb{R}^{2}\times\mathbb{N}\rightarrow[0,\infty)$ is a density
that weights the points on the FoV, assigning them relative importance (see
Sec. IV-C). As reported in Sec. VI, the solution to the above problem also
maximizes the coverage, defined as the ratio between the points with positive
density that lie inside the UAVs’ FoV and the total number of points with
positive density:

$\displaystyle{\rm
Coverage}^{(n)}=\frac{\int_{\mathcal{F}}\mathbb{I}\\!\left\\{\psi(\boldsymbol{v},n)>0\right\\}\mathbb{I}\\!\left\\{\exists\mathcal{B}_{m}^{(n)}:\boldsymbol{v}\in\mathcal{B}_{m}^{(n)}\right\\}d\boldsymbol{v}}{\int_{\mathcal{F}}\mathbb{I}\\!\left\\{\psi(\boldsymbol{v},n)>0\right\\}d\boldsymbol{v}},$
(41)

where $\mathbb{I}\\{\cdot\\}$ is an indicator function whose value is one when
the argument’s condition is true and zero otherwise.

###  IV-B Charging Stage

During the charging stage, a UAV tries to reach one of the $C$ charging
stations. For the $m$th UAV, the closest one is

$\displaystyle
c_{m}=\arg\min_{c}\|\boldsymbol{q}_{m}^{(n_{m})}-\boldsymbol{c}_{c}\|,$ (42)

where $n_{m}$, recall, is the time slot in which the $m$th UAV switches from
tracking to charging. In an unconstrained setup, trivial straight-line
trajectories would be optimal, but the actual solution is subject to
constraints. In particular, choosing to minimize the energy consumption, the
UAV trajectories are optimized via

$\displaystyle\min_{\boldsymbol{q}_{m}^{(n)},N_{c_{m}}}$
$\displaystyle\mspace{4.0mu}\sum\limits_{n=n_{m}}^{n_{m}+N_{c_{m}}}P_{m}^{(n)}\delta,$
(43) $\displaystyle\textrm{s.t.}\quad\eqref{eq:linloc},\eqref{eq:linvel},$
$\displaystyle\eqref{ct:energy}\text{-}\eqref{ct:collision_av}$

where $N_{c_{m}}$ is introduced in Sec. II-D, $P_{m}^{(n)}$ is defined in
(II-D), and the constraints ensure sufficient energy to reach a charging
point, the avoidance of collisions, and the respect of altitude and velocity
limits. The switching between tracking and charging, recall, is controlled by
the remaining energy at the UAV, as furthered in Sec. VI.

###  IV-C Fire Dynamics

The density function $\psi(\boldsymbol{v},n)$ depends on the event being
tracked and can be generated either synthetically through a model or
constructed from actual data. While the solution presented in the next section
is valid for any generic density function, for wildfires this function should
encode certain features that model the spread of a fire. Wildfire simulation
is an active field of research in itself [6]. FARSITE, an established model
employed by government agencies, is adopted here [7]. It dictates that the
spread of every ignition originally follows an ellipse. Subsequently, the
points on the ellipse serve as new fronts as per Huygens principle [7], each
such front point growing as a new ellipse. The set of new front points is
defined by the combination of the local ellipses, i.e., their convex hull [7,
Fig. 1]. Because of the different local weather conditions, fuels, or terrain,
the dimensions of the local ellipses can be different and result in a new
front that is no longer elliptic. Let us define $U^{(n)}$ [m/s] and
$\theta_{\text{wind}}^{(n)}$ as the midflame wind’s speed and direction,
respectively, modelled as

$\displaystyle U^{(n)}$ $\displaystyle\sim|\mathcal{N}(U_{0},\sigma_{U})|$
(44) $\displaystyle\theta_{\text{wind}}^{(n)}$
$\displaystyle\sim\mathcal{N}(\overline{\theta}_{\text{wind}},\sigma_{\text{wind}}).$
(45)

Generating the ellipse corresponding to each front point requires its minor
and major axes, $2a^{(n)}$ and $2b^{(n)}$, respectively, with

$\displaystyle a^{(n)}$
$\displaystyle=\frac{1}{2\,\mathsf{LB}^{(n)}}\left(\mathsf{R}+\frac{\mathsf{R}}{\mathsf{HB}^{(n)}}\right)$
(46) $\displaystyle b^{(n)}$
$\displaystyle=\frac{1}{2}\left(\mathsf{R}+\frac{\mathsf{R}}{\mathsf{HB}^{(n)}}\right),$
(47)

where $\mathsf{R}$ [m/min] is the fire’s steady-state spreading rate and

$\displaystyle\mathsf{LB}^{(n)}$
$\displaystyle=0.936\,e^{0.2566U^{(n)}}+0.461\,e^{–0.1548U^{(n)}}–0.397$ (48)
$\displaystyle\mathsf{HB}^{(n)}$
$\displaystyle=\frac{\mathsf{LB}^{(n)}+\sqrt{[\mathsf{LB}^{(n)}]^{2}–1}}{\mathsf{LB}^{(n)}\,–\,\sqrt{[\mathsf{LB}^{(n)}]^{2}–1}}.$
(49)

Hence, at time $n+1$, the $i$th front point, denoted by
$\boldsymbol{z}_{i}^{(n)}\in\mathbb{R}^{2}$, generates the ellipse

$\displaystyle\boldsymbol{z}_{i}^{(n)}+\delta\begin{bmatrix}c_{x}^{(n)}\sin\theta_{\text{wind}}^{(n)}+a^{(n)}\cos\omega\\\
c_{y}^{(n)}\cos\theta_{\text{wind}}^{(n)}+b^{(n)}\sin\omega\end{bmatrix},$
(50)

where $0\leq\omega\leq 2\pi$ whereas $c_{x}^{(n)}$ and $c_{y}^{(n)}$ represent
the fire spreading gradients, respectively. As the exact FARSITE
implementation requires detailed information about the environment, and with a
view to a general formulation, we consider a simplified set of FARSITE
parameters that has been widely used in the literature [13, 14]. In this
simplified model [7]

$\displaystyle
c_{x}^{(n)}=c_{y}^{(n)}=\frac{\mathsf{R}}{2}\left(1–\frac{1}{\mathsf{HB}^{(n)}}\right).$
(51)

With that, a 2D histogram can be produced to serve as the density function
$\psi(\boldsymbol{v},n)$, assigning a nonzero weight to the points of the
fire’s perimeter being tracked. Fig. 2 illustrates the expansion of a wildfire
from an ignition point (black dot) using FARSITE. Assuming wind in the
direction of the red arrow, an ellipse with minor and major axes, $2a^{(n)}$
and $2b^{(n)}$, respectively, represents the next fire front as shown in Fig.
2b. Then, each blue point in the new fire front acts like a new ignition point
and the convex hull of the corresponding ellipses generates the new fire front
in Fig. 2d. Here, the density function $\psi(\boldsymbol{v},n)$ corresponds to
the perimeter, represented by the red curve.

Figure 2: Summary of FARSITE fire propagation model.

##  V Proposed Solution

###  V-A DQL Fundamentals

Solving (40) and (43) is challenging due to (_i_) the variety of constraints,
(_ii_) the lack of convexity with respect to the optimization variables, and
(_iii_) the time-varying nature of the problem. However, the optimization can
be formulated as a Markov decision process (MDP) problem, for which RL is an
appropriate solution. In general, an MDP is defined by the tuple
$(\mathcal{S},\mathcal{A},\mathcal{P},R,\gamma)$, where $\mathcal{S}$ and
$\mathcal{A}$ are the state and action spaces, respectively. In turn,
$\mathcal{P}$ denotes the state transition probability with
$P\\!\left(s^{(n+1)}|s^{(n)},a^{(n)}\right)$ being the probability of
transitioning from state $s^{(n)}$ to $s^{(n+1)}$ after action $a^{(n)}$.
Finally, $R$ is the reward function and $\gamma\in[0,1]$ denotes the discount
factor in such function. According to the MDP,

$\displaystyle
P\\!\left(s^{(n+1)}|s^{(0)},a^{(0)},\ldots,s^{(n)},a^{(n)}\right)=P\\!\left(s^{(n+1)}|s^{(n)},a^{(n)}\right)$

and

$\displaystyle
R\\!\left(s^{(n+1)}|s^{(0)},a^{(0)},\dots,s^{(n)},a^{(n)}\right)$
$\displaystyle=R\\!\left(s^{(n+1)}|s^{(n)},a^{(n)}\right)$
$\displaystyle=r\\!\left(s^{(n)},a^{(n)}\right),$ (52)

where $r(\cdot)$ is introduced as shorthand notation. The goal in RL is to
find a policy $\pi:\mathcal{S}\rightarrow\mathcal{A}$ that maximizes the
average discounted reward,

$\displaystyle\pi^{*}=\arg\max_{\pi}\mathbb{E}\\!\left\\{\sum\limits_{j=0}^{\infty}\gamma^{j}r\\!\left(s^{(j)},a^{(j)}\right)\bigg{|}\pi\right\\},$
(53)

with expectation over $s^{(0)}\sim P\\!\left(s^{(0)}\right)$,
$a^{(n)}\sim\pi\\!\left(\cdot|s^{(n)}\right)$ and $s^{(n+1)}\sim
P\\!\left(s^{(n+1)}|s^{(n)},a^{(n)}\right)$. In RL, the Q-function measures
the expected future reward when the system is in state $s$, performs action
$a$, and follows policy $\pi$, namely

$\displaystyle
Q^{\pi}(s,a)=\mathbb{E}_{\pi}\\!\left\\{\sum\limits_{j=0}^{\infty}\gamma^{j}r\left(s^{(j^{\prime})},a^{(j^{\prime})}\right)\bigg{|}s^{(n)}=s,a^{(n)}=a\right\\}.$
(54)

where $j^{\prime}=j+n+1$. An interesting property of the Q-function is that
every optimum policy $\pi^{*}$ achieves the optimum Q-function, i.e.,

$\displaystyle Q^{\pi^{*}\\!}(s,a)=Q^{*}(s,a),$ (55)

where $*$ denotes optimality. Hence, one can identify the optimum policy
through the Q-function. In fact, through dynamic programming and the Bellman
equation, learning the Q-function is possible. Under the temporal difference
(TD) learning solution, which combines the current estimate of the Q-function
with samples obtained from the environment, the estimates can be iteratively
updated following [55]

$\displaystyle
Q(s,a)\leftarrow(1-\alpha)Q(s,a)+\alpha\left[r(s,a)+\gamma\max_{a^{\prime}}Q(s^{\prime},a^{\prime})\right].$
(56)

Although (56) provides a straightforward rule to learn the Q-function, the
corresponding learning process is highly dependent on the state and action
space dimensions, as every pair $(s,a)$ needs to be explored many times.
Multi-UAV communication problems, with their complex and dynamic environments,
and high-dimensional state and action spaces, motivate the use of DQL. This
leverages deep neural networks to approximate the Q-value function, with a
slight abuse of the notation, as $Q(s,a)\approx Q(s,a;\boldsymbol{\theta})$
where $\boldsymbol{\theta}$ are the neural network parameters tuned to
minimize

$\displaystyle
L(\boldsymbol{\theta})=\mathbb{E}\\!\left\\{|y(s,a)-Q(s,a;\boldsymbol{\theta})|^{2}\right\\},$
(57)

where $y(s,a)$ is the target value defined as

$\displaystyle
y(s,a)=r(s,a)+\gamma\max_{a^{\prime}}Q(s^{\prime},a^{\prime};\boldsymbol{\theta}).$
(58)

Still, DQL suffers from a major limitation as it considers a discrete action
space. To consider continuous actions, we adopt TD3, a policy-based model-free
algorithm [67]. TD3 leverages an actor-critic architecture, with both modules
composed by neural networks. On the one hand, the critic aims at learning the
Q-function, which is exploited by the actor to know how beneficial an action
is. On the other hand, the actor’s role is to learn a policy that, as opposed
to traditional Q-learning, outputs continuous action values. Such policy is
parameterized by $\boldsymbol{\phi}$, the actor’s neural network parameters,
that are updated by the gradient computed as

$\displaystyle\nabla
J(\boldsymbol{\phi})=\mathbb{E}\\!\left\\{\nabla_{\boldsymbol{\phi}}\pi_{\boldsymbol{\phi}}(s)\
\nabla_{a}Q(s,a;\boldsymbol{\theta})|_{a=\pi_{\boldsymbol{\phi}}(s)}\right\\},$
(59)

where

$\displaystyle
J(\boldsymbol{\phi})=\mathbb{E}_{\pi}\\!\left\\{\sum\limits_{j=0}^{\infty}\gamma^{j}r\\!\left(s^{(j)},a^{(j)}\right)\right\\}.$
(60)

In addition, DQL algorithms typically include the so-called target networks.
These are copies of the original networks whose parameters remain frozen over
a number of iterations. Then, a soft-update rule is applied. Concretely, we
denote by $\boldsymbol{\theta}^{\prime}$ and $\boldsymbol{\phi}^{\prime}$ the
critic and actor neural network target parameters, respectively, with update
rules

$\displaystyle\boldsymbol{\theta}^{{}^{\prime}}\leftarrow(1-\tau_{\rm
T})\boldsymbol{\theta}^{{}^{\prime}}+\tau_{\rm T}\boldsymbol{\theta},$ (61)
$\displaystyle\boldsymbol{\phi}^{{}^{\prime}}\leftarrow(1-\tau_{\rm
T})\boldsymbol{\phi}^{{}^{\prime}}+\tau_{\rm T}\boldsymbol{\phi},$ (62)

where $\tau_{\rm T}$ controls the memory of the updates. The inclusion of
target networks addresses the moving target issue that arises when Q-values
and target values are obtained using the same neural network and are
constantly changing, creating a moving target for the learning; the
algorithm’s goal then keeps varying as the agent learns, hampering the finding
of a solution.

Relative to previous actor-critic methods, TD3 introduces improvements that
enhance the stability during learning:

  1. 1.

It features two critics, and respective target critic networks. This addresses
the overestimation of Q-values caused by the use in previous methods of the
maximum action value to approximate the maximum expected action value, as per
(58). In TD3, the agent selects the Q-target network with smaller Q-value to
construct the target value as

$\displaystyle
y(s,a)=r(s,a)+\gamma\min_{i=1,2}Q(s^{\prime},a^{\prime};\boldsymbol{\theta}_{i}^{\prime}).$
(63)

  2. 2.

To reduce the accumulation of residual errors, the policy network parameters
are updated less frequently than the Q-function parameters. Precisely, [67]
suggests a policy update for every two Q-function updates.

  3. 3.

To prevent overfitting to specific estimated values, TD3 adds clipped random
noise to the target actions,

$\displaystyle\hat{a}=\pi_{\boldsymbol{\phi}^{\prime}}(s)+\hat{\epsilon},$
(64)

where
$\hat{\epsilon}\sim\text{clip}\big{(}\mathcal{N}(0,\hat{\sigma}_{a}),-\hat{\epsilon}_{\rm
max},\hat{\epsilon}_{\rm max}\big{)}$, i.e., the values outside the interval
$[-\hat{\epsilon}_{\rm max},\hat{\epsilon}_{\rm max}]$ are clipped to the
interval edges.

TD3 is applicable to (40) and (43). Next, we define the states, actions, and
rewards for each problem with one agent coordinating a swarm of $M$ UAVs.

###  V-B TD3: Tracking Stage

To solve (40), single-agent or multi-agent algorithms can be employed. A
single-agent solution would coordinate the entire swarm and present a
significant challenge due to the exponential growth of state and action spaces
with the number of UAVs. Alternatively, large MDPs can be factored into
simpler ones, leading to simpler agents and distributed solutions [78, 79].
Accordingly, we use the multi-agent solution, where each UAV is controlled by
a different agent. The following is needed:

  * •

The state space, which includes information about (_i_) the UAV locations,
(_ii_) the portion of the image bit content that has not yet been transmitted,
and (_iii_) the state of the fire coverage. To quantify the latter, let us
compute the center of mass of the area per pixel as per (65).

$\displaystyle\big{(}{\rm x}_{\rm c}^{(n)},{\rm y}_{\rm
c}^{(n)}\big{)}=\frac{\int_{\mathcal{F}}\boldsymbol{v}\tilde{f}\Big{(}\boldsymbol{q}_{m}^{(n)}:m\in\mathcal{S}_{1}^{(n)},\boldsymbol{v}\Big{)}\psi(\boldsymbol{v},n)d\boldsymbol{v}}{\int_{\mathcal{F}}\tilde{f}\Big{(}\boldsymbol{q}_{m}^{(n)}:m\in\mathcal{S}_{1}^{(n)},\boldsymbol{v}\Big{)}\psi(\boldsymbol{v},n)d\boldsymbol{v}}$
(65)

If important points remain uncovered, the agent learns to move in that
direction. The untransmitted portion of the image bit content is included in
the state space as

$\displaystyle i_{m}^{(n)}=\max\bigg{\\{}0,\bigg{(}1-\frac{\delta
W\mathrm{\overline{SE}}_{m}^{(n)}}{B/N}\bigg{)}\bigg{\\}}.$ (66)

Altogether, the state of the $m$th UAV at time $n$ is

$\displaystyle s_{m}^{(n)}=\left\\{\frac{x_{m}^{(n)}}{\rm
S},\frac{y_{m}^{(n)}}{\rm S},\frac{h_{m}^{(n)}}{h_{\rm
max}},\frac{\boldsymbol{v}_{m}^{(n)}}{v_{\rm max}},\frac{{\rm x}_{\rm
c}^{(n)}}{\rm S},\frac{{\rm y}_{\rm c}^{(n)}}{\rm S},\frac{{\rm
d}_{m,j}^{(n)}}{C_{\rm d}},i_{m}^{(n)}\right\\},$ (67)

  * where ${\rm S}$ and $C_{\rm d}$ are positive normalization constants, $\boldsymbol{v}_{m}^{(n)}$ is the UAV’s speed and ${\rm d}_{m,j}^{(n)}$ is the distance between the $m$th and the $j$th UAVs.

  * •

The action space defines the set of actions that the UAVs can take and
consists of variations in the 3D location and transmit power, such that
$a_{m}^{(n)}\in\mathbb{R}^{4\times 1}$. Given the relationships established by
(15) and (16), the problem is solved with respect to the acceleration. The
first three actions relate to the acceleration variations in 3D space, each
constrained by a maximum and minimum value $[-A_{\rm max},A_{\rm max}]$. The
last action component sets the transmit power, adjusted to any value
satisfying $0\leq p_{m}^{(n)}\leq p_{\rm max}$.

  * •

The reward function evaluates the quality of the actions and motivates the
agent to take those leading to desirable outcomes. Not only does it capture
the cost function, but it satisfies the constraints. As per the reward shaping
technique, the reward observed by the $m$th UAV is [80]

$\displaystyle
r_{m}\\!\left(s_{m}^{(n)},a_{m}^{(n)}\right)=\sum\limits_{w=1}^{5}r_{m,w}\\!\left(s_{m}^{(n)},a_{m}^{(n)}\right),$
(68)

whose terms are detailed next. In particular,
$r_{m,1}\\!\big{(}s_{m}^{(n)},a_{m}^{(n)}\big{)}$ contains the variation in
the cost function. After performing action $a_{m}^{(n)}$, the state
transitions to $s_{m}^{(n+1)}$, whose cost is defined as

$\displaystyle c(s_{m}^{(n+1)})$
$\displaystyle=\int_{\mathcal{F}}\tilde{f}\Big{(}\boldsymbol{q}_{m}^{(n+1)}:m\in\mathcal{S}_{1}^{(n+1)},\boldsymbol{v}\Big{)}$
$\displaystyle\quad\cdot\psi(\boldsymbol{v},n+1)d\boldsymbol{v}.$ (69)

Thus, to motivate the UAVs to reduce the cost function,
$r_{m,1}\\!\left(s_{m}^{(n)},a_{m}^{(n)}\right)$ is defined as

$\displaystyle r_{m,1}\\!\left(s_{m}^{(n)},a_{m}^{(n)}\right)=K_{\rm
c}\left(1-\frac{c(s_{m}^{(n+1)})}{I_{\rm r}^{(n)}}\right)$ (70)

where $K_{\rm c}>0$ and $I_{\rm
r}^{(n)}=\frac{1}{\Delta}\int_{\mathcal{F}}\psi(\boldsymbol{v},n+1)d\boldsymbol{v}$.
Moreover, to avoid collisions among UAVs, we define

$\displaystyle
r_{m,2}\left(s_{m}^{(n)},a_{m}^{(n)}\right)=\left\\{\begin{array}[]{ll}-K_{\rm
coll}&\text{if $\exists$ }{\rm d}_{m,j}^{(n)}\leq D_{\rm safe}\\\
0&\text{otherwise}\\\ \end{array},\right.$ (73)

where $K_{\rm coll}>0$. Similarly, we define the reward associated with (38)
as

$\displaystyle
r_{m,3}\\!\left(s_{m}^{(n)},a_{m}^{(n)}\right)=\left\\{\begin{array}[]{ll}-K_{\rm
f}&\text{if }\frac{B}{N}\geq\delta W\mathrm{\overline{SE}}_{m}^{(n)}\\\
0&\text{otherwise}\\\ \end{array},\right.$ (76)

where $K_{\rm f}>0$. A similar expression is used to define $r_{m,4}(\cdot)$
and $r_{m,5}(\cdot)$, which penalizes actions taking the UAV out of limits and
exceeding a maximum velocity, respectively, with respective penalties of
$-K_{\rm h}$ and $-K_{\rm v}$.

The multi-agent TD3 algorithm used to train a single UAV during the tracking
stage is summarized in Alg. 1 whereas Fig. 3 presents its block diagram. Note
that the factorization of the larger problem into $M$ single-agent problems
allows for parallelization, leading to faster convergence and more efficient
learning. Given that the UAVs in tracking mode have the same objective, their
training will converge to the same models. Therefore, training can be
conducted on one UAV, while the others can either follow the learned policy up
to that point in time, perform random movements, or remain static. Besides the
initial network parameters, Alg. 1 accepts the number of episodes
$\mathrm{E}$, which refers to the number of wildfire realizations the system
will observe, as input. Additionally, a replay buffer of size $|\mathcal{M}|$
stores transitions of the type
$\left\\{s_{m}^{(n)},a_{m}^{(n)},r_{m}\\!\left(s_{m}^{(n)},a_{m}^{(n)}\right),s_{m}^{(n+1)}\right\\}$.
Finally, $F$ refers to the update frequency over the policy and target
networks.

###  V-C TD3: Charging Stage

Since the tracking and charging stages are optimized separately, two
independent models can be trained. The same ideas and algorithms described for
tracking can be applied for charging, with only some modifications needed to
address (43). Concretely, during the charging stage, only one UAV needs to be
coordinated. Hence, the states consist of

$\displaystyle s_{m}^{(n)}=\left\\{\frac{{\rm x}_{u_{m}}-x_{m}^{(n)}}{\rm
S},\frac{{\rm y}_{u_{m}}-y_{m}^{(n)}}{\rm S},\frac{{\rm
h}_{u_{m}}-h_{m}^{(n)}}{h_{\rm max}},{\rm
d}_{m,j}^{(n)},\frac{E_{m}^{(n)}}{E_{0}}\right\\},$ (77)

where, recall, the distance between the $m$th and $j$th UAVs is ${\rm
d}_{m,j}^{(n)}$ while $E_{0}$ is a positive constant normalizing the remaining
energy. While $r_{m,2}(\cdot)$, $r_{m,4}(\cdot)$ and $r_{m,5}(\cdot)$ remain
the same, $r_{m,3}\\!\left(s_{m}^{(n)},a_{m}^{(n)}\right)=0$ as the UAVs have
no data to forward. The main difference, though, is in the definition of
$r_{m,1}\\!\left(s^{(n)},a^{(n)}\right)$, in (78) atop the next page, which
includes the variation in the cost function:

  * •

When the charging point is reached, a positive reward $K_{\rm fin}$ is given.

  * •

The term $-K_{\rm e}P_{m}^{(n)}\delta$ aims at minimizing the energy whereas
the second component motivates the UAV to minimize the distance with the
charging point. If it is approached, the UAV receives a positive reward;
otherwise, the reward is negative. Finally, a negative reward of $K_{\rm en}$
is incurred if the remaining energy is exhausted.

$\displaystyle r_{m,1}(s_{m}^{(n)},a_{m}^{(n)})=\begin{cases}K_{\rm
fin}&\text{if $\boldsymbol{q}_{m}^{(n)}=\boldsymbol{c}_{c_{m}}$}\\\ -K_{\rm
e}P_{m}^{(n)}\delta+K_{\rm
d}\big{(}\|\boldsymbol{c}_{c_{m}}-\boldsymbol{q}_{m}^{(n)}\|-\|\boldsymbol{c}_{c_{m}}-\boldsymbol{q}_{m}^{(n+1)}\|\big{)}&\text{otherwise}\end{cases}.$
(78)

Alg. 1 applies to the charging stage with the above modifications in the
rewards. In addition, whereas the tracking stage does not have a final state,
the training of the charging model does reach a final state whenever
$\boldsymbol{q}_{m}^{(n)}=\boldsymbol{c}_{c_{m}}$. Thus, the while condition
in the loop can be changed from $n<N_{e}$ to
$\boldsymbol{q}_{m}^{(n)}\neq\boldsymbol{c}_{c_{m}}$.

Figure 3: TD3 block diagram.

Input: No. of episodes $\mathrm{E}$, initial policy parameters
${\boldsymbol{\phi}}$, critic parameters $\boldsymbol{\theta}_{1}$,
$\boldsymbol{\theta}_{2}$, empty memory replay buffer $\mathcal{M}$ and
parameter update frequency $F$.

Initialize
$\boldsymbol{\theta}_{1}^{{}^{\prime}}\leftarrow\boldsymbol{\theta}_{1}$,
$\boldsymbol{\theta}_{2}^{{}^{\prime}}\leftarrow\boldsymbol{\theta}_{2}$, and
$\boldsymbol{\phi}^{{}^{\prime}}\leftarrow\boldsymbol{\phi}$,

for  _$e=1,\dots,\mathrm{E}$_ do

Set $n=0$ and initialize UAV, AP and ignition point locations.

while  _ $n<N_{e}$_ do

Obtain the current state $s^{(n)}$.

Perform action
$a_{m}^{(n)}=\text{clip}\left(\pi_{\boldsymbol{\phi}}\big{(}s_{m}^{(n)}\big{)}+\epsilon\right)$,
with $\epsilon\sim\mathcal{N}(0,\sigma_{a})$, and observe $s_{m}^{(n+1)}$.

Compute the reward $r\\!\left(s_{m}^{(n)},a_{m}^{(n)}\right)$ as in (68).

Store
$\left\\{s_{m}^{(n)},a_{m}^{(n)},r_{m}\\!\left(s_{m}^{(n)},a_{m}^{(n)}\right),s_{m}^{(n+1)}\right\\}$
in $\mathcal{M}$.

Select a minibatch of $N_{\rm mem}$ experiences
$\left\\{s_{m}^{(n)},a_{m}^{(n)},r_{m}\\!\left(s_{m}^{(n)},a_{m}^{(n)}\right),s_{m}^{(n+1)}\right\\}$
from $\mathcal{M}$.

Compute target actions $\hat{a}_{m}^{(n)}$ as in (64).

Compute targets $y(s,a)$ as in (58).

Update the critics by minimizing
$\mathbb{E}\\{|y(s,a)-Q(s,a;\boldsymbol{\theta}_{i})|^{2}\\}$ for $i=1,2$.

if  _$n$ mod $F$_ then

Update $\boldsymbol{\phi}$ via (59).

Update target networks $\boldsymbol{\theta}_{1}^{{}^{\prime}}$,
$\boldsymbol{\theta}_{2}^{{}^{\prime}}$, and $\boldsymbol{\phi}^{{}^{\prime}}$
as in (61) and (62).

end if

$n=n+1$.

Optional: move the rest of UAVs.

end while

end for

Algorithm 1 TD3 - Tracking stage

##  VI Numerical Results

For the purpose of performance evaluation, a $300$-m2 environment is
considered with the simulation parameters in Table I. The camera parameters
are set based on [72, 53], emulating real devices, while the UAV and cell-free
parameters are borrowed from [24, 36, 81]. Aiming at simulating real
wildfires, the fire propagation components are chosen from [7, 13, 14] whereas
the TD3 parameters are configured as suggested in [68, 69, 70, 71].
Particularly, the variance of the noise added to the actions is
${\sigma}_{a}=\hat{\sigma}_{a}=0.1$, $\hat{\epsilon}_{\rm max}=0.5$, and the
parameter soft updates are handled with $\tau_{\rm T}=0.01$. Also, the
learning rate for actor and critics is set to $5\cdot 10^{-4}$ and $5\cdot
10^{-3}$, respectively. The normalization constants in the states and rewards
are obtained via cross-validation and summarized in Table II; other values
could work as well. The initial UAV and AP locations, charging point
coordinates, and fire ignition point are randomly generated. At each time
slot, the fire perimeter is updated according to Sec. IV-C whereas the UAV
locations are updated following the output of the actor block in the trained
TD3 algorithm. Finally, given that we solve two separate problems, we first
provide training and evaluation results assuming unlimited energy at UAVs.
Then, after discussing the charging performance, the combined results of
tracking and charging are presented.

TABLE I: Simulation parameters Description | Parameter | Value | Description | Parameter | Value
---|---|---|---|---|---
Camera half-view angles |  $\alpha_{1}$, $\alpha_{2}$ | 17.5∘, 13.125∘ | Maximum power | $p_{\rm max}$ | 100 mW
Camera parameters |  $a$, $b$ |  $10^{-6}$, 10 | Pathloss at $1$ m | $\beta_{0}$ | -30 dB
Compression factor | $\rho$ | 0.4 | Noise power | $\sigma^{2}$ | -96 dBm
No. of APs | $L$ | $10$ | Dense urban param. |  $A_{1}$, $A_{2}$ | 0, 6.4 dB
No. of charging points | $C$ | $4$ | Regularizer in (39) | $\Delta$ | $10^{-5}$
Carrier frequency | $f_{\rm c}$ | 2.4 GHz | Time slot | $\delta$ | 0.5 s
Maximum acceleration | $A_{\rm max}$ | 1 | Transmission bandwidth | $B$ | 10 MHz
Maximum speed | $V_{\rm max}$ | 20 | Learning rate | $\gamma$ | 0.85
Minimum and maximum altitude |  $H_{\rm min}$, $H_{\rm max}$ | 100 m, 150 m | Image forwarding period | $N$ | 2
Pilot sequence length | $\tau$ | 200 | No. of channel uses | $\tau_{\rm c}$ | 6250
Mean midflame wind speed | $U_{0}$ | 5 m/s | Minimum safety dist. | $D_{\rm safe}$ | 4 m
Midflame wind speed variance | $\sigma_{U}$ | 1 m/s | Mean wind direction | $\hat{\theta}_{\rm wind}$ | $\mathcal{U}[0,2\pi]$
Fire spreading rate | $\mathsf{R}$ | 35 m/min | Wind direction variance | $\sigma_{\rm wind}$ | 0.1
TABLE II: State and reward parameters Description | Parameter | Value
---|---|---
Normalization in (67), (77) | ${\rm S}$ | 300
Constant in (70) | $K_{\rm c}$ | 50
Penalty in (73) | $K_{\rm coll}$ | 100
Penalty in (76) | $K_{\rm f}$ | 15
Flying-out-of-limits penalty | $K_{\rm h}$ | 60
Normalization in (77) | $E_{0}$ | 12000
Rewards and constants in (78) |  $K_{\rm fin}$,$K_{\rm e}$, $K_{\rm d}$ | 200, 1, 0.1
Out-of-energy penalty | $K_{\rm en}$ | 100

For starters, let us validate the cell-free asymptotic derivations, as they
are leveraged in (40). To that end, Fig. 4(a) plots the average user spectral
efficiency for different $M$ and $L$; as one would expect, a smaller $M/L$,
i.e., more APs per UAV, yields a better spectral efficiency. In turn, Fig.
4(b) indicates that Thm. 1 provides an excellent approximation even in low-
dimensional systems, say $M=3$ and $L=6$.

(a)

(b)

Figure 4: (a) Spectral efficiency for different $M,L$; (b) validation of
Theorem 1.

Next, the attention shifts to the tracking stage. Fig. 5 assesses the TD3
algorithm’s training performance for $M=1,\dots,4$, where a 3-layer feed-
forward neural network is used, each layer containing 256 neurons. In the
initial 750 episodes, the replay buffer accumulates transitions to ensure
that, upon training commencement, every minibatch comprises samples with
minimal correlation, enhancing the learning process. Hence, learning truly
begins after episode 750, leading to an overall increase in the average reward
for all models. Once 3000 episodes are reached, rewards stabilize. For
subsequent evaluations, the models saved at the 6000-episode mark are
utilized. Fig. 5(b) evaluates the trained models over a single realization of
a wildfire with the cost function in (40). Despite the random initialization
of the UAV locations, the swarm efficiently repositions itself, ensuring full
perimeter coverage by $n=100$. However, as $n$ reaches 300, the perimeter
dimensions become too extensive to be managed by $M=1$ UAV, resulting in a
significant cost increase. A similar trend is observed for $M=2$ UAVs, but
$M=3$ and $M=4$ UAVs manage to maintain a satisfactory performance. Note how,
due to the regularization parameter $\Delta$, the cost function exhibits a
high variance. When all points are covered, the function yields a small value.
However, if even a single perimeter point remains uncovered, the function
experiences sharp spikes. The coverage metric defined in (41) then comes
handy, as a smoother complement to the cost function, and hence both values
are reported in the sequel.

(a)

(b)

Figure 5: For the tracking problem (a) average training reward; (b) cost
function value as a function of $n$.

Next, we measure the cumulative distribution function (CDF) of the cost
function and fire coverage for $H_{\rm min}=125$ and $H_{\rm max}=150$ at
$n=350$ over 1,000 independent wildfire events. As one would expect,
increasing the number of UAVs results in an increased coverage and a smaller
cost function. Particularly, the case with $M=4$ UAVs achieves a 90% coverage
with very high probability. Also, note that both figures consider $n=350$,
i.e., fires that have already expanded and are large. Although not shown for
the sake of conciseness, for $n=100$, meaning for fires that are still small,
90% coverage is achieved regardless of the number of UAVs.

(a)

(b)

Figure 6: At $n=350$ for $H_{\rm min}=125$ and $H_{\rm max}=150$ (a) cost
function; (b) coverage.

To quantify the impact of the flying altitude on the problem, we include Fig.
7. These figures measure the swarm’s area per pixel (blue) and coverage (red)
over 400 time steps and average the results across 1,000 different fire
realizations, considering only the tracking stage where UAVs have unlimited
energy. Fig. 7(a) corresponds to $H_{\rm min}=125$ and $H_{\rm max}=150$,
while Fig. 7(b) is for $H_{\rm min}=100$ and $H_{\rm max}=125$. As expected,
reducing the allowed flying altitudes decreases the FoV, leading to a
degradation in both the cost function and coverage. This degradation is
especially notable when the fire is large, say $n>300$, whereas, at the onset
of the event, lower flying altitudes suffice.

(a)

(b)

Figure 7: Cost and coverage for tracking and $(H_{\rm min},H_{\rm max})$: (a)
(125,150) (b) (100,125).

Next, our RL solution is tested against two benchmarks. The first one
distributes the UAVs uniformly at random throughout the region and is denoted
by U. The second one, denoted by G, places the UAVs around the ignition point
following a 2D Gaussian random distribution with covariance
$10\boldsymbol{I}$. Against these benchmarks, two initializations are
considered for the RL approach, again uniformly at random or Gaussian, denoted
respectively by U+RL and G+RL. Fig. 8 presents results for $M=3$ and $M=4$
with $H_{\rm min}=125$ and $H_{\rm max}=150$. The following is observed:

  * •

The U method performs poorly. An approximate value for the coverage achieved
by this method can be gauged from the expected coverage with $M$ UAVs,
computed as the ratio between the area covered by the UAV’s FoV and the total
area of the region, namely

$M\frac{2^{2}\,\mathbb{E}\\{h_{m}^{(n)}\\}^{2}\tan\alpha_{1}\tan\alpha_{2}}{S^{2}}.$
(79)

For $\mathbb{E}\\{h_{m}^{(n)}\\}=(H_{\rm max}+H_{\rm min})/2$, the values
obtained for $M=3$ and $M=4$ are, respectively, 0.199 and 0.265, consistent
with the coverage curves for the U plots.

  * •

The G benchmark performs well during the initial stages. However, as the fire
grows large, say $n>150$, both the area per pixel and coverage degrade
significantly. In contrast, our method manages to provide solid coverage under
both initializations. Starting with the UAVs around the ignition point,
corresponding to G+RL, ensures better coverage at first, but the U+RL method
rapidly relocates the UAVs to cover the fire effectively too.

(a)

(b)

Figure 8: Cost and coverage for tracking under different solutions: (a) $M=3$
(b) $M=4$.

Before merging the tracking and charging stages, we validate the TD3 solution
applied to the charging stage via (43), where the same three-layer feedforward
neural network is used. For training, episodes are limited to a maximum length
of 200. The maximum number of UAVs coincides with that of the tracking stage,
i.e., $M=4$. UAV locations are randomly initialized. During each episode, one
UAV intends to reach a charging point, while the remaining UAVs perform random
movements. Fig. 9(a) plots the average training reward versus the episode
number for various values of $M$ with $\gamma=0.85$. After 3,000 episodes, the
agent is fully trained, achieving an average reward of around 200—the value of
$r_{\rm fin}$—thus resulting in
$\boldsymbol{q}_{m}^{(n)}=\boldsymbol{c}_{c_{m}}$. Once the agent is trained,
the energy threshold at which the UAV switches from tracking to charging needs
to be set. For that purpose, a 2D histogram is produced relating the distance
between the UAV and the charging point when switching to the charging stage
(x-axis), and the amount of energy that the UAV needs to reach the charging
point (y-axis). Fig. 9(b) shows such a histogram after evaluating the trained
agent over 50,000 episodes. Since the UAV should reach the charging point with
very high probability, a conservative threshold should be chosen. For example,
the UAV switches to charging when, at a certain distance, its level of energy
falls below $1.2$ times the maximum energy at that distance in Fig. 9(b).

(a)

(b)

Figure 9: For the charging problem (a) average training reward; (b) 2D
histogram for the distance vs energy in a trained agent.

Finally, for the combined tracking and charging results, the UAV energies are
initialized at random from $\mathcal{U}[3.5,E_{\rm max}]$ kJ for different
$E_{\rm max}$. This is justified by the UAVs’ previous surveillance resulting
in different energy values. The minimum energy level for UAVs is set to 3.5
kJ, ensuring that UAVs can reach a charging point in case their initial energy
levels are low. We measure the swarm’s coverage over 400 time steps and
average the results across 1,000 different fire realizations. Fig. 10 plots
the cost function (blue) and coverage (red) given $E_{\rm max}=125$ kJ for
various $M$ and $(H_{\rm min},H_{\rm max})$. While the initial coverage is
small, UAVs rapidly adjust their locations to track the varying perimeter.
Although, at each time instant, we aim at minimizing the cost function, the
corresponding curves are increasing. This is due to the fire perimeter
expanding with time, with a higher cost even if full coverage is achieved.
Also, note that the values attained in Fig. 10 are similar to those in Fig. 7,
where infinite energy is assumed, only with a 3-5% degradation.

To further assess the effects of finite energy batteries, Fig. 11 presents
results for lower values of $E_{\rm max}$, precisely for $E_{\rm max}=125$ kJ,
$E_{\rm max}=100$ kJ and $E_{\rm max}=80$ kJ, for a variety of $M$ with
$H_{\rm min}=125$ and $H_{\rm max}=150$. Interestingly, since the UAVs switch
to charging more often when their energy levels are lower, the performance
worsens. There is a 10% degradation in coverage for $M=3$ while, for $M=4$,
that degradation is smaller, especially when the fire is extensive.
Conversely, when a fire is in its initial stage, coverages above $90$% can be
achieved at all energy levels.

(a)

(b)

Figure 10: Cost and coverage for $E_{\rm max}=125$ kJ and $(H_{\rm min},H_{\rm
max})$: (a) (125,150) (b) (100,125).

(a)

(b)

Figure 11: Cost and coverage for different $E_{\rm max}$ (a) $M=3$, (b) $M=4$.

Fig. 12 illustrates the system’s operation for $M=3$ UAVs and a given wildfire
realization. Precisely, Figs. 12(a) and 12(b) depict the UAV trajectories
until $n=150$ and $n=400$, respectively, alongside the wildfire perimeter
represented by means of a 2D histogram, where the value of each point
represents its importance. For that same realization, the cost function and
coverage are depicted in Fig. 12(c) whereas Fig. 12(d) shows the excess bit
rate percentage in the cell-free connectivity constraint given by (38). At the
beginning of the mission, UAV-1 (red) needs to recharge its battery and
therefore heads to the closest charging point with coordinates ($225,80,10$).
In parallel, UAV-3 (blue) approaches the wildfire and, by $n=36$, full
coverage is achieved; see Fig. 12(c). However, the energy levels of UAV-3 are
low, and at $n=74$ it switches to charging mode and heads to the charging
point located at ($230,225,10$), causing a reduction in coverage. At that
time, UAV-1 (red) and UAV-2 (green) have enough energy and are in tracking
mode, repositioning themselves closer to the wildfire. At $n=84$, coverage
begins to increase, reaching full coverage by $n=140$. Afterwards, the three
UAVs continue adjusting their locations to cover the fire perimeter,
maintaining the coverage well above $80$% for the remainder of the mission.
Finally, Fig. 12(d) presents the excess bit rate in (38) as a percentage of
$\frac{B}{N}$. (A value of 0 is assigned during the charging periods as UAVs
do not transmit.). The combination of trajectory and power optimization yields
a positive gap, indicating that the images are correctly relayed to the
network.

(a)

(b)

(c)

(d)

Figure 12: For a given wildfire realization: (a) wildfire and UAV trajectories
until $n=150$; (b) wildfire and UAV trajectories until $n=400$; (c) cost
function and coverage over the entire mission; (d) excess bit rate in (38).

##  VII Conclusion

This paper has considered a cell-free UAV network whose aim is to track a
wildfire while satisfying a set of mechanical, energetic, and communication
constraints. Two complex nonconvex optimization problems have been formulated,
for tracking and charging, and a reinforcement learning framework has been
applied to tackle them. Particularly, the TD3 algorithm has been used.
Extensive results have shown that a small swarm of UAVs can reliably provide
coverage. Concretely, if the energy levels and the flying altitudes are
moderately high, an average coverage of more than 90% can be achieved with
only a few UAVs, with that coverage shrinking with the charging level and the
altitude. Altogether, the tradeoff among the number of UAVs, energy, and
flying altitude, has been established.

## Appendix A

######  Theorem 2.

([82, Thm. 1]) Let $\boldsymbol{{D}}\in\mathbb{C}^{M\times M}$ and
$\boldsymbol{{S}}\in\mathbb{C}^{M\times M}$ be Hermitian nonnegative-definite
while $\boldsymbol{{H}}\in\mathbb{C}^{M\times K}$ is a random matrix with
zero-mean independent vectors $\boldsymbol{h}_{k}$, each with covariance
$\frac{1}{M}\boldsymbol{R}_{k}$. In turn, $\boldsymbol{{D}}$ and
$\boldsymbol{{R}}_{k}$ have uniformly bounded spectral norm w.r.t. $M$. For
$z>0$ and $M,K\to\infty$,

$\frac{1}{M}\,\mathrm{tr}\\!\left[\boldsymbol{{D}}\big{(}\boldsymbol{{H}}\boldsymbol{{H}}^{*}+\boldsymbol{{S}}+z\boldsymbol{{I}}_{M})^{-1}\right]-\frac{1}{M}\,\mathrm{tr}[\boldsymbol{{D}}\boldsymbol{{T}}]\stackrel{{\scriptstyle\text{a.s.}}}{{\to}}0,$

where

$\boldsymbol{{T}}=\bigg{(}\frac{1}{M}\sum\limits_{j=1}^{K}\frac{\boldsymbol{{R}}_{j}}{1+e_{j}}+\boldsymbol{{S}}+z\boldsymbol{{I}}_{M}\bigg{)}^{\\!-1}$
(80)

with coefficients $e_{k}=\text{lim}_{n\xrightarrow{}\infty}e_{k}^{(n)}$ for

$e_{k}^{(n)}=\frac{1}{M}\,\mathrm{tr}\\!\left[\boldsymbol{{R}}_{k}\bigg{(}\frac{1}{M}\sum\limits_{j=1}^{K}\frac{\boldsymbol{{R}}_{j}}{1+e_{j}^{(n-1)}}+\boldsymbol{{S}}+z\boldsymbol{{I}}_{M}\bigg{)}^{\\!-1},\right]$
(81)

with initial values $e_{k}^{(0)}=M$.

## Appendix B

Let us drop the time index and define the matrix

$\boldsymbol{\Omega}_{m}=\left(\boldsymbol{{\hat{G}}}\boldsymbol{{P}}{}\boldsymbol{{\hat{G}}}^{*}-\boldsymbol{\hat{g}}_{m}\boldsymbol{\hat{g}}_{m}^{*}p_{m}+\boldsymbol{\Sigma}\right)^{\\!-1},$
(82)

with
$\boldsymbol{P}=\mathrm{diag}\\{p_{m}\mspace{4.0mu}\forall\mspace{4.0mu}m\\}$
and $\boldsymbol{\Omega}_{m}^{\prime}=L\boldsymbol{\Omega}_{m}$. Then, (30)
can be written as

$\displaystyle\mathrm{SINR}_{m}$
$\displaystyle=\boldsymbol{\hat{g}}_{m}^{*}\boldsymbol{\Omega}_{m}\boldsymbol{\hat{g}}_{m}\,p_{m}$
(83)
$\displaystyle=\frac{p_{m}}{L}\,\mathrm{tr}\\!\left[\boldsymbol{\hat{g}}_{m}\boldsymbol{\hat{g}}_{m}^{*}\boldsymbol{\Omega}_{m}^{\prime}\right].$
(84)

For $M$,$L$ $\xrightarrow{}\infty$, using [82, Lemma 4] and Theorem 2,

$\frac{p_{m}}{L}\,\mathrm{tr}\\!\left[\boldsymbol{\hat{g}}_{m}\boldsymbol{\hat{g}}_{m}^{*}\boldsymbol{\Omega}_{m}^{\prime}\right]-\frac{p_{m}}{L}\mathrm{tr}[\boldsymbol{\Gamma}_{m}\boldsymbol{{T}}_{m}]\stackrel{{\scriptstyle\text{a.s.}}}{{\to}}0.$
(85)

In our case, the role of
$\big{(}\boldsymbol{{H}}\boldsymbol{{H}}^{*}+\boldsymbol{{S}}+z\boldsymbol{{I}}_{M})^{\\!-1}$
in Theorem 2 is played by $\boldsymbol{\Omega}_{m}^{\prime}$. There is a
direct mapping between the terms in the aforementioned theorem and our
problem, namely (_i_) $\boldsymbol{{D}}=\boldsymbol{\Gamma}_{m}\,p_{m}$,
(_ii_) $\boldsymbol{{R}}_{j}=\boldsymbol{\Gamma}_{j}\,p_{j}$, and (_iii_)
$\boldsymbol{{S}}+z\boldsymbol{{I}}_{M}=\frac{1}{L}\boldsymbol{\Sigma}$ with
$\boldsymbol{{T}}_{m}$ following the structure of $\boldsymbol{{T}}$ in
Theorem 2,

$\boldsymbol{{T}}_{m}=\bigg{(}\frac{1}{L}\sum\limits_{j\neq
m}^{M}\frac{\boldsymbol{\Gamma}_{j}}{1+e_{j}}\,p_{j}+\frac{1}{L}\boldsymbol{\Sigma}\bigg{)}^{\\!-1}.$
(86)

The coefficients can be calculated as $e_{j}=\lim_{t\to\infty}e_{j}^{(t)}$
with

$\displaystyle e_{j}^{(t)}$
$\displaystyle=p_{j}\,\mathrm{tr}\Bigg{[}\boldsymbol{\Gamma}_{j}\bigg{(}\sum\limits_{i\neq
j}^{M}\frac{\boldsymbol{\Gamma}_{i}}{1+e_{i}^{(t-1)}}\,p_{i}+\boldsymbol{\Sigma}\bigg{)}^{\\!-1}\Bigg{]}.$
(87)

The fixed-point algorithm can be used to compute $e_{j}$ and has been proved
to converge [82]. Finally, since matrices $\boldsymbol{\Gamma}_{m}$ and
$\boldsymbol{{T}}_{m}$ are diagonal, (85) can be written as

$\displaystyle\frac{p_{m}}{L}\mathrm{tr}[\boldsymbol{\Gamma}_{m}\boldsymbol{{T}}_{m}]={p_{m}}\,\mathrm{tr}\\!\left[\boldsymbol{\Gamma}_{m}\bigg{(}\sum\limits_{i\neq
m}^{M}\frac{\boldsymbol{\Gamma}_{i}}{1+e_{i}}\,p_{i}+\boldsymbol{\Sigma}\bigg{)}^{\\!-1}\right],$
(88)

and, after some algebra, Theorem 1 is proved.

## References

  * [1] J. T. Abatzoglou and A. P. Williams, “Impact of anthropogenic climate change on wildfire across western US forests,” Proc. of the National Academy of Sciences, vol. 113, no. 42, pp. 11770–11775, 2016.
  * [2] E. Chuvieco et al., “Development of a framework for fire risk assessment using remote sensing and geographic information system technologies,” Ecological modelling, vol. 221, pp. 46–58, Jan. 2010.
  * [3] M. P. Thompson and D. E. Calkin, “Uncertainty and risk in wildland fire management: A review,” Journal of Environmental Management, vol. 92, pp. 1895–1909, Aug. 2011.
  * [4] J. Cohen, “The wildland-urban interface fire problem,” Fremontia, vol. 38, no. 2, p. 3, 2010.
  * [5] M. A. Finney et al., “A simulation of probabilistic wildfire risk components for the continental united states,” Stochastic Environmental Research and Risk Assessment, vol. 25, pp. 973–1000, Mar. 2011.
  * [6] G. Perry, “Current approaches to modelling the spread of wildland fire: a review,” Progress in Physical Geography: Earth and Environment, vol. 22, no. 2, p. 222–245, 1998.
  * [7] M. A. Finney, “FARSITE: Fire area simulator-model development and evaluation,” tech. rep., 1998.
  * [8] H. Jafarkhani, “Taking to the air to help on the ground: How UAVs can help fight wildfires,” in IEEE ComSoc Technology News, Oct. 2022.
  * [9] M. C. P. Santos et al., “A novel null-space-based UAV trajectory tracking controller with collision avoidance,” IEEE/ASME Trans. Mechatronics, vol. 22, pp. 2543–2553, Sept. 2017.
  * [10] K. Meng et al., “UAV trajectory and beamforming optimization for integrated periodic sensing and communication,” IEEE Wireless Commun. Lett., vol. 11, pp. 1211–1215, Mar. 2022.
  * [11] C. Zhan et al., “Energy-efficient data collection in UAV enabled wireless sensor network,” vol. 7, pp. 328–331, Nov. 2018.
  * [12] Y.-J. Chen et al., “Autonomous tracking using a swarm of UAVs: A constrained multi-agent reinforcement learning approach,” IEEE Trans. Vehicular Technology, vol. 69, pp. 13702–13717, Sept. 2020.
  * [13] E. Seraj and M. Gombolay, “Coordinated control of UAVs for human-centered active sensing of wildfires,” in 2020 American Control Conf. (ACC), pp. 1845–1852, Jul. 2020.
  * [14] K. A. Ghamry et al., “Cooperative control of multiple UAVs for forest fire monitoring and detection,” in IEEE Int’l Conf. on Mechatronic and Embedded Systems and Applications, Oct. 2016.
  * [15] K. A. Ghamry and Y. Zhang, “Fault-tolerant cooperative control of multiple UAVs for forest fire detection and tracking mission,” in Conf. on Control and Fault-Tolerant Systems (SysTol), pp. 133–138, 2016.
  * [16] E. Seraj, A. Silva, and M. Gombolay, “Multi-UAV planning for cooperative wildfire coverage and tracking with quality-of-service guarantees,” Autonomous Agents and Multi-Agent Systems, vol. 36, no. 2, p. 39, 2022.
  * [17] C. Yuan, Z. Liu, and Y. Zhang, “Fire detection using infrared images for UAV-based forest fire surveillance,” in 2017 Int’l Conf. on Unmanned Aircraft Systems (ICUAS), pp. 567–572, IEEE, 2017.
  * [18] C. Yuan et al., “UAV-based forest fire detection and tracking using image processing techniques,” in 2015 Int’l Conf. on Unmanned Aircraft Systems (ICUAS), pp. 639–643, IEEE, 2015.
  * [19] B. Alzahrani, O. S. Oubbati, A. Barnawi, M. Atiquzzaman, and D. Alghazzawi, “UAV assistance paradigm: State-of-the-art in applications and challenges,” Journal of Network and Computer Applications, vol. 166, p. 102706, 2020.
  * [20] O. S. Oubbati, H. Badis, A. Rachedi, A. Lakas, and P. Lorenz, “Multi-UAV assisted network coverage optimization for rescue operations using reinforcement learning,” in IEEE Consumer Communications & Networking Conference (CCNC), pp. 1003–1008, 2023.
  * [21] O. S. Oubbati, M. Atiquzzaman, A. Baz, H. Alhakami, and J. Ben-Othman, “Dispatch of UAVs for urban vehicular networks: A deep reinforcement learning approach,” IEEE Trans. Vehic. Techn., vol. 70, no. 12, pp. 13174–13189, 2021.
  * [22] S. Venkatesan, A. Lozano, and R. Valenzuela, “Network MIMO: Overcoming intercell interference in indoor wireless systems,” in Asilomar Conf. on Signals, Systems and Computers, pp. 83–87, 2007.
  * [23] H. Q. Ngo et al., “Cell-free massive MIMO versus small cells,” IEEE Trans. Wireless Commun., vol. 16, pp. 1834–1850, Mar. 2017.
  * [24] H. Q. Ngo et al., “On the total energy efficiency of cell-free massive MIMO,” IEEE Trans. Green Commun. and Net., vol. 2, pp. 25–39, Nov. 2018.
  * [25] E. Björnson and L. Sanguinetti, “Making cell-free massive MIMO competitive with MMSE processing and centralized implementation,” IEEE Trans. Wireless Commun., vol. 19, pp. 77–90, Jan. 2020.
  * [26] M. Attarifar, A. Abbasfar, and A. Lozano, “Subset MMSE receivers for cell-free networks,” IEEE Trans. Wireless Commun., vol. 19, pp. 4183–4194, Jun. 2020.
  * [27] R. Mosayebi, M. M. Mojahedian, and A. Lozano, “Linear interference cancellation for the cell-free C-RAN uplink,” IEEE Trans. Wireless Commun., vol. 20, pp. 1544–1556, Nov. 2021.
  * [28] E. Björnson and L. Sanguinetti, “Scalable cell-free massive MIMO systems,” IEEE Trans. Commun., vol. 68, pp. 4247–4261, Jul. 2020.
  * [29] C. D’Andrea et al., “Analysis of UAV communications in cell-free massive MIMO systems,” IEEE Open Journal of the Commun. Soc., vol. 1, pp. 133–147, Jan. 2020.
  * [30] Ö. T. Demir and E. Björnson, “Joint power control and LSFD for wireless-powered cell-free massive MIMO,” IEEE Trans. Wireless Commun., vol. 20, pp. 1756–1769, Mar. 2020.
  * [31] V. Tentu, E. Sharma, D. N. Amudala, and R. Budhiraja, “UAV-enabled hardware-impaired spatially correlated cell-free massive MIMO systems: Analysis and energy efficiency optimization,” IEEE Trans. Commun, vol. 70, no. 4, pp. 2722–2741, 2022.
  * [32] C. Diaz-Vilor, A. Lozano, and H. Jafarkhani, “Cell-free UAV networks: Asymptotic analysis and deployment optimization,” IEEE Trans. Wireless Commun., vol. 22, no. 5, pp. 3055–3070, 2023.
  * [33] C. Diaz-Vilor, A. Lozano, and H. Jafarkhani, “On the deployment problem in cell-free UAV networks,” in 2021 IEEE Global Commun. Conf., Dec. 2021.
  * [34] J. Zheng, J. Zhang, and B. Ai, “UAV communications with WPT-aided cell-free massive MIMO systems,” IEEE J. Sel. Areas Commun., pp. 3114–3128, Jan. 2021.
  * [35] J. An and F. Zhao, “Trajectory optimization and power allocation algorithm in MBS-assisted cell-free massive MIMO systems,” IEEE Access, vol. 9, pp. 30417–30425, 2021.
  * [36] C. Diaz-Vilor, A. Lozano, and H. Jafarkhani, “Cell-free UAV networks with wireless fronthaul: Analysis and optimization,” IEEE Trans. Wireless Commun., pp. 1–1, 2023.
  * [37] J. Guo, P. Walk, and H. Jafarkhani, “Optimal deployments of UAVs with directional antennas for a power-efficient coverage,” IEEE Trans. Commun., vol. 68, pp. 5159–5174, Aug. 2020.
  * [38] E. Koyuncu et al., “Deployment and trajectory optimization of UAVs: A quantization theory approach,” IEEE Trans. Wireless Commun., vol. 17, pp. 8531–8546, Dec. 2018.
  * [39] C. Diaz-Vilor and H. Jafarkhani, “Optimal 3D-UAV trajectory and resource allocation of DL UAV-GE links with directional antennas,” 2020 IEEE Global Commun. Conf., pp. 1–6, Dec. 2020.
  * [40] J. Guo and H. Jafarkhani, “Movement-efficient sensor deployment in wireless sensor networks with limited communication range,” IEEE Trans. Wireless Commun., vol. 18, pp. 3469–3484, July 2019.
  * [41] Y. Zeng, X. Xu, and R. Zhang, “Trajectory design for completion time minimization in UAV-enabled multicasting,” IEEE Trans. Wireless Commun., vol. 17, pp. 2233–2246, Apr. 2018.
  * [42] Q. Wu, Y. Zeng, and R. Zhang, “Joint trajectory and communication design for multi-UAV enabled wireless networks,” IEEE Trans. Wireless Commun., vol. 17, pp. 2109–2121, Mar. 2018.
  * [43] J. Guo and H. Jafarkhani, “Sensor deployment with limited communication range in homogeneous and heterogeneous wireless sensor networks,” IEEE Trans. Wireless Commun., vol. 15, pp. 6771–6784, Oct. 2016.
  * [44] C. You and R. Zhang, “3D trajectory optimization in Rician for UAV-enabled data harvesting,” IEEE Trans. Wireless Commun., vol. 18, pp. 3192–3207, Jun. 2019.
  * [45] J. Guo, E. Koyuncu, and H. Jafarkhani, “A source coding perspective on node deployment in two-tier networks,” IEEE Trans. Commun., vol. 66, pp. 3035–3049, Jul. 2018.
  * [46] F. Cheng et al., “UAV trajectory optimization for data offloading at the edge of multiple cells,” IEEE Trans. Vehicular Technology, vol. 67, pp. 6732–6736, Jul. 2018.
  * [47] S. Karimi-Bidhendi, J. Guo, and H. Jafarkhani, “Energy-efficient node deployment in heterogeneous two-tier wireless sensor networks with limited communication range,” IEEE Trans. Wireless Commun., vol. 20, no. 1, pp. 40–55, 2021.
  * [48] Y. Zeng and R. Zhang, “Energy-efficient UAV communication with trajectory optimization,” IEEE Trans. Wireless Commun., vol. 16, pp. 3747–3760, Jun. 2017.
  * [49] C. Sun, W. Ni, and X. Wang, “Joint computation offloading and trajectory planning for uav-assisted edge computing,” IEEE Trans. Wireless Commun., vol. 20, pp. 5343–5358, Mar. 2021.
  * [50] Y. Zeng, J. Xu, and R. Zhang, “Energy minimization for wireless communication with rotary-wing UAV,” IEEE Trans. Wireless Commun., vol. 18, pp. 2329–2345, Apr. 2019.
  * [51] S. Karimi-Bidhendi et al., “Energy-efficient deployment in static and mobile heterogeneous multi-hop wireless sensor networks,” IEEE Trans. Wireless Commun., vol. 21, no. 7, pp. 4973–4988, 2022.
  * [52] C. Diaz-Vilor et al., “Sensing and communication in UAV cellular networks: Design and optimization,” IEEE Trans. Wireless Commun., pp. 1–1, 2023.
  * [53] N. V. Cuong, Y.-W. P. Hong, and J.-P. Sheu, “UAV trajectory optimization for joint relay communication and image surveillance,” IEEE Trans. Wireless Commun., vol. 21, pp. 10177–10192, Jun. 2022.
  * [54] H. Ucgun et al., “A review on applications of rotary-wing unmanned aerial vehicle charging stations,” Int’l Journal of Advanced Robotic Systems, vol. 18, p. 17298814211015863, May 2021.
  * [55] R. S. Sutton and A. G. Barto, Reinforcement learning: An introduction, vol. 2.  MIT press, Nov. 2018.
  * [56] K. Arulkumaran, M. P. Deisenroth, M. Brundage, and A. A. Bharath, “Deep reinforcement learning: A brief survey,” IEEE Signal Processing Magazine, vol. 34, pp. 26–38, Nov. 2017.
  * [57] V. Mnih et al., “Playing atari with deep reinforcement learning,” arXiv preprint arXiv:1312.5602, 2013.
  * [58] R. S. Sutton, D. McAllester, S. Singh, and Y. Mansour, “Policy gradient methods for reinforcement learning with function approximation,” Advances in Neural Information Proc. Systems, vol. 12, 1999.
  * [59] V. Konda and J. Tsitsiklis, “Actor-critic algorithms,” Advances in Neural Information Proc. Systems, vol. 12, 1999.
  * [60] V. Mnih et al., “Asynchronous methods for deep reinforcement learning,” in Int’l Conf. on Machine Learning, pp. 1928–1937, PMLR, 2016.
  * [61] T. P. Lillicrap et al., “Continuous control with deep reinforcement learning,” arXiv preprint arXiv:1509.02971, 2015.
  * [62] H. Yang et al., “Privacy-preserving federated learning for UAV-enabled networks: Learning-based joint scheduling and resource management,” IEEE J. on Sel. Areas Commun., vol. 39, pp. 3144–3159, Jun. 2021.
  * [63] Y. Yu et al., “Multi-objective optimization for UAV-assisted wireless powered IoT networks based on extended DDPG algorithm,” IEEE Trans. Commun., vol. 69, pp. 6361–6374, Jun. 2021.
  * [64] R. Ding, F. Gao, and X. S. Shen, “3D UAV trajectory design and frequency band allocation for energy-efficient and fair communication: A deep reinforcement learning approach,” IEEE Trans. Wireless Commun., vol. 19, pp. 7796–7809, Aug. 2020.
  * [65] Z. Xia et al., “Multi-agent reinforcement learning aided intelligent UAV swarm for target tracking,” IEEE Trans. Vehicular Tech., vol. 71, pp. 931–945, Nov. 2022.
  * [66] R. Ding et al., “Trajectory design and access control for air–ground coordinated communications system with multiagent deep reinforcement learning,” IEEE IoT Journal, vol. 9, pp. 5785–5798, Feb. 2022.
  * [67] S. Fujimoto, H. Hoof, and D. Meger, “Addressing function approximation error in actor-critic methods,” in Int’l Conf. on Machine Learning, pp. 1587–1596, PMLR, 2018.
  * [68] M. Sun, X. Xu, X. Qin, and P. Zhang, “AoI-energy-aware UAV-assisted data collection for IoT networks: A deep reinforcement learning method,” IEEE IoT Journal, vol. 8, pp. 17275–17289, May 2021.
  * [69] D. Hong et al., “Energy-efficient online path planning of multiple drones using reinforcement learning,” IEEE Trans. Vehicular Tech., vol. 70, pp. 9725–9740, Aug. 2021.
  * [70] Y. Wang et al., “Trajectory design for UAV-based internet of things data collection: A deep reinforcement learning approach,” IEEE IoT Journal, vol. 9, pp. 3899–3912, Aug. 2022.
  * [71] Y. Li and A. H. Aghvami, “Radio resource management for cellular-connected UAV: A learning approach,” IEEE Trans. Commun., vol. 71, pp. 2784–2800, Mar. 2023.
  * [72] M. Schwager, B. J. Julian, M. Angermann, and D. Rus, “Eyes in the sky: Decentralized control for the deployment of robotic camera networks,” Proceedings of the IEEE, vol. 99, pp. 1541–1561, July 2011.
  * [73] R. W. Heath Jr. and A. Lozano, Foundations of MIMO communication.  Cambridge University Press, 2019.
  * [74] S. Shimamoto and Iskandar, “Channel characterization and performance evaluation of mobile communication employing stratospheric platforms,” IEICE Trans. Commun., vol. E89-B, pp. 937–944, Mar. 2006.
  * [75] A. Lozano, “Interplay of spectral efficiency, power and Doppler spectrum for reference-signal-assisted wireless communication,” IEEE Trans. Wireless Commun., vol. 7, no. 12, pp. 5020–5029, 2008.
  * [76] Q. Sun, D. Cox, A. Lozano, and H. Huang, “Training-based channel estimation for continuous flat fading BLAST,” in IEEE Int’l Conf. Commun. (ICC), vol. 1, pp. 325–329, 2002.
  * [77] H. B. Mann and A. Wald, “On stochastic limit and order relationships,” Annals of Mathematical Statistics, vol. 14, pp. 217–226, 1943.
  * [78] C. Guestrin, D. Koller, and R. Parr, “Multiagent planning with factored mdps,” in Advances in Neural Information Processing Systems (T. Dietterich, S. Becker, and Z. Ghahramani, eds.), vol. 14, MIT Press, 2001.
  * [79] T. Degris, O. Sigaud, and P.-H. Wuillemin, “Learning the structure of factored markov decision processes in reinforcement learning problems,” in Proc. of the 23rd Int’l Conf. on Machine Learning, (New York, NY, USA), p. 257–264, Association for Computing Machinery, 2006.
  * [80] A. Y. Ng et al., “Policy invariance under reward transformations: Theory and application to reward shaping,” in Proc. of the Sixteenth Int’l Conf. on Machine Learning, ICML ’99, (San Francisco, CA, USA), p. 278–287, Morgan Kaufmann Publishers Inc., 1999.
  * [81] “Enhanced LTE support for aerial vehicles,” Tech. Rep. 36.777, 3GPP, Dec. 2017.
  * [82] S. Wagner et al., “Large system analysis of linear precoding in correlated MISO broadcast channels under limited feedback,” IEEE Trans. Inf. Th., vol. 58, pp. 4509–4537, Mar. 2012.


# Validation of tight-binding model in system of two semiconductor single-
electron lines coupled electrostatically by Schrödinger formalism
*

Krzysztof Pomorski1,2 1: Institute of Physics, Lodz University of Technology,
Lodz, Poland
2: Quantum Hardware Systems, Lodz, Poland

###### Abstract

Validation of tight-binding model is given basing on two single-electron lines
coupled perturbatively electrostatically in Schrödinger formalism. Scheme for
conversion of quantum information from eigenenergy qubits to position based
qubits is given. The procedure for determination of system ground is
presented. Additional extension schemes of two coupled single-electron lines
in case of Rabi oscillations occurring in each position based qubits are
given.

###### Index Terms:

single-electron line, position-based qubits, quantum swap gate, Rabi
oscillations

##  I Introduction to position-based qubits

Currently there is worldwide race in development of quantum technologies
leading to new quantum computation platforms, quantum communication and
quantum metrology systems [1] , [2]. Due to integration possibility the
dominant technology platforms are mesoscopic superconducting structures and
nano-size semiconductor electronic elements as coupled semiconductor quantum
dots or miniaturized cryogenic CMOS circuits. Starting from eigenergy qubit as
given by IBM Q-Experience or many other technologies
($|\gamma_{E1}(t)|^{2}+|\gamma_{E2}(t)|^{2}=1$) we arrive into

$\displaystyle\hat{H}(t)_{E}\ket{\psi(t)}_{E}=\begin{pmatrix}E_{1\rightarrow
1}(t)&E_{2\rightarrow 1}(t)\\\ E_{1\rightarrow 2}(t)&E_{2\rightarrow
2}(t)\end{pmatrix}_{E}\begin{pmatrix}\gamma_{E1}(t)\\\
\gamma_{E2}(t)\end{pmatrix}_{E}=$
$\displaystyle=i\hbar\frac{d}{dt}\begin{pmatrix}\gamma_{E1}(t)\\\
\gamma_{E2}(t)\end{pmatrix}_{E}=i\hbar\frac{d}{dt}(\ket{\psi(t)})_{E}=$
$\displaystyle=i\hbar\frac{d}{dt}(\gamma_{E1}(t)\ket{E_{1}(t)}+\gamma_{E2}(t)\ket{E_{2}(t)})=$
$\displaystyle=\begin{pmatrix}E_{1\rightarrow 1}(t)&E_{2\rightarrow 1}(t)\\\
E_{1\rightarrow 2}(t)&E_{2\rightarrow
2}(t)\end{pmatrix}_{E}\begin{pmatrix}+\cos(\Theta)&-\sin(\Theta)\\\
+\sin(\Theta)&+\cos(\Theta)\end{pmatrix}\times$
$\displaystyle\times\begin{pmatrix}+\cos(\Theta)&+\sin(\Theta)\\\
-\sin(\Theta)&+\cos(\Theta)\end{pmatrix}\begin{pmatrix}\gamma_{E1}(t)\\\
\gamma_{E2}(t)\end{pmatrix}_{E}=$ $\displaystyle
i\hbar\frac{d}{dt}\Bigg{[}\begin{pmatrix}+\cos(\Theta)&-\sin(\Theta)\\\
+\sin(\Theta)&+\cos(\Theta)\end{pmatrix}\times$
$\displaystyle\times\begin{pmatrix}+\cos(\Theta)&+\sin(\Theta)\\\
-\sin(\Theta)&+\cos(\Theta)\end{pmatrix}\begin{pmatrix}\gamma_{E1}(t)\\\
\gamma_{E2}(t)\end{pmatrix}_{E}\Bigg{]}$ (1)

where we have $E_{s\rightarrow k}(t)=\bra{E_{k}(t)}\hat{H}(t)\ket{E_{s}(t)}$
since we preassume $\hat{H}(t)=E_{1\rightarrow
1}(t)\ket{E_{1}(t)}\bra{E_{1}(t)}+E_{2\rightarrow
2}(t)\ket{E_{2}(t)}\bra{E_{2}(t)}+E_{1\rightarrow
2}(t)\ket{E_{2}(t)}\bra{E_{1}(t)}+E_{2\rightarrow
1}(t)\ket{E_{1}(t)}\bra{E_{2}(t)}$. Here $|\gamma_{E1}(t)|^{2}$ and
$|\gamma_{E2}(t)|^{2}$ accounts for probability of occupancy of $E_{1}$ and
$E_{2}$ energetic states. Maintenance of quantum information is about keeping
quantum information encoded in Bloch sphere by two parameters ($\Theta(t)$,
$\phi(t)$) and given in reference to quantum state by
$\ket{\psi}_{E}=\gamma_{E1}(t)\ket{E_{1}(t)}+\gamma_{E2}(t)\ket{E_{2}(t)}=\\\
=e^{i\phi_{0}(t)}(cos(\Theta(t))\ket{E_{1}(t)}+e^{i\phi(t)}sin(\Theta(t))\ket{E_{2}(t)})$.
Due to the fact of compact size of semiconductor quantum dots as depicted in
Fig.1. it is vital to express position based qubit by stating its Hamiltonian
in so called tight-binding model, where quantum state is linear superposition
of two Wannier wave functions denoting maximal localization on the left or
right quantum dot ($|\alpha_{q}(t)|^{2}+|\beta_{q}(t)|^{2}=1$):

$\displaystyle\hat{H}(t)_{W}\ket{\psi(t)}_{W}=\begin{pmatrix}E_{p1}(t)&t_{s12}(t)\\\
t_{s12}^{*}(t)&E_{p2}(t)\end{pmatrix}_{W}\begin{pmatrix}\alpha_{q}(t)\\\
\beta_{q}(t)\end{pmatrix}_{W}=$ $\displaystyle
i\hbar\frac{d}{dt}(\alpha_{q}(t)\ket{w_{L}(t)}+\beta_{q}(t)\ket{w_{R}(t)})_{W}=i\hbar\frac{d}{dt}\ket{\psi(t)}_{W},$
(2)

where $\ket{w_{L}(t)}$ and $\ket{w_{R}(t)}$ are maximum localized Wannier wave
functions
($\bra{w_{L}}\ket{w_{R}}=0$,$\bra{w_{L}}\ket{w_{L}}=1=\bra{w_{R}}\ket{w_{R}}$
) respectively on left and right side, while $|\alpha_{q}(t)|^{2}$ and
$|\beta_{q}(t)|^{2}$ gives the probabilities for particle (single electron or
hole) to stay on the left and right quantum dot. Encoding quantum information
by electron or hole position brings compactification of quantum chips. Quantum
state in Wannier representation $\ket{\psi(t)}_{W}$ can be mapped to quantum
state in energetic representation given as $\ket{\psi(t)}_{E}$ or reversely
and those two ways of description are fully equivalent.

Following [4], [5] it will be useful to introduce parameter

$\displaystyle
r=\frac{\int_{-\infty}^{0}dx(\psi_{E1}(x)\psi_{E2}^{{\dagger}}(x)+\psi_{E1}(x)^{{\dagger}}\psi_{E2}(x))}{\int_{-\infty}^{0}dx(|\psi_{E1}(x)|^{2}-|\psi_{E2}(x)|^{2})}.$
(3)

The precondition for occurrence of position-based qubit is occupancy at least
two energetic levels, so neither $|\gamma_{E1}|^{2}$ nor $|\gamma_{E2}|^{2}$
is zero. Finally in case of static Hamiltonian from equation 1 also given by
[4] the tight-binding Hamiltonian [5] will appear to be

$\displaystyle\begin{pmatrix}+\cos(\Theta)&\sin(\Theta)\\\
-\sin(\Theta)&\cos(\Theta)\end{pmatrix}\begin{pmatrix}E_{1}&0\\\
0&E_{2}\end{pmatrix}_{E}\begin{pmatrix}+\cos(\Theta)&-\sin(\Theta)\\\
-\sin(\Theta)&+\cos(\Theta)\end{pmatrix}=\begin{pmatrix}E_{p1}&t_{s21}\\\
t_{s12}&E_{p2}\end{pmatrix}_{W}=$
$\displaystyle=E_{p1}\ket{w_{L}}\bra{w_{L}}+E_{p2}\ket{w_{R}}\bra{w_{R}}+t_{s21}\ket{w_{R}}\bra{w_{L}}+t_{s12}\ket{w_{L}}\bra{w_{R}}=\hat{H}_{W}=$
$\displaystyle\begin{pmatrix}E_{1}+|sin(\frac{1}{2}ArcTan(r))|^{2}(E_{2}-E_{1})&(E_{2}-E_{1})\frac{1}{2}sin(ArcTan(r))\\\
(E_{2}-E_{1})\frac{1}{2}sin(ArcTan(r))&E_{1}+|cos(\frac{1}{2}ArcTan(r))|^{2}(E_{2}-E_{1})\end{pmatrix},$
(4)

We propose maximum localized orthonormal Wannier functions (so occupancy of
$w_{L}$ is maximized from minus infinity to zero) of the form

$\displaystyle\Theta=\frac{1}{2}ArcTan(r),\ket{\psi}_{W}=\begin{pmatrix}w_{L}(x)\\\
w_{R}(x)\end{pmatrix}_{W}=\begin{pmatrix}+\cos(\Theta)&\sin(\Theta)\\\
-\sin(\Theta)&\cos(\Theta)\end{pmatrix}\begin{pmatrix}\psi_{L}(x)\\\
\psi_{R}(x)\end{pmatrix}_{E}=\begin{pmatrix}+\cos(\Theta)&\sin(\Theta)\\\
-\sin(\Theta)&\cos(\Theta)\end{pmatrix}\ket{\psi}_{E},$
$\displaystyle\begin{pmatrix}+cos(\frac{1}{2}ArcTan[\frac{\int_{-\infty}^{0}dx(\psi_{E1}(x)\psi_{E2}^{{\dagger}}(x)+\psi_{E1}(x)^{{\dagger}}\psi_{E2}(x))}{\int_{-\infty}^{0}dx(|\psi_{E1}(x)|^{2}-|\psi_{E2}(x)|^{2})}])&sin(\frac{1}{2}ArcTan[\frac{\int_{-\infty}^{0}dx(\psi_{E1}(x)\psi_{E2}^{{\dagger}}(x)+\psi_{E1}(x)^{{\dagger}}\psi_{E2}(x))}{\int_{-\infty}^{0}dx(|\psi_{E1}(x)|^{2}-|\psi_{E2}(x)|^{2})}])\\\
-sin(\frac{1}{2}ArcTan[\frac{\int_{-\infty}^{0}dx(\psi_{E1}(x)\psi_{E2}^{{\dagger}}(x)+\psi_{E1}(x)^{{\dagger}}\psi_{E2}(x))}{\int_{-\infty}^{0}dx(|\psi_{E1}(x)|^{2}-|\psi_{E2}(x)|^{2})}])&cos(\frac{1}{2}ArcTan[\frac{\int_{-\infty}^{0}dx(\psi_{E1}(x)\psi_{E2}^{{\dagger}}(x)+\psi_{E1}(x)^{{\dagger}}\psi_{E2}(x))}{\int_{-\infty}^{0}dx(|\psi_{E1}(x)|^{2}-|\psi_{E2}(x)|^{2})}])\end{pmatrix}\begin{pmatrix}\psi_{E1}(x)\\\
\psi_{E2}(x)\end{pmatrix}_{E}.$

##  II Two interacting single-electron lines described in tight-binding model
validated by perturbative Schrödinger model

We consider the system of two position based qubits interacting
perturbatively. Two electrons on two separated single-electron lines tends to
be in anti-correlated positions due to repulsion force as depicted in Fig.1.
In such a way we encounter the situation of two electrons placed on two
different left (right) quantum dots with probabilities $|\xi_{1}(t)|^{2}$
($|\xi_{4}(t)|^{2}$) or the case left (right) occupancy of first dot and right
(left) occupancy of second dot with probability $|\xi_{2}(t)|^{2}$
($|\xi_{3}(t)|^{2}$) that specifies normalized quantum state
($|\xi_{1}(t)|^{2}+..+|\xi_{4}(t)|^{2}=1$) to be in the form as given by [3]:

$\displaystyle\ket{\psi}=\xi_{1}(t)\ket{0}_{A}\ket{0}_{B}+\xi_{2}(t)\ket{0}_{A}\ket{1}_{B}+\xi_{3}(t)\ket{1}_{A}\ket{0}_{B}+\xi_{4}(t)\ket{1}_{A}\ket{1}_{B}.$
(6)

In case of 2 Single-Electron Lines we can write the effective Hamiltonian to
be of the form

$\displaystyle\hat{H}_{A-B}=\hat{H}_{A}\times\hat{I}_{B}+\hat{I}_{A}\times\hat{H}_{B}+\hat{H}_{Coulomb:A-B}=$
$\displaystyle=(E_{p1A}\ket{w_{L,A}}\bra{w_{L,A}}+E_{p2A}\ket{w_{R,A}}\bra{w_{R,A}}+t_{s12A}\ket{w_{R,A}}\bra{w_{L,A}}+t_{s21A}\ket{w_{L,A}}\bra{w_{R,A}})\times\hat{I}_{B}+$
$\displaystyle+\hat{I}_{A}\times(E_{p1B}\ket{w_{L,B}}\bra{w_{L,B}}+E_{p2B}\ket{w_{R,B}}\bra{w_{R,B}}+t_{s12B}\ket{w_{R,B}}\bra{w_{L,B}}+t_{s21B}\ket{w_{L,B}}\bra{w_{R,B}})+$
$\displaystyle+\Bigg{[}q_{1}\ket{w_{L,A}}\ket{w_{L,B}}\bra{w_{L,A}}\bra{w_{L,B}}+q_{2}\ket{w_{L,A}}\ket{w_{L,B}}\bra{w_{R,A}}\bra{w_{R,B}}+q_{3}\ket{w_{L,A}}\ket{w_{L,B}}\bra{w_{R,A}}\bra{w_{R,B}}+$
$\displaystyle+q_{4}\ket{w_{R,A}}\ket{w_{R,B}}\bra{w_{R,A}}\bra{w_{R,B}}\Bigg{]}_{Coulomb}$
(7)

Finally we obtain validation of tight-binding model in terms of Schrödinger
equation, so we end up with following effective Hamiltonian

$\displaystyle\hat{H}=$
$\displaystyle\begin{pmatrix}E_{1a}cos(\Theta_{A})^{2}+E_{2a}sin(\Theta_{A})^{2}&(E_{2b}-E_{1b})\frac{1}{2}sin(2\Theta_{B})&(E_{2a}-E_{1a})\frac{1}{2}sin(2\Theta_{A})&0\\\
(E_{2b}-E_{1b})\frac{1}{2}sin(2\Theta_{B})&E_{1a}cos(\Theta_{A})^{2}+E_{2a}sin(\Theta_{A})^{2}&0&(E_{2a}-E_{1a})\frac{1}{2}sin(2\Theta_{A})\\\
(E_{2a}-E_{1a})\frac{1}{2}sin(2\Theta_{A})&0&E_{1a}sin(\Theta_{A})^{2}+E_{2a}cos(\Theta_{A})^{2}&(E_{2b}-E_{1b})\frac{1}{2}sin(2\Theta_{B})\\\
0&(E_{2a}-E_{1a})\frac{1}{2}sin(2\Theta_{A})&(E_{2b}-E_{1b})\frac{1}{2}sin(2\Theta_{B})&E_{1a}sin(\Theta_{A})^{2}+E_{2a}cos(\Theta_{A})^{2}\\\
\end{pmatrix}+$
$\displaystyle\begin{pmatrix}E_{1b}cos(\Theta_{B})^{2}+E_{2b}sin(\Theta_{B})^{2}&0&0&0\\\
0&E_{1b}sin(\Theta_{B})^{2}+E_{2b}cos(\Theta_{B})^{2}&0&0\\\
0&0&E_{1b}cos(\Theta_{B})^{2}+E_{2b}sin(\Theta_{B})^{2}&0\\\
0&0&0&E_{1b}sin(\Theta_{B})^{2}+E_{2b}cos(\Theta_{B})^{2}\end{pmatrix}+$
$\displaystyle+\begin{pmatrix}q_{1}&0&0&0\\\ 0&q_{2}&0&0\\\ 0&0&q_{3}&0\\\
0&0&0&q_{4}\end{pmatrix}$ (8)

where 4 terms are responsible for 2 qubits being at a distance d mutual
Coulomb interaction (as depicted in Fig.1.) and given as

$\displaystyle
q_{1}=\int_{-\infty}^{+\infty}\int_{-\infty}^{+\infty}\frac{dx_{A}dx_{B}e^{2}}{4\pi\epsilon_{0}\sqrt{d^{2}+(x_{A}-x_{B})^{2}}}|cos(\Theta_{A})\psi_{E1}(x_{A})+sin(\Theta_{A})\psi_{E2}(x_{A})|^{2}*|cos(\Theta_{B})\psi_{E1}(x_{B})+sin(\Theta_{B})\psi_{E2}(x_{B})|^{2}$
(9) $\displaystyle
q_{2}=\int_{-\infty}^{+\infty}\int_{-\infty}^{+\infty}\frac{dx_{A}dx_{B}e^{2}}{4\pi\epsilon_{0}\sqrt{d^{2}+(x_{A}-x_{B})^{2}}}|cos(\Theta_{A})\psi_{E1}(x_{A})+sin(\Theta_{A})\psi_{E2}(x_{A})|^{2}*|-sin(\Theta_{B})\psi_{E1}(x_{B})+cos(\Theta_{B})\psi_{E2}(x_{B})|^{2}$
(10) $\displaystyle
q_{3}=\int_{-\infty}^{+\infty}\int_{-\infty}^{+\infty}\frac{dx_{A}dx_{B}e^{2}}{4\pi\epsilon_{0}\sqrt{d^{2}+(x_{A}-x_{B})^{2}}}|-sin(\Theta_{A})\psi_{E1}(x_{A})+cos(\Theta_{A})\psi_{E2}(x_{A})|^{2}*|cos(\Theta_{B})\psi_{E1}(x_{B})+sin(\Theta_{B})\psi_{E2}(x_{B})|^{2}$
(11) $\displaystyle
q_{4}=\int_{-\infty}^{+\infty}\int_{-\infty}^{+\infty}\frac{dx_{A}dx_{B}e^{2}}{4\pi\epsilon_{0}\sqrt{d^{2}+(x_{A}-x_{B})^{2}}}|-sin(\Theta_{A})\psi_{E1}(x_{A})+cos(\Theta_{A})\psi_{E2}(x_{A})|^{2}|-sin(\Theta_{B})\psi_{E1}(x_{B})+cos(\Theta_{B})\psi_{E2}(x_{B})|^{2}.$

In first approximation we consider two-position based qubits as non-
interacting and slightly perturbed by electrostatic single-electron to single-
electron interaction what naturally leads to first approximation of system
eigenstate as tensor product of two Hilbert spaces belonging to each qubits in
the form as given below:

$\displaystyle\ket{\psi}_{Energy,A-B}=\begin{pmatrix}\gamma_{E1A}\\\
\gamma_{E2A}\\\ \end{pmatrix}\times\begin{pmatrix}\gamma_{E1B}\\\
\gamma_{E2B}\\\
\end{pmatrix}=\begin{pmatrix}\gamma_{E1A}\gamma_{E1B}=cos(\xi_{A})e^{i\phi_{E1A}}cos(\xi_{B})e^{i\phi_{E1B}}\\\
\gamma_{E1A}\gamma_{E2B}=cos(\xi_{A})e^{i\phi_{E1A}}sin(\xi_{B})e^{i\phi_{E2B}}\\\
\gamma_{E2A}\gamma_{E1B}=sin(\xi_{A})e^{i\phi_{E2A}}cos(\xi_{B})e^{i\phi_{E1B}}\\\
\gamma_{E2A}\gamma_{E2B}=sin(\xi_{A})e^{i\phi_{E2A}}sin(\xi_{B})e^{i\phi_{E2B}}\\\
\end{pmatrix}.$ (13)

Search for system ground state naturally implies minimization of Coulomb
energy as criteria for lowest energy state that can be expressed by proper
self-alignment of electrons in anticorrelated way and thus leading to
minimization of Coulomb interaction energy. Coulomb energy is given as

$\displaystyle E_{c}(\xi_{A},\xi_{B})=\bra{\psi}\hat{H}_{c}\ket{\psi}=$
$\displaystyle=q_{1}cos(\xi_{A})^{2}cos(\xi_{B})^{2}+q_{2}cos(\xi_{A})^{2}sin(\xi_{B})^{2}+q_{3}sin(\xi_{A})^{2}cos(\xi_{B})^{2}+q_{4}sin(\xi_{A})^{2}sin(\xi_{B})^{2}.$
(14)

Coulomb energy minimization condition brings following constrains

$\displaystyle\frac{d}{d\xi_{A}}E_{c}(\xi_{A},\xi_{B})=0,\frac{d}{d\xi_{B}}E_{c}(\xi_{A},\xi_{B})=0,$
(15)

and hence we have

$\displaystyle
0=\frac{d}{d\xi_{A}}E_{c}(\xi_{A},\xi_{B})=2(-sin(\xi_{A})[q_{1}cos(\xi_{B})^{2}+q_{2}sin(\xi_{B})^{2}]+$
$\displaystyle+cos(\xi_{A})[q_{3}cos(\xi_{B})^{2}+q_{4}sin(\xi_{B})^{2}]),$
(16)

that leads to

$\displaystyle
tan(\xi_{A})=\frac{[q_{3}cos(\xi_{B})^{2}+q_{4}sin(\xi_{B})^{2}]}{[q_{1}cos(\xi_{B})^{2}+q_{2}sin(\xi_{B})^{2}]}=\frac{[q_{3}+q_{4}tan(\xi_{B})^{2}]}{[q_{1}+q_{2}tan(\xi_{B})^{2}]}.$
(17)

Now we exercise second constrain on energy minimization what brings

$\displaystyle
0=\frac{d}{d\xi_{B}}E_{c}(\xi_{A},\xi_{B})=2(+cos(\xi_{A})^{2}[-q_{1}sin(\xi_{B})+q_{2}cos(\xi_{B})]+$
$\displaystyle+sin(\xi_{A})^{2}[-q_{3}sin(\xi_{B})+q_{4}cos(\xi_{B})]),$ (18)

and we obtain

$\displaystyle
tan(\xi_{A})^{2}=\frac{[-q_{1}sin(\xi_{B})+q_{2}cos(\xi_{B})]}{[q_{3}sin(\xi_{B})-q_{4}cos(\xi_{B})]}=\frac{[q_{3}cos(\xi_{B})^{2}+q_{4}sin(\xi_{B})^{2}]^{2}}{[q_{1}cos(\xi_{B})^{2}+q_{2}sin(\xi_{B})^{2}]^{2}}.$
(19)

Under assumption that qubit A is occupying 2 energy levels $E_{1A}$ and
$E_{2A}$, while qubit B is occupying 2 energy levels $E_{1B}$ and $E_{2B}$ we
have

$\displaystyle
tan(\xi_{A})^{2}=\frac{[-q_{1}tan(\xi_{B})+q_{2}]}{[q_{3}tan(\xi_{B})-q_{4}]}=\frac{[q_{3}+q_{4}tan(\xi_{B})^{2}]^{2}}{[q_{1}+q_{2}tan(\xi_{B})^{2}]^{2}},$
(20)

$|\gamma_{E1A}|^{2}=cos(\xi_{A})^{2}=p_{E1A}$,
$|\gamma_{E2A}|^{2}=cos(\xi_{A})^{2}=p_{E1A}$,
$|\gamma_{E1B}|^{2}=cos(\xi_{B})^{2}=p_{E1B}$ and
$|\gamma_{E2B}|^{2}=sin(\xi_{B})^{2}=p_{E2B}$. Geometry of two
electrostatically interacting structures is encoded in $q_{1}$, $q_{2}$,
$q_{3}$ and $q_{4}$ coefficients. We compute $(r_{A}$,$r_{A}$) in accordance
with formula 3 for qubits (A,B) and thus consequently we obtain
$\Theta_{A},\Theta_{B}$ coefficients. By having values of $q_{1}$ .. $q_{4}$
and consequently initial probability of occupancy eigenergies in case of
separated, but perturbatively interacting qubits $(E_{1A},E_{2A})$ and
$(E_{1B},E_{2B})$ we can determine quantum dot occupancy during next time
steps. Such procedure can be attempted for N interacting qubits under
condition of occupancy of 2 eigenenergy levels.

Figure 1: Case of position based qubit encoded in system of two coupled
quantum dot controlled electrostatically and case of simple Quantum Swap Gate
[6].

##  III Conclusions

It is quite straightforward to generalize the system under the consideration
of the case of N perturbatively electrostatically interacting position based
qubits described by Wannier functions. Those qubits can be further subjected
to Rabi oscillations by local microwave fields acting on qubits in their
direct proximity which might imply evolution of the system ground state with
time.

$\displaystyle\hat{H}_{A-B}=\hat{H}_{A}\times\hat{I}_{B}+\hat{I}_{A}\times\hat{H}_{B}+\hat{H}_{Coulomb:A-B}=$
$\displaystyle=(E_{p1A}\ket{w_{L,A}}\bra{w_{L,A}}+E_{p2A}\ket{w_{R,A}}\bra{w_{R,A}}+t_{s12A}\ket{w_{R,A}}\bra{w_{L,A}}+t_{s21A}\ket{w_{L,A}}\bra{w_{R,A}})\times\hat{I}_{B}+$
$\displaystyle+\hat{I}_{A}\times(E_{p1B}\ket{w_{L,B}}\bra{w_{L,B}}+E_{p2B}\ket{w_{R,B}}\bra{w_{R,B}}+t_{s12B}\ket{w_{R,B}}\bra{w_{L,B}}+t_{s21B}\ket{w_{L,B}}\bra{w_{R,B}})+$
$\displaystyle+\Bigg{[}q_{1}\ket{w_{L,A}}\ket{w_{L,B}}\bra{w_{L,A}}\bra{w_{L,B}}+q_{2}\ket{w_{L,A}}\ket{w_{L,B}}\bra{w_{R,A}}\bra{w_{R,B}}+q_{3}\ket{w_{L,A}}\ket{w_{L,B}}\bra{w_{R,A}}\bra{w_{R,B}}+$
$\displaystyle+q_{4}\ket{w_{R,A}}\ket{w_{R,B}}\bra{w_{R,A}}\bra{w_{R,B}}\Bigg{]}_{Coulomb}+$
$\displaystyle+\sum_{r,s,u,d}(\ket{w_{L},A}\bra{w_{L},A}+\ket{w_{R},A}\bra{w_{R},A})(\ket{w_{L},B}\bra{w_{L},b}+\ket{w_{R},B}\bra{w_{R},B})$
$\displaystyle
g(t,r,s,u,d)(\ket{E_{L,r}}\ket{E_{R,s}})(\bra{E_{L,u}}\bra{E_{R,d}})(\ket{w_{L},A}\bra{w_{L},A}+\ket{w_{R},A}\bra{w_{R},A})(\ket{w_{L},B}\bra{w_{L},b}+\ket{w_{R},B}\bra{w_{R},B}).$
(21)

It is quite straightforward to apply the presented scheme for the case of N
position based qubits implementing various quantum gates [7],[8]. Obtained
results have its significance in quantum communication as well [9].

## References

  * [1] A.R.Mills, D.M.Zajac, M.J.Gullans et al. , _Shuttling a single charge across a one-dimensional array of silicon quantum dots_ , Nature Communication, Vol. 10, 2019.
  * [2] T. Fujisawa, T. Hayashi, Y. Hirayama, H. D. Cheong, Y. H. Jeong, _Electron counting of single-electron tunneling current_ , Applied Physical Letters, Vol. 84, 2004
  * [3] K.Pomorski, P.Giounanlis, E.Blokhina, D.Leipold, R.B.Staszewski, _Analytic view on coupled single-electron lines_ , Semiconductor Science and Technology, Vol. 34, Nr 12, 2019
  * [4] K.Pomorski, _Electrostatically interacting Wannier qubits in curved space_ , Arxiv: 2112.13470v1, 2019.
  * [5] K.Pomorski, _Equivalence between finite state stochastic machine, non-dissipative and dissipative tight-binding and Schrödinger model_ , Mathematics and Computers in Simulation, Vol. 209, 2023, Elsevier.
  * [6] K.Pomorski, P.Peczkowski, R.B.Staszewski, _Analytical solutions for N interacting electron system confined in graph of coupled electrostatic semiconductor and superconducting quantum dots in tight-binding model_ , Cryogenics, Vol. 109, 2020.
  * [7] K.Pomorski, _Analytic view on N body interaction in electrostatic quantum gates and decoherence effects in tight-binding model_ , Vol. 19, Issue 04, 2021, International Journal of Quantum Information
  * [8] K.Pomorski, _Analytical solutions for N-Electron interacting system confined in graph of coupled electrostatic semiconductor and superconducting quantum dots in tight-binding model with focus on quantum information processing_ , Springer Proceedings in Physics, Vol. 279, 2023
  * [9] K.Pomorski, R.B.Staszewski, _Towards quantum internet and non-local communication in position-based qubits_ , AIP Conf. Proc., 2241, 2020