# Scheduling Quantum Annealing for Active User Detection in a NOMA Network

Romain Piron, Claire Goursaud

► To cite this version:

Romain Piron, Claire Goursaud. Scheduling Quantum Annealing for Active User Detection in a NOMA Network. Fifth IEEE International Conference on Quantum Computing and Engineering (QCE 2024), IEEE, Sep 2024, Montréal (Québec), Canada. hal-04664779

HAL Id: hal-04664779

https://hal.science/hal-04664779v1

Submitted on 30 Jul 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
# Scheduling Quantum Annealing for Active User Detection in a NOMA Network

Romain Piron

Université de Lyon, INSA Lyon, INRIA
CITI EA 3720

69621 Villeurbanne, France
romain.piron@insa-lyon.fr

Claire Goursaud

Université de Lyon, INSA Lyon, INRIA
CITI EA 3720

69621 Villeurbanne, France
claire.goursaud@insa-lyon.fr

***Abstract***—Active user detection in a non-orthogonal multiple access (NOMA) network is a major challenge for 5G/6G applications. However, classical algorithms that can perform this task suffer either from complexity or reduced performances. This work aims at proposing a quantum annealing approach to overcome this trade-off. Firstly, we show that the maximum a posteriori decoder of the activity pattern of the network can be seen as the ground state of an Ising Hamiltonian. For N users in a network with perfect channels, we propose a universal control function to schedule the annealing process. Our approach avoids to continuously compute the optimal control function but still ensures high success probability while demanding a lower annealing time than a linear control function. This advantage holds even in the presence of imperfections in the network.

***Index Terms***—Quantum annealing, local adiabaticity, control function, non-orthogonal multiple access, active user detection

## I. INTRODUCTION

One of the major challenges of 5G wireless networks is to achieve ultra-reliable low latency communications (uRLLC) requirements to serve applications as the Industrial Internet of Things [1]. However, the allocation of a dedicated time slot or frequency prior to transmission for each user is no longer adapted for typical uRLLC networks because of their growing size [2], [3]. Instead, each time that the users wish to transmit a signal, they must complete a random access procedure [4].

This procedure begins by an activity notification sent by the active users in the uplink channel. To do so, the code division multiple access (CDMA) approach can be used [5]. An identification sequence is attributed to each node of the network and all active users simply transmit their code. This approach permits sharing a common time-frequency resource to all the users without interference.

It recently appeared that leveraging non orthogonal multiple access (NOMA) in CDMA systems is promising [6]. In practice, it allows to attribute non orthogonal identification sequences to the nodes. Thus, one can reduce the size of these sequences for a given size of the users' set.

Detecting the active users — or the activity pattern — accurately from the received signal is called active user detection (AUD) [4]. It poses a significant challenge in NOMA systems since the cross-correlations between the sequences do not vanish. Thus, the users cannot be detected separately. While methods like the zero-forcing technique [7] or traditional cor-

relation receiver decoders [8] offer computational efficiency, their performance is limited. Conversely, the maximum a posteriori (MAP) Bayesian estimator [9] is the most reliable but computationally demanding, being NP-hard to compute.

In AUD scenarios in which the channel state is known, one can show that the computation of the MAP estimator is a quadratic unconstrained binary optimization (QUBO) [10] problem. The equivalence between QUBO problems and Ising Hamiltonians allows one to formulate the AUD problem as a ground state searching problem. However, an exhaustive search over the space of configurations to find the ground state has an exponential cost with the system size.

Several methods have been investigated to avoid this exhaustive search. Quantum annealing (QA) originally introduced by Kadowaki and Nishimori as a quantum analogue of thermal annealing in [11] is considered today as a promising quantum computing-based strategy for the ground state searching of several classes of Hamiltonians [10], [12]. Like many authors today, we use the term quantum annealing to refer to QA-inspired adiabatic quantum computation [10], [13]. Much work remains to be done regarding execution time, but QA already outperformed classical thermal annealing in reaching the ground state of some classes of Hamiltonians [10].

The control function that introduces transverse quantum fluctuations in QA usually needs to be adapted according to the parameters of the target Hamiltonian [14]. Unfortunately, in typical AUD models several variables are not deterministic but rather random variables [4]. This would require to compute the control function for each received signal. This is why, in this paper, we propose a strategy to build a generic control function to parameterize a QA process for AUD given the size *N* of the wireless network.

The rest of the paper is organized as follows. Sec. II defines our model and build the Ising Hamiltonian associated with the MAP estimator for AUD. In Sec. III, we briefly review the adiabatic theorem and show how to derive a control function for QA that ensures adiabaticity locally during the annealing process. Then, in Sec. IV we adopt a simplified scenario of our model to propose a universal control function to schedule a QA approach for AUD once the size of the network has been fixed. Finally, Sec. V takes into account the channel state information in our AUD model.
## Notations

The set of matrices over a field $\mathbb{K}$ of size $N \times M$ is denoted $\mathbb{K}^{N \times M}$ with the usual identification $\mathbb{K}^{N \times 1} \equiv \mathbb{K}^N$. We respectively denote a scalar, vector, matrix by $x$, $\mathbf{x}$, $\mathbf{X}$.

The diagonal matrix whose diagonal elements correspond to the components of a vector $\mathbf{x}$ is denoted $\text{diag}(\mathbf{x})$. We also introduce $\mathbf{1}_N = [1, \dots, 1] \in \mathbb{K}^N$, finally the identity matrix of size $M \times M$ corresponds to $I_M = \text{diag}(\mathbf{1}_M)$.

The multivariate normal distribution with the mean vector $\mu$ and the covariance matrix $\Sigma$ is denoted $\mathcal{N}(\mu, \Sigma)$. The Dirac distribution picked at $x$ is denoted $\text{Dirac}(x)$. The uniform distribution over a set $\Omega$ is denoted $\mathcal{U}(\Omega)$.

## II. QA FORMULATION FOR AUD

### A. Context and problem definition

We consider a network with $N$ users. Each user $\#i$ is equipped with a single antenna and assigned with a bipolar identification sequence [15] – also called code – of length $M$ denoted $c_i$. This code constitutes a signal that is transmitted by the $i$-th user. Its power is given by:

$$P_{\text{signal},i} = \|c_i\|^2. \quad (1)$$

It is often convenient to work with normalized signals such that $P_{\text{signal},i} = 1$. Hence, the bipolar codes are chosen in the set:

$$c_i \in \left\{ -\frac{1}{\sqrt{M}}, \frac{1}{\sqrt{M}} \right\}^M. \quad (2)$$

These vectors must also satisfy a quasi-linearly independence to ensure the decodability, which is discussed in Appendix A.

Let us assume that a random subset $\mathcal{N}_A$ of these users wish to transmit their code to a base station (BS) to notify their activity on the wireless channel. This subset constitutes the activity pattern of the network encoded in the bit-string $b^{(0)} \in \{0, 1\}^N$ so that one can write:

$$\mathcal{N}_A = \{i \in \{1, \dots, N\} \text{ s.t. } b_i^{(0)} = 1\}. \quad (3)$$

Each propagation path between the $i$-th user and the BS might suffer from attenuation due to the Rayleigh fading coming from the isotropic scattering assumption [16]. It is modelled by $N$ fading coefficients $(w_i)_{i=1,...,N}$ that obey a different probability law depending on the quality of the channels:

$$\mathbf{w} \sim \begin{cases} \text{Dirac}(\mathbf{1}_N) & \text{for perfect channels} \\ \mathcal{N}(\mathbf{0}, I_N) & \text{for imperfect channels} \end{cases} \quad (4)$$

where $\mathbf{w} = (w_i)_{i=1,...,N}$. At the end, one measures the following signal at the receiver side:

$$\mathbf{y} = \sum_{i=1}^{N} w_i b_i^{(0)} \mathbf{c}_i + \mathbf{n}, \quad (5)$$

where $\mathbf{n} \sim \mathcal{N}(\mathbf{0}, \xi^2 I_M)$ is an additive white Gaussian noise that models the thermal noise at the BS.

The AUD problem consists in retrieving the activity pattern $b^{(0)}$ of the network with the knowledge of the received

signal $\mathbf{y}$ and the fading coefficients $\mathbf{w}$. As mentioned in the introduction, the MAP estimator is the most reliable decoder for such task. It is defined by:

$$\hat{\mathbf{b}}^{\text{MAP}} = \arg \max_{\mathbf{b} \in \{0,1\}^N} \mathbb{P}_{\mathbf{b}|\mathbf{y},\mathbf{w}}(\mathbf{b}|\mathbf{y}, \mathbf{w}) \quad (6)$$

In the next part, we express the probability mass function $\mathbb{P}_{\mathbf{b}|\mathbf{y},\mathbf{w}}(\mathbf{b}|\mathbf{y}, \mathbf{w})$ in terms of the variables $b_i$.

### B. QUBO formulation for the MAP detector

Eq. 5 implies that $f_{\mathbf{y}|\mathbf{b},\mathbf{w}}$ simply corresponds to the probability density function of a Gaussian law. In addition, the Bayes' formula states that [17]:

$$\mathbb{P}_{\mathbf{b}|\mathbf{y},\mathbf{w}}(\mathbf{b}|\mathbf{y}, \mathbf{w}) f_{\mathbf{y}|\mathbf{w}}(\mathbf{y}|\mathbf{w}) = f_{\mathbf{y}|\mathbf{b},\mathbf{w}}(\mathbf{y}|\mathbf{b}, \mathbf{w}) \mathbb{P}_{\mathbf{b}}(\mathbf{b}). \quad (7)$$

From this point, a simple hypothesis consists in assuming a uniform distribution over the activity patterns. In this scenario, the prior distribution $\mathbb{P}_b$ is constant which yields:

$$\mathbb{P}_{\mathbf{b}|\mathbf{y},\mathbf{w}}(\mathbf{b}|\mathbf{y}, \mathbf{w}) \propto \exp \left( - \frac{\|\mathbf{y} - \mathbf{C} \cdot \text{diag}(\mathbf{w}) \cdot \mathbf{b}\|^2}{2\xi^2} \right), \quad (8)$$

where the sequences have been concatenated in the matrix $\mathbf{C} = [c_1, \dots, c_N] \in \mathbb{C}^{M \times N}$. Hence the MAP detector reads:

$$\hat{\mathbf{b}}^{\text{MAP}} = \arg \min_{\mathbf{b} \in \{0,1\}^N} \left\{ \|\mathbf{y} - \mathbf{C} \cdot \text{diag}(\mathbf{w}) \cdot \mathbf{b}\|^2 \right\}. \quad (9)$$

The computation of $\hat{\mathbf{b}}^{\text{MAP}}$ aims at recovering $b^{(0)}$ with $(\mathbf{y}, \mathbf{w})$ as inputs. Developing this objective function yields to the expression:

$$\|\mathbf{y} - \mathbf{C} \cdot \text{diag}(\mathbf{w}) \cdot \mathbf{b}\|^2 = \|\mathbf{y}\|^2 - 2 \sum_{i=1}^{N} (\mathbf{y} \cdot w_i \mathbf{c}_i) b_i + \sum_{i,j=1}^{N} (w_i \mathbf{c}_i \cdot w_j \mathbf{c}_j) b_i b_j. \quad (10)$$

Once one throws away the irrelevant constant $\|\mathbf{y}\|^2$, our minimization problem corresponds to the general form

$$\hat{\mathbf{b}}^{\text{MAP}} = \arg \min_{\mathbf{b} \in \{0,1\}^N} \left\{ \sum_{i,j} Q_{ij} b_i b_j + \sum_i a_i b_i \right\}. \quad (11)$$

Thus, the MAP detector can be written as a QUBO problem [10] expressed in terms of the binary variables $b_i$.

### C. Building the Ising Hamiltonian

As underlined by the authors of [10], a QUBO problem can be mapped to an Ising Hamiltonian by the change of variable $\sigma_i = 1 - 2b_i$. We defer the details to appendix B and simply give here the final expression of the problem Hamiltonian:

$$H_P = - \sum_i h_i \sigma_i - \sum_{i<j} J_{ij} \sigma_i \sigma_j, \quad (12)$$
expressed with classical Ising spins $\sigma_i \in \{+1 \equiv \uparrow, -1 \equiv \downarrow\}$ and the parameters:

$$ \left\{ \begin{aligned} h_i &= -w_i (\mathbf{c}_i \cdot \tilde{\mathbf{y}}) \\ J_{ij} &= -\frac{1}{2}(w_i w_j) (\mathbf{c}_i \cdot \mathbf{c}_j), \end{aligned} \right. \qquad (13) $$

where the signal $\tilde{\mathbf{y}}$ is defined by:

$$ \tilde{\mathbf{y}} = \mathbf{y} - \frac{1}{2} \sum_{j=1}^{N} w_j \mathbf{c}_j. \qquad (14) $$

In terms of Ising variables, an inactive (resp. active) user corresponds to $\sigma_i = \uparrow$ (resp. $\sigma_i = \downarrow$). The MAP estimator is the ground state of this Hamiltonian. It is the spin configuration $\boldsymbol{\sigma}^{(0)} \in \{\uparrow, \downarrow\}^N$ associated to $\mathbf{b}^{(0)}$ that minimizes the energy of the system. The decodability constraint previously mentioned about the codes $\mathbf{c}_i$ ensures that the ground state is unique in the idealized situation $\mathbf{w} = \mathbf{1}_N$ and $\xi = 0$.

Let us quickly mention a simple sanity check that confirms that this Hamiltonian is well suited to our problem. We adopt the scenario $\mathbf{w} = \mathbf{1}_N$, $\xi = 0$ and relax the NOMA constraint by considering orthogonal identifications sequences, ie. $\mathbf{c}_i \cdot \mathbf{c}_j = 0$ for $i \neq j$. Then, the activity pattern is trivially recovered by:

$$ b_i^{(0)} = |\mathbf{y} \cdot \mathbf{c}_i|. \qquad (15) $$

On the other hand, Eq. 13 in this ideal situation reduces to:

$$ \left\{ \begin{aligned} h_i &= ||\mathbf{c}_i||^2 \sigma_i^{(0)} , \\ J_{ij} &= 0 \end{aligned} \right. , \qquad (16) $$

where $\boldsymbol{\sigma}^{(0)}$ is the spin configuration associated to the initial activity pattern $\mathbf{b}^{(0)}$. The ground state of $H_P$ is trivially $\boldsymbol{\sigma}^{(0)}$, as one might reasonably expect.

The ground state of an Ising Hamiltonian can be found with perfect accuracy through an exhaustive search in the space of configurations $\{\uparrow, \downarrow\}^N$. However, such task would require $O(2^N)$ iterations which is why it is known to be NP-hard. In order to perform ground state searching with QA, we promote the Ising spins $\sigma_i$ to quantum spins $\hat{\sigma}_i^z$. The problem Hamiltonian is thus encoded in the operator:

$$ \hat{H}_P = - \sum_i h_i \hat{\sigma}_i^z - \sum_{i<j} J_{ij} \hat{\sigma}_i^z \hat{\sigma}_j^z. \qquad (17) $$

### III. STATE OF THE ART FOR QA SCHEDULING

#### A. Basic principle of QA

Given a quantum Hamiltonian $\hat{H}_P$ whose ground state is hard to find, one uses a control Hamiltonian $\hat{H}_C$ with a simple ground state to initialize the system in this known configuration. In the initial formulation of QA by Kadowaki and Nishimori [11], $\hat{H}_C$ is chosen as the uniform transverse field Hamiltonian:

$$ \hat{H}_C = - \sum_{i=1}^{N} \hat{\sigma}_i^x, \qquad (18) $$

whose ground state is the uniform superposition $|+\rangle^{\otimes N}$. The QA process slowly evolves the system from an initial dynamic dictated by $\hat{H}_C$ towards a final dynamic encoded by $\hat{H}_P$ through the global Hamiltonian:

$$ \hat{H}(t) = (1 - u(t)) \hat{H}_P + u(t) \hat{H}_C. \qquad (19) $$

The function $u(t)$ is a control function obeying the initial condition $u(0) = 1$ so that the system starts in the ground state $|+\rangle^{\otimes N}$. If the control function decreases slowly towards 0 during a period $T$, the final state of the system is expected to be the ground state of $\hat{H}_P$ with high probability. This period $T$ is called the annealing time and its duration can be estimated thanks to the adiabatic theorem.

#### B. Adiabatic approximation

The so-called adiabatic theorem is a mathematical statement on the time evolution of the global Hamiltonian $\hat{H}(t)$ to ensure that the instantaneous state $|\psi(t)\rangle$ of the system remains close to the instantaneous ground state $|\phi_0(t)\rangle$. Several formulations exist with different precision levels reviewed in [18], [19]. If one denotes $\{\varepsilon_n(t), |\phi_n(t)\rangle\}$ the instantaneous eigensystem of $\hat{H}(t)$ and $\Delta(t) = \varepsilon_1(t) - \varepsilon_0(t)$ is the spectral gap between the ground state and the first excited state, the most commonly used adiabatic approximation reads [10]:

$$ \frac{|\langle\phi_1(t)|\partial_t \hat{H}(t)|\phi_0(t)\rangle|}{\Delta(t)^2} \ll 1 \quad \forall t \in [0, T] \Rightarrow P_{QA}(T) \approx 1, \qquad (20) $$

where we introduced the notation used by [11]:

$$ P_{QA}(t) = |\langle\phi_0(T)|\psi(t)\rangle|^2. \qquad (21) $$

The success probability of QA naturally corresponds to $P_{QA}(T)$. In our case, the time dependence of $\hat{H}(t)$ is fully contained in the control function $u$. Using Eq. 19, the adiabatic approximation reads:

$$ \left| \frac{du}{dt} \right| \left| \frac{\langle \hat{H}_C - \hat{H}_P \rangle_{1,0}(u)}{\Delta(u)^2} \right| \ll 1 \quad \forall t \in [0, T] \Rightarrow P_{QA}(T) \approx 1. \qquad (22) $$

which constraints the control function $u$. It is quite easy to show that the numerator $|\langle \hat{H}_C - \hat{H}_P \rangle_{1,0}(u)|$ scales at most as $O(\text{Poly}(N))$ (see [19] for instance). Thus, a qualitative discussion of the adiabatic approximation as a function of the system size is often focused on the spectral gap $\Delta^2$. At the end, one usually uses the following constraint on the control function $u$:

$$ \left| \frac{du}{dt} \right| \Delta(u)^{-2} \ll 1 \quad \forall t \in [0, T]. \qquad (23) $$

#### C. Linear scheduling

Several theoretical analysis of QA use a linear control function [20] to parameterize the annealing process:

$$ u_{\text{lin}}(t) = 1 - \frac{t}{T_{\text{lin}}}. \qquad (24) $$

<footnote n="1">¹Recall that $\hat{\sigma}_i^\alpha$ denotes the Pauli matrices acting on the i-th qubit, where $\alpha = x, y, z$</footnote>
<figure><img src="image_2.png" /><figcaption>Figure 1: Two first eigenvalues of $\hat{H}(u)$ for some values of the control function $u \in [0, 1]$ (a) and corresponding evolution of the spectral gap (b). It allows to solve numerically Eq. 28 to obtain $u_{(y,w)}$ (c), which yields the following evolution for the overlap $P_{QA}(t)$ (d)</figcaption></figure>

In this case, $T_{\text{lin}}$ naturally corresponds to the annealing time. However, the first derivative of $u$ is constant which heavily constraints $T_{\text{lin}}$. Indeed, Eq. 23 immediately imposes:

$$T_{\text{lin}} \gg \Delta_{\min}^{-2}, \qquad (25)$$

where $\Delta_{\min}^2 = \min_{0 \le u \le 1} \Delta^2(u)$. However, if the system encounters a first order quantum phase transition during the annealing process, it is known that the minimum of the gap closes exponentially with the system size [10]. Concretely, it means that $T_{\text{lin}} \sim e^{\alpha N}$ with $\alpha > 1$ so QA would offer no advantage over an exhaustive search to compute $\hat{b}^{\text{MAP}}$.

The work [14] underlines that in this case one can enhance the annealing time by looking for a control function with a non constant derivative.

### D. Optimized control function

Recall that $u$ is a decreasing function from $u(0) = 1$ to $u(T) = 0$, which means that $du/dt < 0$. As it is done in [14], we can fix an arbitrary small constant $\epsilon > 0$ and impose:

$$\left\{ \begin{array}{l} \frac{du}{dt} = -\epsilon \Delta^2(u(t)), \\ u(0) = 1 \end{array} \right. \qquad (26)$$

which directly comes from Eq. 23. We refer to the solution of ordinary differential equation (ODE) as the optimized control function. Qualitatively, this equation tells that the larger the gap, the faster the control function can decrease. Integrating both sides of this ODE from $t=0$ to $t=T$ yields the same expression than [21] for the annealing time:

$$T = \frac{1}{\epsilon} \int_{0}^{1} \frac{du}{\Delta^{2}(u)}. \qquad (27)$$

These two equations 26, 27 underline that once the gap is known, one can fully determine the annealing time and the control function associated to a problem encoded in $\hat{H}_P$.

### E. An example

Going back to AUD, one can build the quantum Ising Hamiltonian 17 once the received signal $y$ and the channel coefficients $w$ are known. Thus, we call a pair $(y, w)$ a problem instance that allows to compute the Ising parameters

of Eq. 13. Since the time dependence of the full Hamiltonian $\hat{H}$ is entirely contained in the control function, one can take $u \in [0, 1]$ as an affine parameter to evaluate the spectral gap as a function of $u$ before even knowing the time dependence $u(t)$. Thus, we denote $\Delta_{(y,w)}^2(u)$ the spectral gap associated to the problem instance $(y, w)$. From Eq. 26 and 27, one defines the associated control function:

$$\left\{ \begin{array}{l} \frac{du_{(y,w)}}{dt} = -\epsilon \Delta_{(y,w)}^2(u_{(y,w)}(t)), \\ u_{(y,w)}(0) = 1 \end{array} \right. \qquad (28)$$

which is the optimal one for $(y, w)$. The associated annealing time reads:

$$T_{(y,w)} = \frac{1}{\epsilon} \int_{0}^{1} \frac{du}{\Delta_{(y,w)}^{2}(u)}. \qquad (29)$$

As an example, we considered a network of $N=8$ users and generated a problem instance $(y, w)$ in the following scenario:

$$\left\{ \begin{array}{l} \text{Activity pattern:} \quad \mathbf{b}^{(0)} = (1, 0, 0, 0, 0, 0, 0, 0) \\ \text{Perfect channels:} \quad \mathbf{w} \sim \text{Dirac}(\mathbf{1}_N) \\ \text{No additive noise:} \quad \xi = 0 \end{array} \right. .$$

We show on Fig. 1a the two first energy levels of the global Hamiltonian $\hat{H}(u)$. Fig. 1b gives the corresponding spectral gap $\Delta_{(y,w)}^2$. These results have been obtained by exact diagonalization of the Hamiltonian using the QuTiP<fnref n="2" /> library. Then we numerically evaluate the annealing time $T_{(y,w)}$ and solve Eq. 28 on the interval $[0, T_{(y,w)}]$ to obtain the optimized control function (Fig. 1c). To do so, we fixed our precision parameter to $\epsilon = 0.1$. The evolution of the overlap during the annealing process scheduled with this control function is finally shown on Fig. 1d. One can see that the success probability $P_{QA}(T_{(y,w)})$ reached at the end of the evolution is very close to 1, as expected. However, it would be very long and inconvenient to require to compute $T_{(y,w)}$ and $u_{(y,w)}$ for each instance $(y, w)$ of our problem. The goal of the following sections is to emphasize the possibility to build a generic approach to parameterize a QA process for AUD.

<footnote n="2">²https://qutip.org</footnote>
<figure><img src="image_3.png" /><figcaption>Figure 2: Several gaps $\Delta^2_{(y,w)}$ associated to random samples of $(y, w = 1_N)$ and their mean value (in red) for $N = 8$.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Figure 3: Shape of the mean control functions obtained for $N = 6, \dots, 9$ with respect to $s = t/T_{(\text{mean}),N}$</figcaption></figure>

## IV. PROPOSITION: GENERIC CONTROL FUNCTION WITH PERFECT CHANNELS

We would like to check if, despite the randomness of the parameters $(y, w)$, we can estimate a generic control function for a given size of the network $N$. We first consider a simple scenario with perfect channels.

### A. Mean gap in the absence of additive noise

Let us begin with a scenario where the size of the network $N$ is fixed and where the additive Gaussian noise is neglected, which amounts to fix $\xi = 0$. We randomly generated 30 samples of $(y, w)$ to check how the associated gaps $\Delta^2_{(y,w)}$ behave. The signals $y$ are generated from uniform samples of activity patterns $b^{(0)}$ and the channel coefficients are given by $w = 1_N$. Interestingly, Fig. 2 shows that a pattern emerges for the shape of the gap against $u$ despite the dependence of the local magnetic fields $h_i$ in $y$.

These results motivate the definition of the mean gap:

$$ \Delta^2_{(\text{mean}),N}(u) = \mathbb{E}_{(y,w)}(\Delta^2). \quad (30) $$

In our specific case, this quantity can be exactly computed. Indeed, the channel coefficients are fixed to $w = 1_N$ and the received signal $y$ is fully determined from the initial activity pattern $b^{(0)}$. Using that $b^{(0)}$ is uniformly distributed, one has:

$$ \Delta^2_{(\text{mean}),N} = \frac{1}{2^N} \sum_{b^{(0)} \in \{0,1\}^N} \Delta^2_{(y,w=1_N)}(u). \quad (31) $$

Then, we can use this mean gap to compute the associated annealing time $T_{(\text{mean}),N}$ with Eq. 27:

$$ T_{(\text{mean}),N} = \frac{1}{\epsilon} \int_{0}^{1} \frac{du}{\Delta^2_{(\text{mean}),N}(u)}, \quad (32) $$

and the resulting control function $u_{(\text{mean}),N}(t)$ from Eq. 26:

$$ \left\{ \begin{array}{l} \frac{du_{(\text{mean}),N}}{dt} = -\epsilon \Delta^2_{(\text{mean}),N} (u_{(\text{mean}),N}(t)). \quad (33) \\ u_{(\text{mean}),N}(0) = 1 \end{array} \right. $$

Once for all, we adopt the precision level $\epsilon = 0.1$. We reported the shape obtained for $u_{(\text{mean}),N}(s)$ for $N = 6, \dots, 9$ on Fig. 3 with respect to the normalized time $s = t/T_{(\text{mean}),N}$.

Then, we must test for each value of $N$ whether these functions are good candidates to schedule a QA process on the period $t \in [0, T_{(\text{mean}),N}]$ to estimate $\hat{b}^{\text{MAP}}$.

### B. Definition of the metrics

We will compare for several $(y, w)$ the success probability of a QA process scheduled with our approach with the one that would been obtained with an optimal scheduling for each problem instance. To do so, one need to define appropriate metrics. We first generate a pair $(y, w)$ and simulate a QA process parameterized with $u_{(\text{mean}),N}$ by numerically solving the Schrödinger's equation:

$$ i \partial_t |\psi(t)\rangle = \hat{H}(u_{(\text{mean}),N}(t)) |\psi(t)\rangle, \quad t \in [0, T_{(\text{mean}),N}], \quad (34) $$

with the initial condition $|\psi(0)\rangle = |+\rangle^{\otimes N}$ (we work with natural units where $\hbar = 1$) using again the QuTiP library. The obtained success probability is denoted with the short-hand notation:

$$ P_{\text{QA}}(T_{(\text{mean}),N}; (y, w)) \equiv P_{\text{QA}}^{(\text{mean})}(y, w). \quad (35) $$

On the other hand, the reference success probability is the one corresponding to a QA process scheduled with the well-suited optimized control function $u_{(y,w)}$ (Eq. 26) and the corresponding annealing time $T_{(y,w)}$ (Eq. 27). We adopt also a short hand notation:

$$ P_{\text{QA}}(T_{(y,w)}; (y, w)) \equiv P_{\text{QA}}^{(\text{opti})}(y, w). \quad (36) $$

Since the instances $(y, w)$ are randomly distributed, $P_{\text{QA}}^{(\text{mean})}$ and $P_{\text{QA}}^{(\text{opti})}$ are random variables. Thus, we propose to use their expectation values to quantify the reliability of each scheduling strategy. We turn on the additive Gaussian noise to realize the performance analysis, hence the expectation values are given by:

$$ \mathbb{E}_{(y,w)} (P_{\text{QA}}^{(\alpha)}) = \sum_{b^{(0)} \in \{0,1\}^N} \int dn f(b^{(0)}, n) P_{\text{QA}}^{(\alpha)}(y, w), \quad (37) $$

for $\alpha \in \{\text{mean, opti}\}$. The function $f(b^{(0)}, n)$ is the probability density function associated to the instances $(y, w)$.
<figure><img src="image_5.png" /><figcaption>Figure 4: Distributions of $P_{QA}^{(\text{mean})}$ (in blue) and $P_{QA}^{(\text{opti})}$ (in red) for several network sizes $N$. The noise level corresponds to SNR = 20dB We performed $N_{\text{samples}} = 500$ samples for each histogram.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Figure 5: Same distributions than Fig. 4 but with SNR = 15dB.</figcaption></figure>

The channel coefficients are still fixed to $\mathbf{w} = \mathbf{1}_N$ but the signal $\mathbf{y}$ is now generated from the sample of an activity pattern $\mathbf{b}^{(0)}$ and a Gaussian vector $\mathbf{n}$. Those two parameters are independently distributed, respectively following the laws $U(\{0, 1\}^N)$ and $N(0, \xi^2)$. Thus, the weights of Eq. 37 are given by:

$$f(\mathbf{b}^{(0)}, \mathbf{n}) = \frac{1}{2^N} \times \frac{1}{(2\pi\xi)^{M/2}} \exp\left(-\frac{\|\mathbf{n}\|^2}{2\xi^2}\right). \quad (38)$$

We also need a metric to quantify the intensity of the Gaussian noise. A well adapted one is the signal-to-noise ratio (SNR) expressed in dB. Since the power of the signal sent individually by each user is normalized to $P_{\text{signal},i} = 1$, the SNR roughly reads:

$$\text{SNR} = -20 \log_{10}(\xi). \quad (39)$$

### C. Reliability of our strategy

The instances $(\mathbf{y}, \mathbf{w})$ have no longer a finite support. Thus, we generated $N_{\text{samples}} = 500$ pairs $(\mathbf{y}, \mathbf{w} = \mathbf{1}_N)$ according to the weights of Eq.38 for $N = 6, \dots, 9$ to estimate the distributions of $P_{QA}^{(\text{mean})}$ and $P_{QA}^{(\text{opti})}$. The expectation values for $\alpha \in \{\text{mean}, \text{opti}\}$ are estimated by:

$$\mathbb{E}_{(\mathbf{y}, \mathbf{w})} (P_{QA}^{(\alpha)}) = \frac{1}{N_{\text{samples}}} \sum_{(\mathbf{y}, \mathbf{w}) \in \text{samples}} P_{QA}^{(\alpha)}(\mathbf{y}, \mathbf{w}). \quad (40)$$

We reported on Fig. 4 the obtained histograms for a noise level corresponding to SNR = 20dB. We can observe that the values of $P_{QA}^{(\text{mean})}$ are more spread out along the axis

than the values of $P_{QA}^{(\text{opti})}$ which is not surprising. Indeed, the control function $u_{\text{mean},N}$ cannot be well suited for all possible outcomes $(\mathbf{y}, \mathbf{w})$. Nevertheless, the expectation value $\mathbb{E}_{(\mathbf{y}, \mathbf{w})} (P_{QA}^{(\text{mean})})$ is always above 0.8 which is reasonably below $\mathbb{E}_{(\mathbf{y}, \mathbf{w})} (P_{QA}^{(\text{opti})})$.

We also reported the obtained distributions for a level of noise corresponding to SNR = 15dB. The gap between $\mathbb{E}_{(\mathbf{y}, \mathbf{w})} (P_{QA}^{(\text{mean})})$ and $\mathbb{E}_{(\mathbf{y}, \mathbf{w})} (P_{QA}^{(\text{opti})})$ is higher but our approach can still offer a mean success probability above 0.7.

Thus, our generic control function appears as a good compromise to schedule a QA approach for any set of parameters $(\mathbf{y}, \mathbf{w})$ given a size of the network $N$.

### D. Advantage over a linear control function

Once the annealing time $T_{(\text{mean}),N}$ is known thanks to Eq. 32, one might question the necessity to use the control function $u_{(\text{mean}),N}$ on the annealing period $[0, T_{(\text{mean}),N}]$. Indeed, the differential equation 26 to obtain the optimized control function associated to a given gap results from several approximations. Thus, one might be tempted to simply use a linear control function with an annealing period $T_{\text{lin},N} = T_{(\text{mean}),N}$ and expect the resulting success probability to be not too low even if the condition of Eq. 25 is not satisfied. To answer this question, we define the success probability of a QA process planned according to this strategy:

$$P_{QA} (T_{\text{lin}}, N = T_{(\text{mean}),N}; (\mathbf{y}, \mathbf{w})) \equiv P_{QA}^{(\text{lin})} (\mathbf{y}, \mathbf{w}). \quad (41)$$
<figure><img src="image_7.png" /><figcaption>Figure 6: Distributions of $P_{QA}^{(\text{mean})}$ (in blue) and $P_{QA}^{(\text{lin})}$ (in green) for several network sizes $N$. The noise level corresponds to SNR = 20dB We performed $N_{\text{samples}} = 500$ samples for each histogram.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Figure 7: Same distributions than Fig. 6 but with SNR = 15dB.</figcaption></figure>

Following the same approach than above, we reported on Fig. 6 the distributions of $P_{QA}^{(\text{lin})}$ and $P_{QA}^{(\text{mean})}$ for a level of noise corresponding to SNR = 20dB. The expectation value $\mathbb{E}_{(y,w)}(P_{QA}^{(\text{mean})})$ is higher than $\mathbb{E}_{(y,w)}(P_{QA}^{(\text{lin})})$ which confirms the advantage of using $u_{\text{mean},N}$.

As the noise increases, Fig. 7 shows that the gap between the two expectation values $\mathbb{E}_{(y,w)}(P_{QA}^{(\text{mean})})$ and $\mathbb{E}_{(y,w)}(P_{QA}^{(\text{lin})})$ is reduced. Nevertheless at SNR = 15dB, our control function is still expected to yield a success probability higher than the one that would been obtained with a linear control function

V. CONTROL FUNCTION IN THE GENERAL CASE

In the previous section we proposed a generic method to evaluate an appropriate control function and demonstrated its validity on random problem instances with thermal noise at the BS. However in more realistic scenarios, the received signal can also be impacted by the channel coefficients. In this section, we take into account the channels between the users and the BS. From now on, the channel coefficients are normal variables to account for Rayleigh fading:

$$w \sim \mathcal{N}(0, I_N). \qquad (42)$$

A. Behavior of the gap

First of all, let us go back to our scenario of $N = 8$ users with an initial activity pattern $b^{(0)} = (1, 0, 0, 0, 0, 0, 0, 0)$. We generate one random problem instance $(y, w)$ with no additional Gaussian noise by independently sampling the $N$

<figure><img src="image_9.png" /><figcaption>Figure 8: Two first eigenvalues of $\hat{H}(u)$ in a scenario with attenuation (a) and corresponding spectral gap against $u$ (b).</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>Figure 9: Shape of the mean control functions with fading coefficients in the network obtained for $N = 6, \dots, 9$ with respect to $s = t/T_{(\text{mean}),N}$</figcaption></figure>
<figure><img src="image_11.png" /><figcaption>Figure 10: Distributions of $P_{QA}^{(\text{mean})}$ (in blue) and $P_{QA}^{(\text{opti})}$ (in red) for several network sizes $N$ with channel imperfections taken into account. The noise level corresponds to SNR = 15dB. We performed $N_{\text{samples}} = 500$ samples for each histogram.</figcaption></figure>

channel coefficients according to the standard normal distribution.

As done before, we reported on Fig. 8 the eigenvalues $\epsilon_{0,1}$ of the global Hamiltonian $\hat{H}$ and the corresponding spectral gap. This time, the minimum of the gap is reached at the end of the annealing process when the control function vanishes.

This new behavior can be explained by the modifications induced in the spectrum of the Ising Hamiltonian $\hat{H}_P$. In the case where a user $i$ is strongly attenuated by a coefficient $w_i$ close from 0, the spin flip $\sigma_i \to -\sigma_i$ will not significantly modify the associated eigenenergy. Thus, taking into account the attenuation on the propagation paths between the users and the access point strongly reduces the gaps between the eigenvalues of the associated Ising Hamiltonian $\hat{H}_P$. It suggests that the mean gap should be estimated again in this new scenario.

### B. Mean control function with fading coefficients

As done previously, we first neglect the additive Gaussian noise in order to evaluate the mean gap from several problem instances. From this quantity, we derive the associated mean control function $u_{(\text{mean}),N}$ defined on $[0, T_{(\text{mean}),N}]$. Since the channel coefficients $w$ are now also randomly distributed, Eq. 31 is generalized to:

$$ \Delta_{(\text{mean}),N}^2(u) = \sum_{\mathbf{b}^{(0)} \in \{0,1\}^N} \int d\mathbf{w} f(\mathbf{b}^{(0)}, \mathbf{w}) \Delta_{(y,\mathbf{w})}^2(u). \quad (43) $$

Since $\mathbf{b}^{(0)}$ and $\mathbf{w}$ are independently distributed, the above probability density function is simply given by:

$$ f(\mathbf{b}^{(0)}, \mathbf{w}) = \frac{1}{2^N} \times \frac{1}{(2\pi)^{N/2}} \exp\left(-\frac{\|\mathbf{w}\|^2}{2}\right). \quad (44) $$

Contrarily to the scenario with perfect channels, the mean gap cannot be evaluated exactly through Eq. 43. Thus, we sample the probability distribution $f(\mathbf{b}^{(0)}, \mathbf{w})$ with $N_{\text{samples, gap}} = 2000$ samples and evaluate the mean gap with the estimator:

$$ \Delta_{(\text{mean}),N}^2(u) = \frac{1}{N_{\text{samples, gap}}} \sum_{(\mathbf{y},\mathbf{w}) \in \text{samples}} \Delta_{(\mathbf{y},\mathbf{w})}^2(u). \quad (45) $$

As done previously, we reported on Fig. 9 the shape of the control functions against the normalized time $s = t/T_{(\text{mean}),N}$ for $N = 6, \dots, 9$. These shapes indicate that the first derivative is monotonically decreasing over the annealing period, which is expected regarding the previous comments about the behavior of the gaps $\Delta_{(y,w)}^2$. Let us now check whether this approach still offers a good compromise in terms of success probability.

### C. Mean success probability

As we did with the perfect channel scenario, we would like to check that the mean control function still offers a good success probability when the noise is taken into account. We use again the metrics $P_{QA}^{(\text{mean})}$ and $P_{QA}^{(\text{opti})}$. The expression of their expectation values is modified according to the distribution of the instances $(\mathbf{y}, \mathbf{w})$:

$$ \mathbb{E}_{(\mathbf{y}, \mathbf{w})} (P_{\text{QA}}^{(\alpha)}) = \sum_{\mathbf{b}^{(0)} \in \{0,1\}^N} \int d\mathbf{n} d\mathbf{w} f(\mathbf{b}^{(0)}, \mathbf{n}, \mathbf{w}) P_{\text{QA}}^{(\alpha)}(\mathbf{y}, \mathbf{w}). \quad (46) $$

for $\alpha \in \{\text{mean}, \text{opti}\}$. Using again the independence of the parameters used to generate a pair $(\mathbf{y}, \mathbf{w})$, the above probability distribution is given by:

$$ f(\mathbf{b}^{(0)}, \mathbf{n}, \mathbf{w}) = \frac{1}{2^N} \times \frac{1}{(2\pi\xi)^{M/2}} \exp\left(-\frac{\|\mathbf{n}\|^2}{2\xi^2}\right) \times \frac{1}{(2\pi)^{N/2}} \exp\left(-\frac{\|\mathbf{w}\|^2}{2}\right). \quad (47) $$

An estimation of the distributions of $P_{QA}^{(\text{mean})}$ and $P_{QA}^{(\text{opti})}$ at SNR = 15dB is shown on the histograms of Fig. 10. These results have been obtained by performing again $N_{\text{samples}} = 500$ samples of problem instances to evaluate $P_{QA}^{(\text{mean})}(\mathbf{y}, \mathbf{w})$ and $P_{QA}^{(\text{opti})}(\mathbf{y}, \mathbf{w})$. Unfortunately, the mean success probability $\mathbb{E}_{(\mathbf{y}, \mathbf{w})} (P_{QA}^{(\text{mean})})$ is around 0.65 while a well-suited control function for each problem instance still ensures a mean success probability $\mathbb{E}_{(\mathbf{y}, \mathbf{w})} (P_{QA}^{(\text{opti})})$ above 0.9. We must then investigate a way to improve the mean control function in this more complex AUD scenario.
<figure><img src="image_12.png" /><figcaption>Figure 11: Distributions $P_{QA}^{(\text{mean})}$ (in blue) and $P_{QA}^{(\text{opti})}$ (in red) for $N = 8$ and SNR = 15dB with different dilatation factors $\lambda$. The annealing period $[0, T_{(\text{mean})}]$ is extended to $[0, \lambda T_{(\text{mean})}]$ and $[0, T_{(y,w)}]$ is left unchanged.</figcaption></figure>

<figure><img src="image_13.png" /><figcaption>Figure 12: Distributions $P_{QA}^{(\text{mean})}$ (in blue) and $P_{QA}^{(\text{lin})}$ (in green) for $N = 8$ and SNR = 15dB with different dilatation factors $\lambda$.</figcaption></figure>

### D. Time dilatation

The reduction of the mean success probability $E_{(y,w)}(P_{QA}^{(\text{mean})})$ is likely due to the fact that $u_{(\text{mean}),N}$ decreases too fast at the end of the annealing process. Indeed, the example previously shown on Fig. 8 suggests that when the control function approaches 0, the gap closes which means that $u$ must decrease slowly enough to offer a good reliability. In order to force $u_{(\text{mean}),N}$ to decrease slower, we propose to apply a time dilatation factor $\lambda$ to its associated annealing period $[0, T_{(\text{mean}),N}]$. In practice, it amounts to use the time-dilated mean control function:

$$u_{(\text{mean}),\lambda,N}(t) = u_{(\text{mean}),N}\left(\frac{t}{\lambda}\right). \quad (48)$$

Consequently, the annealing period is extended to $[0, \lambda T_{(\text{mean}),N}]$. Of course, we aim at keeping $\lambda$ reasonably low to preserve the order of magnitude of the annealing time.

Nevertheless, we show on Fig. 11 that time dilatation factors of $\lambda = 2, 3, 4$ already allow to boost the mean success probability $E_{(y,w)}(P_{QA}^{(\text{mean})})$. These results correspond to a network with $N = 8$ users and a level of additive noise fixed to SNR = 15dB. For instance, picking $\lambda = 2$ allows to reach a mean success probability above 0.75 which is significantly better than the results obtained in the previous paragraph.

We also compared again our proposal to a linear control function parameterized over the same annealing period $[0, \lambda T_{(\text{mean}),N}]$. The results on Fig. 12 confirm that

$u_{(\text{mean}),\lambda,N}$ offers significantly a higher mean success probability $E_{(y,w)}(P_{QA}^{(\text{mean})})$ than a linear control function.

Thus, our strategy is also promising in this more complex AUD scenario. A time dilatation can be applied to the mean control function to improve the performances.

## VI. CONCLUSION

### A. Contributions

In this work, we proposed a generic strategy to parame-terize a QA process with the usual transverse Hamiltonian to compute the MAP estimator for AUD in a NOMA network. Despite the randomness of several parameters involved in the couplings and the local fields of the corresponding targeted Ising Hamiltonian, the evolution of the spectral gap of the global Hamiltonian against $u(t)$ follows a recurrent shape. We have shown that our mean control function ensures a reasonable expectation value for the success probability of QA after (if wanted) applying a time dilatation operation. This function has the advantage to be universal once the size of the network $N$ is fixed and permits to reach better performances than the simple linear control function.

### B. Scalability with $N$

The work conducted in this paper is limited to small network sizes. For larger values of $N$, the computation of the quantum gaps like those involved in Eq. 45 would be too long with the exact diagonalization method of the full Hamiltonian used in
this work. A further improvement would be to use Quantum Monte Carlo methods as done by [22] to extract quantum gaps involved in the mean gap estimator in order to compute it for $N \ge 9$.

### C. Knowledge of the channel coefficients

Regarding the AUD model, another improvement would be to relax some assumptions we made. In fact, the knowledge of the channel coefficients $w$ at the access point is a strong assumption that does not hold true in practice. Some modern formulations of AUD assume that the statistics of the channel vector $w$ (which is a Gaussian law in this work) is known by the access point but not the realizations. A further perspective to this work could be to adapt the signal processing model to such case.

## ACKNOWLEDGMENT

The authors would like to thank Lélio Chetot and Bruno Fedrici for helpful discussions.

This work was supported by the ANR under the France 2030 program, grant "NF-PERSEUS : ANR-22-PEFT-0004".

## APPENDIX A
### CHOICE OF THE SEQUENCES

The identification sequences of the users must be chosen such that each activity pattern $b^{(0)}$ yields a different signal $y$ in the absence of imperfections. It means that in the case $\xi = 0$ and $w = 1_N$, one requires:

$$ \left\{ \begin{aligned} \mathbf{y}(\mathbf{b}^{(0)}, \mathbf{n}=\mathbf{0}, \mathbf{w}=\mathbf{1}_N) &= \mathbf{y}(\mathbf{b}^{(1)}, \mathbf{n}=\mathbf{0}, \mathbf{w}=\mathbf{1}_N) \\ \Rightarrow \mathbf{b}^{(0)} &= \mathbf{b}^{(1)} \end{aligned} \right\}. \quad (49) $$

One can show that this condition is ensured if and only if:

$$ \left\{ \exists \boldsymbol{\lambda} \in \{-1, 0, 1\}^N \text{ s.t. } \sum_{i=1}^{N} \lambda_i \mathbf{c}_i = \mathbf{0} \right\} \Rightarrow \{\boldsymbol{\lambda} = \mathbf{0}\}. \quad (50) $$

Throughout this work, we used coding schemes constructed by a random selection of $N$ codes of size $M$ followed by an exhaustive verification of the condition 50. We give below the matrices $C = [c_1, \dots, c_N] \in \mathbb{C}^{M \times N}$ corresponding to the codes we used for $N = 6, \dots, 9$.

$$ C_6 = \frac{1}{\sqrt{5}} \begin{pmatrix} -1 & -1 & 1 & 1 & -1 & -1 \\ 1 & 1 & -1 & 1 & 1 & -1 \\ -1 & -1 & 1 & 1 & 1 & 1 \\ 1 & -1 & 1 & 1 & -1 & -1 \\ -1 & 1 & 1 & 1 & -1 & 1 \end{pmatrix}, $$

$$ C_7 = \frac{1}{\sqrt{6}} \begin{pmatrix} 1 & 1 & 1 & 1 & -1 & -1 & -1 \\ 1 & 1 & 1 & 1 & 1 & -1 & 1 \\ 1 & -1 & -1 & 1 & -1 & 1 & 1 \\ -1 & 1 & 1 & 1 & 1 & -1 & -1 \\ 1 & 1 & 1 & -1 & -1 & 1 & -1 \\ -1 & -1 & 1 & -1 & -1 & -1 & 1 \end{pmatrix}, $$

$$ C_8 = \frac{1}{\sqrt{6}} \begin{pmatrix} 1 & 1 & -1 & 1 & 1 & -1 & 1 & -1 \\ 1 & -1 & -1 & -1 & 1 & -1 & 1 & 1 \\ -1 & 1 & 1 & -1 & 1 & -1 & -1 & 1 \\ 1 & 1 & -1 & -1 & 1 & 1 & -1 & -1 \\ 1 & 1 & -1 & 1 & 1 & 1 & -1 & 1 \\ -1 & 1 & -1 & 1 & -1 & -1 & -1 & -1 \end{pmatrix}, $$

$$ C_9 = \frac{1}{\sqrt{7}} \begin{pmatrix} -1 & -1 & -1 & -1 & -1 & -1 & -1 & -1 \\ -1 & -1 & 1 & -1 & -1 & -1 & 1 & -1 \\ -1 & -1 & -1 & 1 & -1 & 1 & -1 & 1 \\ -1 & 1 & 1 & -1 & 1 & -1 & -1 & 1 \\ -1 & -1 & 1 & 1 & 1 & -1 & 1 & -1 \\ -1 & 1 & 1 & 1 & 1 & 1 & 1 & 1 \\ -1 & -1 & 1 & -1 & -1 & 1 & 1 & 1 \\ -1 & -1 & 1 & -1 & -1 & 1 & 1 & 1 \end{pmatrix}. $$

## APPENDIX B
### QUBO - ISING MAPPING

One uses the change of variable $\sigma_i = 1 - 2b_i$ in the expression of the objective function 10 without taking into account the irrelevant constant $||\mathbf{y}||^2$:

$$ \|\mathbf{y} - \mathbf{C} \cdot \operatorname{diag}(\mathbf{w}) \cdot \mathbf{b}\|^2 \sim -2 \sum_{i=1}^{N} (\mathbf{y} \cdot w_i \mathbf{c}_i) \left(\frac{1}{2}(1 - \sigma_i)\right) + \sum_{i,j=1}^{N} (w_i \mathbf{c}_i \cdot w_j \mathbf{c}_j) \left(\frac{1}{2}(1 - \sigma_i)\right) \left(\frac{1}{2}(1 - \sigma_j)\right). \quad (51) $$

We throw away the constant terms (ie. not depending on the $\sigma$'s) to obtain:

$$ \|\mathbf{y} - \mathbf{C} \cdot \operatorname{diag}(\mathbf{w}) \cdot \mathbf{b}\|^2 \sim \sum_{i=1}^{N} \sigma_i \left( w_i \mathbf{y} \cdot \mathbf{c}_i - \frac{1}{2} w_i \mathbf{c}_i \cdot \sum_j w_j \mathbf{c}_j \right) + \frac{1}{4} \sum_{i,j=1}^{N} (w_i w_j \mathbf{c}_i \cdot \mathbf{c}_j) \sigma_i \sigma_j. \quad (52) $$

Using that $\sigma_i^2 = 1$, the diagonal contribution of the double sum is also an irrelevant constant. Since the sum is symmetric under $i \leftrightarrow j$, we can use $\sum_{i \neq j} \equiv 2 \sum_{i<j}$. Hence the expression of the couplings in the main text:

$$ J_{ij} = -\frac{1}{2}(w_i w_j)(\mathbf{c}_i \cdot \mathbf{c}_j). \quad (53) $$

As for the local fields, one can directly identify:

$$ h_i = -w_i \mathbf{c}_i \cdot \left( \mathbf{y} - \frac{1}{2} \sum_j w_j \mathbf{c}_j \right). \quad (54) $$

It matches the expression of the main text:

$$ h_i = -w_i \mathbf{c}_i \cdot \tilde{\mathbf{y}}, \quad (55) $$

with the signal $\tilde{\mathbf{y}}$ expressed as in Eq. 14
REFERENCES

[1] P. Popovski, C. Stefanovic, J. J. Nielsen, et al., "Wireless Access in Ultra-Reliable Low-Latency Communication (URLLC)", en, *IEEE Transactions on Communications*, vol. 67, no. 8, pp. 5783-5801, Aug. 2019. [Online]. Available: https://ieeexplore.ieee.org/document/8705373/.

[2] G. A. Akpakwu, B. J. Silva, G. P. Hancke, and A. M. Abu-Mahfouz, "A Survey on 5G Networks for the Internet of Things: Communication Technologies and Challenges", en, *IEEE Access*, vol. 6, pp. 3619-3647, 2018. [Online]. Available: http://ieeexplore.ieee.org/document/8141874/.

[3] Y. L. Lee, D. Qin, L.-C. Wang, and G. H. Sim, "6G Massive Radio Access Networks: Key Applications, Requirements and Challenges", en, *IEEE Open Journal of Vehicular Technology*, vol. 2, pp. 54-66, 2021. [Online]. Available: https://ieeexplore.ieee.org/document/9295376/.

[4] L. Chetot, M. Egan, and J.-M. Gorce, "Active User Detection and Channel Estimation for Grant-Free Random Access with Gaussian Correlated Activity", en, in *2023 IEEE 97th Vehicular Technology Conference (VTC2023-Spring)*, Florence, Italy: IEEE, Jun. 2023, pp. 1-6. [Online]. Available: https://ieeexplore.ieee.org/document/10199877/.

[5] A. K. Fletcher, S. Rangan, and V. K. Goyal, *On-Off Random Access Channels: A Compressed Sensing Framework*, arXiv:0903.1022 [cs, math] version: 2, Mar. 2009. [Online]. Available: http://arxiv.org/abs/0903.1022.

[6] L. Dai, B. Wang, Z. Ding, Z. Wang, S. Chen, and L. Hanzo, "A Survey of Non-Orthogonal Multiple Access for 5G", en, *IEEE Communications Surveys & Tutorials*, vol. 20, no. 3, pp. 2294-2323, 2018. [Online]. Available: https://ieeexplore.ieee.org/document/8357810/.

[7] J. Xiao, X. Ma, and S. W. McLaughlin, "Quantifying information rate losses with zero-forcing and maximum-likelihood detectors", en, in *2010 IEEE International Conference on Acoustics, Speech and Signal Processing*, Dallas, TX, USA: IEEE, 2010, pp. 3342-3345. [Online]. Available: http://ieeexplore.ieee.org/document/5496010/.

[8] C. Goursaud, A. Julien-Vergonjanne, C. Aupetit-Berthelemot, J.-P. Cances, and J.-M. Dumas, "DS-OCDMA receivers based on parallel interference cancellation and hard limiters", en, *IEEE Transactions on Communications*, vol. 54, no. 9, pp. 1663-1671, Sep. 2006. [Online]. Available: http://ieeexplore.ieee.org/document/1703824/.

[9] L. Zdeborová and F. Krzakala, "Statistical physics of inference: Thresholds and algorithms", *Advances in Physics*, vol. 65, no. 5, pp. 453-552, Sep. 2016, arXiv:1511.02476 [cond-mat, stat]. [Online]. Available: http://arxiv.org/abs/1511.02476.

[10] P. Hauke, H. G. Katzgraber, W. Lechner, H. Nishimori, and W. D. Oliver, "Perspectives of quantum annealing: Methods and implementations", *Reports on Progress in Physics*, vol. 83, no. 5, p. 054401, May 2020, arXiv:1903.06559 [quant-ph]. [Online]. Available: http://arxiv.org/abs/1903.06559.

[11] T. Kadowaki and H. Nishimori, "Quantum annealing in the transverse Ising model", en, *Physical Review E*, vol. 58, no. 5, pp. 5355-5363, Nov. 1998. [Online]. Available: https://link.aps.org/doi/10.1103/PhysRevE.58.5355.

[12] Y. Seki and H. Nishimori, "Quantum annealing with antiferromagnetic fluctuations", *Physical Review E*, vol. 85, no. 5, p. 051112, May 2012, arXiv:1203.2418 [cond-mat, physics:quant-ph]. [Online]. Available: http://arxiv.org/abs/1203.2418.

[13] L. T. Brady, L. Kocia, P. Bienias, A. Bapat, Y. Kharkov, and A. V. Gorshkov, *Behavior of Analog Quantum Algorithms*,

arXiv:2107.01218 [quant-ph], Jul. 2021. [Online]. Available: http://arxiv.org/abs/2107.01218.

[14] J. Roland and N. J. Cerf, "Quantum Search by Local Adiabatic Evolution", *Physical Review A*, vol. 65, no. 4, p. 042308, Mar. 2002, arXiv:quant-ph/0107015. [Online]. Available: http://arxiv.org/abs/quant-ph/0107015.

[15] N. Karafolas and D. Uttamchandani, "Optical CDMA system using bipolar codes and based on narrow passband optical filtering and direct detection", en, *IEEE Photonics Technology Letters*, vol. 7, no. 9, pp. 1072-1074, Sep. 1995. [Online]. Available: http://ieeexplore.ieee.org/document/414706/.

[16] T. Marzetta, E. Larsson, H. Yang, and H. Ngp, *Fundamentals of Massive MIMO*, Cambridge University Press. 2016.

[17] L. Chetot, "Activity models and Bayesian estimation algorithms for wireless grant-free random access",

[18] T. Albash and D. A. Lidar, "Adiabatic Quantum Computing", *Reviews of Modern Physics*, vol. 90, no. 1, p. 015002, Jan. 2018, arXiv:1611.04471 [quant-ph]. [Online]. Available: http://arxiv.org/abs/1611.04471.

[19] S. Morita and H. Nishimori, "Mathematical Foundation of Quantum Annealing", *Journal of Mathematical Physics*, vol. 49, no. 12, p. 125210, Dec. 2008, arXiv:0806.1859 [cond-mat, physics:quant-ph]. [Online]. Available: http://arxiv.org/abs/0806.1859.

[20] E. Farhi, J. Goldstone, S. Gutmann, and M. Sipser, *Quantum Computation by Adiabatic Evolution*, arXiv:quant-ph/0001106, Jan. 2000. [Online]. Available: http://arxiv.org/abs/quant-ph/0001106.

[21] C. Monroe, W. C. Campbell, L.-M. Duan, et al., "Programmable Quantum Simulations of Spin Systems with Trapped Ions", *Reviews of Modern Physics*, vol. 93, no. 2, p. 025001, Apr. 2021, arXiv:1912.07845 [cond-mat, physics:quant-ph]. [Online]. Available: http://arxiv.org/abs/1912.07845.

[22] S. Boixo, T. F. Rønnow, S. V. Isakov, et al., "Quantum annealing with more than one hundred qubits", *Nature Physics*, vol. 10, no. 3, pp. 218-224, Mar. 2014, arXiv:1304.4595 [quant-ph]. [Online]. Available: http://arxiv.org/abs/1304.4595.

[23] E. Gabdulsattarov, K. Rabie, X. Li, and G. Nauryzbayev, *Towards Quantum Annealing for Multi-user NOMA-based Networks*, arXiv:2301.03872 [cs, eess, math], Jan. 2023. [Online]. Available: http://arxiv.org/abs/2301.03872.

[24] Y. Seki and H. Nishimori, *Quantum Annealing with Antiferromagnetic Transverse Interactions for the Hopfield Model*, arXiv:1410.0450 [cond-mat, physics:quant-ph], Aug. 2015. [Online]. Available: http://arxiv.org/abs/1410.0450.

[25] Y. Kimura and H. Nishimori, "Rigorous convergence condition for quantum annealing", *Journal of Physics A: Mathematical and Theoretical*, vol. 55, no. 43, p. 435302, Oct. 2022, arXiv:2207.12096 [quant-ph]. [Online]. Available: http://arxiv.org/abs/2207.12096.

[26] D. Duchemin, L. Chetot, J.-M. Gorce, and C. Goursaud, "Coded random access for massive MTC under statistical channel knowledge", en, in *2019 IEEE 20th International Workshop on Signal Processing Advances in Wireless Communications (SPAWC)*, Cannes, France: IEEE, Jul. 2019, pp. 1-5. [Online]. Available: https://ieeexplore.ieee.org/document/8815491/.