# Hybrid Beamforming and Combining for Millimeter Wave Full Duplex Massive
MIMO Interference Channel

Chandan Kumar Sheemar and Dirk Slock  Communications Systems Department,
EURECOM, France email<EMAIL_ADDRESS>
###### Abstract

Full-Duplex (FD) communication can revolutionize wireless communications as it
avoids using independent channels for bi-directional communications. In this
work, we generalize the point-to-point FD communication in millimeter wave
(mmWave) band consisting of K-pair of massive MIMO FD nodes operating
simultaneously. To enable the coexistence of massive MIMO FD links cost-
efficiently, we present novel joint hybrid beamforming (HYBF) and combining
scheme for weighted sum-rate (WSR) maximization. The proposed algorithm is
based on alternative optimization based on the minorization-maximization
method. Moreover, we present a novel SI and massive MIMO interference channel
aware power allocation scheme to include the optimal power control. Simulation
results show significant performance improvement compared to a traditional
bidirectional half-duplex communication system.

###### Index Terms:

Hybrid beamforming, Full Duplex, Millimeter Wave, massive MIMO Interference
channel

Full-duplex (FD) has to potential to double the performance of a wireless
communication system as it allows simultaneous transmission and reception in
the same frequency band. It avoids using two independent channels for bi-
directional communication by allowing more flexibility in spectrum
utilization, improving data security, and reducing the air interface latency
and delay issues [1, 2, 3]. To enable FD operation, Self-Interference (SI) is
a major challenge to deal with to achieve FD operation, which could be around
$110$ dB compared to the received signal of interest.

However, continuous advancement in the SI cancellation (SIC) techniques has
made the FD operation feasible. The research towards FD communication in the
mmWave band (from $30$ to $300$ GHz) has recently started, which offers much
wider bandwidths and results to be a vital resource for future wireless
communications. However, shifting the potential of FD to the mmWave is much
more challenging, as the signal can suffer from shadowing effects, higher
Doppler spreads, rapid channel fluctuations, intermittent connectivity and low
signal-to-noise ratio (SNR). To overcome the propagation challenges, FD MIMO
nodes need to equipped with a massive number of antennas [4, 1]. Therefore, to
enable FD communication in mmWave cost-efficiently, we must rely on efficient
hybrid beamforming (HYBF) designs, consisting of large dimensional analog
processing and low dimensional digital processing.

In [5], HYBF for an FD relay assisted mmWave macro-cell scenario is
investigated. In [4], the authors proposed a novel HYBF and combining for a
point-to-point FD communication. In [6], a machine-learning-based HYBF for a
one-way mmWave FD relay is investigated. In [7], transmit beamforming for two
massive MIMO nodes is presented. In [8], HYBF for a bidirectional point-to-
point OFDM FD system is available. In [9], a novel HYBF design for an FD
mmWave MIMO relay is proposed. In [10], a novel HYBF design with one hybrid
uplink and one hybrid downlink multi-antenna half-duplex (HD) user is
proposed. In [11], HYBF and combining for a multi-user MIMO uplink and
downlink users is proposed. In [12], a HYBF for FD integrated access and
backhaul is proposed. Note that the contributions on the point-to-point mmWave
massive MIMO FD are limited only to a single communication link [4, 6, 7, 8].

In this paper, we generalize the point-to-point FD communication to a $K$-pair
of massive MIMO nodes in mmWave. All the nodes are assumed to be equipped with
a massive number of antennas and only a limited number of radio-frequency (RF)
chains. The coexistence of multiple FD nodes in mmWave leads to an FD massive
MIMO interference channel in which each node suffers from SI and interference
from all the other nodes. To enable FD communication in such a challenging
scenario, we present a novel HYBF and combining design based on minorization-
maximization [13] for weighted sum-rate (WSR) maximization. Moreover, we
present an optimal power allocation scheme for FD massive MIMO nodes, which
are both SI and massive MIMO interference channel aware. Simulation results
show that the proposed HYBF and combining design exhibit significant
performance improvement over a fully digital massive MIMO HD communication
system with only a limited number of RF chains 111Notation: Boldface lower and
upper case characters denote vectors and matrices, respectively, and
$E\\{\cdot\\},\mbox{Tr}\\{\cdot\\},(\cdot)^{H},\mathbf{I}$, and
$\mathbf{D}_{d}$ represent expectation, trace, Hermitian transpose, identity
matrix and the $d$ dominant vector selection matrix, respectively. The
operators $\mbox{vec}(\mathbf{X}),\mbox{unvec}(\mathbf{x})$ stacks the column
of $\mathbf{X}$ into a vector $\mathbf{x}$, stacks back the vector
$\mathbf{x}$ into a matrix $\mathbf{X}$ and $\angle\mathbf{X}$ normalize the
amplitude to be unit-modulus. .

In summary the contributions of our work are

  1. 1.

Generalization of the point-to-point massive MIMO FD communication in mmWave
to $K$-pair links.

  2. 2.

Introduction of a WSR maximization problem over a massive MIMO interference
channel for HYBF and combining.

  3. 3.

Novel HYBF and combining design.

  4. 4.

Optimal SI and massive MIMO interference channel aware power allocation
scheme.

##  I SYSTEM MODEL

In this paper, we consider a setup of $K$-pair mmWave massive MIMO FD point-
to-point communication links as shown in Figure 1.

Figure 1: Massive MIMO interference channel in the mmWave consisting of Hybrid
nodes with separate transmit (TX) and receive (RX) array.

We distinguish the total nodes into two sets, left and right, denoted with
$\mathcal{L}=\\{1_{l},...,K_{l}\\}$ and $\mathcal{R}=\\{1_{r},...,K_{r}\\}$,
respectively. Nodes involved in the communication link $1\leq i\leq K$ are
denoted with $i_{a}$, in which the subscript $a\in\mathcal{L}$ or
$\in\mathcal{R}$. We consider a multi-stream approach and let
$\mathbf{s}_{i_{l}}\in\mathbb{C}^{d_{i_{l}}\times 1}$ denote the white and
unitary variance data-streams transmitted from node $i_{l}\in\mathcal{L}$
intended for node $i_{r}\in\mathcal{R}$. Let
$\mathbf{V}_{i_{l}}\in\mathbb{C}^{M_{i_{l}}^{t}\times d_{i_{l}}}$ and
$\mathbf{G}_{i_{l}}\in\mathbb{C}^{N_{i_{l}}^{t}\times M_{i_{l}}^{t}}$ denote
the digital and the analog beamformer at node $i_{l}\in\mathcal{L}$,
respectively. Let $\mathbf{F}_{i_{l}}\in\mathbb{C}^{M_{i_{l}}^{r}\times
N_{i_{l}}^{r}}$ denote the analog combiner at node $i_{l}\in\mathcal{L}$ for
the data streams $\mathbf{s}_{i_{r}}$ transmitted from node
$i_{r}\in\mathcal{R}$. Let $M_{i_{l}}^{t}$ and $M_{i_{l}}^{r}$ denote the
number of transmit and receive RF chains for node ${i_{l}}\in\mathcal{L}$,
respectively. Let $N_{i_{l}}^{t}$ and $N_{i_{l}}^{r}$ denote the total number
of transmit and receive antennas for node $i_{l}\in\mathcal{L}$, respectively.
The signal received at node $i_{l}\in\mathcal{L}$ for the case of massive MIMO
FD interference channel can be written as

$\displaystyle\mathbf{y}_{i_{l}}=\mathbf{F}_{i_{l}}\mathbf{H}_{{i_{l}},{i_{r}}}\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}}\mathbf{s}_{i_{r}}+\mathbf{F}_{i_{l}}\mathbf{H}_{{i_{l}},{i_{l}}}\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}}\mathbf{s}_{i_{l}}+\mathbf{F}_{i_{l}}\mathbf{n}_{i_{l}}+$
(1) $\displaystyle\sum_{\begin{subarray}{c}m_{l}\in\mathcal{L}\\\ m_{l}\neq
i_{l}\end{subarray}}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},m_{l}}\mathbf{G}_{m_{l}}\mathbf{V}_{m_{l}}\mathbf{s}_{m_{l}}+\sum_{\begin{subarray}{c}m_{r}\in\mathcal{R}\\\
m_{r}\neq
i_{r}\end{subarray}}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},m_{r}}\mathbf{G}_{m_{r}}\mathbf{V}_{m_{r}}\mathbf{s}_{m_{r}}$

where $\mathbf{n}_{i_{l}}\sim\mathcal{CN}(0,\sigma_{i_{l}}^{2}\mathbf{I})$
denote the noise vector at node $i$,
$\mathbf{H}_{i_{l},i_{r}}\in\mathbb{C}^{N_{i_{l}}^{r}\times N_{i_{r}}^{t}}$,
and $\mathbf{H}_{i_{l},i_{l}}\in\mathbb{C}^{N_{i_{l}}^{r}\times
N_{i_{l}}^{t}}$ denote the channel between the node $i_{l}\in\mathcal{L}$ and
$i_{r}\in\mathcal{R}$, and the SI channel for node $i_{l}$, respectively. The
matrices $\mathbf{H}_{i_{l},m_{l}}\in\mathbb{C}^{N_{i_{l}}^{r}\times
N_{m_{l}}^{t}}$ and
$\mathbf{H}_{i_{l},m_{r}}\in\mathbb{C}^{N_{i_{l}}^{r}\times N_{m_{r}}^{t}}$
denote the interference channels from node $m_{l}\in\mathcal{L},m_{l}\neq
i_{l}$ and from node $m_{r}\in\mathcal{R},m_{r}\neq i_{r}$, respectively. Let
$\mathbf{T}_{i_{l},i_{r}}=\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}}\mathbf{V}_{i_{r}}^{H}\bm{G}_{i_{r}}^{H}$
denote the transmit covariance matrix from node $i_{r}\in\mathcal{R}$ intended
for node $i_{l}\in\mathcal{L}$. Let ($\mathbf{R}_{i_{l}}$) and
$\mathbf{R}_{\overline{i_{l}}}$ denote the (signal plus) interference and
noise covariance matrix received at node $i_{l}\in\mathcal{L}$. The covariance
matrix $\mathbf{R}_{i_{l}}$ can be written as

$\displaystyle\mathbf{R}_{i_{l}}$
$\displaystyle=\underbrace{\mathbf{F}_{i_{l}}\mathbf{H}_{{i_{l}},{i_{r}}}\mathbf{T}_{i_{l},i_{r}}\mathbf{H}_{{i_{l}},{i_{r}}}^{H}\mathbf{F}_{i_{l}}^{H}}_{\triangleq\mathbf{S}_{i_{l}}}+\mathbf{F}_{i_{l}}\mathbf{H}_{{i_{l}},{i_{l}}}\mathbf{T}_{i_{r},i_{l}}\mathbf{H}_{{i_{l}},{i_{l}}}^{H}\mathbf{F}_{i_{l}}^{H}$
(2)
$\displaystyle+\sigma_{i_{l}}^{2}\mathbf{F}_{i_{l}}\mathbf{F}_{i_{l}}^{H}+\sum_{\begin{subarray}{c}m_{l}\in\mathcal{L},m_{l}\neq
i_{l}\end{subarray}}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},m_{l}}\mathbf{T}_{m_{r},m_{l}}\mathbf{H}_{i_{l},m_{l}}^{H}\mathbf{F}_{i_{l}}^{H}$
$\displaystyle+\sum_{\begin{subarray}{c}m_{r}\in\mathcal{R},m_{r}\neq
i_{r}\end{subarray}}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},m_{r}}\mathbf{T}_{m_{l},m_{r}}\mathbf{H}_{i_{l},m_{r}}^{H}\mathbf{F}_{i_{l}}^{H}$

and $\mathbf{R}_{\overline{i_{l}}}$ can be written as
$\mathbf{R}_{\overline{i_{l}}}=\mathbf{R}_{i_{l}}-\mathbf{S}_{i_{l}}$, where
$\mathbf{S}_{i_{l}}$ denotes the useful signal part for node
$i_{l}\in\mathcal{L}$. In the mmWave, the channel $\mathbf{H}_{i_{l},i_{r}}$
can be modelled as

$\mathbf{H}_{i_{l},i_{r}}=\sqrt{\frac{N_{i_{l}}^{r}N_{i_{r}}^{t}}{N_{c}N_{p}}}\sum_{n_{c}=1}^{N_{c}}\sum_{n_{p}=1}^{N_{p}}\alpha_{i_{l},i_{r}}^{n_{p},n_{c}}\mathbf{a}_{i_{l}}(\phi_{i_{l}}^{n_{p},n_{c}})\mathbf{a}_{i_{r}}^{T}(\theta_{i_{r}}^{n_{p},n_{c}}),$
(3)

where $N_{c}$ and $N_{p}$ denote the number of clusters and number of rays,
respectively, $\alpha_{i_{l},i_{r}}^{n_{p},n_{c}}\sim\mathcal{CN}(0,1)$ is a
complex Gaussian random variable with amplitudes and phases distributed
according to the Rayleigh and uniform distribution, respectively, and
$\mathbf{a}_{i_{l}}(\phi_{i_{l}}^{n_{p},n_{c}})$ and
$\mathbf{a}_{i_{r}}^{T}(\theta_{i_{r}}^{n_{p},n_{c}})$ denote the receive and
transmit antenna array response with angle of arrival (AoA)
$\phi_{i_{l}}^{n_{p},n_{c}}$ and angle of departure (AoD)
$\theta_{i_{r}}^{n_{p},n_{c}}$, respectively. Also, the channel matrices
$\mathbf{H}_{i_{l},m_{r}}$ and $\mathbf{H}_{i_{l},m_{l}}$ can be modelled as
(3). The SI channel can be modelled with the Rician fading channel model given
as

$\mathbf{H}_{i_{l},i_{l}}=\sqrt{\frac{\kappa_{i_{l}}}{\kappa_{i_{l}}+1}}\mathbf{H}_{i_{l}}^{LoS}+\sqrt{\frac{1}{\kappa_{i_{l}}+1}}\mathbf{H}_{i_{l}}^{Ref},$
(4)

where $\kappa_{i_{l}}$, $\mathbf{H}_{i_{l}}^{LoS}$ and
$\mathbf{H}_{i_{l}}^{Ref}$ denote the Rician factor, the line-of-sight (LoS)
component matrix and the matrix of reflected components, respectively, at node
$i_{l}\in\mathcal{L}$. The channel matrix $\mathbf{H}_{i_{l}}^{Ref}$ can be
modelled as in (3). The elements of matrix $\mathbf{H}_{i_{l}}^{LoS}$ at
$m$-th row and $n$-th column can be modelled as

$\mathbf{H}_{i_{l}}^{LoS}(m,n)=\frac{\rho}{r_{m,n}}e^{-j2\pi\frac{r_{m,n}}{\lambda}}$
(5)

where $\rho$ is the power normalization constant which assure that
$\mathbb{E}(||\mathbf{H}_{i_{l}}^{LoS}(m,n)||_{F}^{2})=N_{i_{l}}^{r}N_{i_{l}}^{t}$
and $r_{m,n}$ is the distance parameter between the $m$-th receive and the
$n$-th transmit antenna. The WSR maximization problem for a massive MIMO
interference channel under the total sum-power constraint and the unit-modulus
phase shifters over the mmWave FD massive MIMO interference channel can be
stated as

$\underset{\begin{subarray}{c}\mathbf{G},\mathbf{V},\\\
\mathbf{F}\end{subarray}}{\text{max}}\quad\sum_{i_{l}\in\mathcal{L}}w_{i_{l}}\mbox{ln
det}(\mathbf{R}_{\overline{{i_{l}}}}^{-1}\mathbf{R}_{{i_{l}}})+\sum_{{i_{r}}\in\mathcal{R}}w_{i_{r}}\mbox{ln
det}(\mathbf{R}_{\overline{{i_{r}}}}^{-1}\mathbf{R}_{{i_{r}}}).$ (6a)
$\text{s.t.}\quad\mbox{Tr}(\mathbf{G}_{a}\mathbf{V}_{a}\mathbf{V}_{a}^{H}\mathbf{G}_{a}^{H})\preceq
p_{a},\quad\forall a\in\mathcal{L}\;\mbox{or}\;\in\mathcal{R},$ (6b)
$\quad\quad\quad|\mathbf{G}_{a}(m,n)|^{2}=1,\quad\forall\;m,n,\quad\forall
a\in\mathcal{L}\;\mbox{or}\;\in\mathcal{R},$ (6c)
$\quad\quad\quad|\mathbf{F}_{a}(m,n)|^{2}=1,\quad\forall\;m,n,\quad\forall
a\in\mathcal{L}\;\mbox{or}\;\in\mathcal{R}.$ (6d)

The scalars $w_{i_{l}}$ and $w_{i_{r}}$ denote the rate weights for node
$i_{l}\in\mathcal{L}$ and $i_{r}\in\mathcal{R}$, respectively, and $p_{a}$
denote the total sum-power constraint for node
$a\in\mathcal{L}\;\mbox{or}\;\in\mathcal{R}$. The constraints $\eqref{c2}$ and
$\eqref{c3}$ denote the unit-modulus phase-shifter constraint on node $a$ for
the analog beamformer and combiner, respectively. In the problem statement,
$\mathbf{G},\mathbf{V},$ and $\mathbf{F}$, denote the collection of the analog
and digital beamformers and the analog combiners, respectively.

_Remark-_ Note that (6), stated as $\mbox{ln det}(\cdot)$, is not affected by
the digital combiners. They can be chosen to be the well known MMSE combiners
and the achieved rate would not be affected (Please see (4)-(9) [14]).

##  II Problem Simplification

The problem (6) is non-concave in terms of the covariance matrices for the
weighted rates of nodes for which the covariance matrices generate
interference. Given a non-concave structure of the problem due to
interference, it makes the searching for global optima extremely challenging.

To render a feasible solution, we leverage the minorization-maximization [13]
approach. Note that the WSR can be written as

$WSR=WSR_{L}+WSR_{R}=\sum_{i_{l}\in\mathcal{L}}WR_{i_{l}}+\sum_{i_{r}\in\mathcal{R}}WR_{i_{r}}$
(7)

where $WR$ denotes the weighted-rate of one particular user, denoted with the
subscript. Note that $WR_{i_{l}}$ is only concave in
$\mathbf{T}_{i_{l},i_{r}}$ (point-to-point link between nodes $i_{l}$ and
$i_{r}$) and the remaining WRs are non-concave in $\mathbf{T}_{i_{l},i_{r}}$.
Let $\overline{i_{l}}$( $\overline{i_{r}}$) denote the summation in set
$\mathcal{L}$ ( $\mathcal{R}$) with the element $i_{l}$ ( $i_{r}$). Since a
linear function is simultaneously convex and concave, difference of convex
(DC) programming introduces the first order Taylor series expansion of
$WSR_{\overline{i_{l}}}$ and $WSR_{R}$ in $\mathbf{T}_{i_{l},i_{r}}$, around
$\hat{\mathbf{T}}_{i_{l},i_{r}}$, (i.e. all $\mathbf{T}_{i_{l},i_{r}}$). Let
$\hat{\mathbf{T}}$ denote the set containing all such
$\hat{\mathbf{T}}_{i_{l},i_{r}}$. The linearized tangent expression for the FD
point-to-point link for the transmit covariance matrices of node
$i_{l}\in\mathcal{L}$ and $i_{r}\in\mathcal{L}$ can be written by computing
the following gradients

$\hat{\mathbf{A}}_{i_{r}}=-\frac{\partial\mbox{WSR}_{\overline{i_{l}}}}{\partial\hat{\mathbf{T}}_{i_{l},i_{r}}}\Bigr{|}_{\hat{\mathbf{T}}},\quad\hat{\mathbf{B}}_{i_{r}}=-\frac{\partial\mbox{WSR}_{R}}{\partial\hat{\mathbf{T}}_{i_{l},i_{r}}}\Bigr{|}_{\hat{\mathbf{T}}},$
(8)
$\hat{\mathbf{A}}_{i_{l}}=-\frac{\partial\mbox{WSR}_{\overline{i_{r}}}}{\partial{i_{r},i_{l}}}\Bigr{|}_{\hat{\mathbf{T}}},\quad\hat{\mathbf{B}}_{i_{l}}=-\frac{\partial\mbox{WSR}_{L}}{\partial\mathbf{T}_{i_{r},i_{l}}}\Bigr{|}_{\hat{\mathbf{T}}}.$
(9)

where $\hat{\mathbf{A}}_{a}$ and $\hat{\mathbf{B}}_{a}$, for
$a\;\in\;\mathcal{L}\;\mbox{or}\;\mathcal{R}$, denote the linearization with
respect to the same set and the other set, respectively. The gradients can be
computed by using the matrix differentiation properties, which leads to the
expressions

$\hat{\mathbf{A}}_{i_{r}}=\sum_{\begin{subarray}{c}m_{l}\in\mathcal{L},m_{l}\neq
i_{l}\end{subarray}}\mathbf{H}_{m_{l},i_{r}}^{H}\mathbf{F}_{m_{l}}^{H}(\mathbf{R}_{\overline{m_{l}}}^{-1})-\mathbf{R}_{m_{l}}^{-1}))\mathbf{F}_{m_{l}}\mathbf{H}_{m_{l},i_{r}}$
(10a)
$\hat{\mathbf{B}}_{i_{r}}=\sum_{\begin{subarray}{c}n_{r}\in\mathcal{R}\end{subarray}}\mathbf{H}_{n_{r},i_{r}}^{H}\mathbf{F}_{n_{r}}^{H}(\mathbf{R}_{\overline{n_{r}}}^{-1})-\mathbf{R}_{n_{r}}^{-1}))\mathbf{F}_{n_{r}}\mathbf{H}_{n_{r},i_{r}},$
(10b)
$\hat{\mathbf{A}}_{i_{l}}=\sum_{\begin{subarray}{c}m_{r}\in\mathcal{R},m_{r}\neq
i_{r}\end{subarray}}\mathbf{H}_{m_{r},i_{l}}^{H}\mathbf{F}_{m_{r}}^{H}(\mathbf{R}_{\overline{m_{r}}}^{-1})-\mathbf{R}_{m_{r}}^{-1}))\mathbf{F}_{m_{r}}\mathbf{H}_{m_{r},i_{r}},$
(10c)
$\hat{\mathbf{B}}_{i_{l}}=\sum_{\begin{subarray}{c}n_{l}\in\mathcal{L}\end{subarray}}\mathbf{H}_{n_{l},i_{l}}^{H}\mathbf{F}_{n_{l}}^{H}(\mathbf{R}_{\overline{n_{l}}}^{-1})-\mathbf{R}_{n_{l}}^{-1}))\mathbf{F}_{n_{l}}\mathbf{H}_{n_{l},i_{l}}.$
(10d)

Note that the rate of node $i_{r}$ depends on the transmit covariance matrix
from node $i_{l}$ and vice-versa, and the gradients $\hat{\mathbf{B}}_{i_{r}}$
and $\hat{\mathbf{B}}_{i_{l}}$ take into account also the SI generated at the
node $i_{l}$ and $i_{r}$, respectively. Let $\lambda_{i_{l}}$ and
$\lambda_{i_{r}}$ denote the Lagrange multipliers associated with the sum-
power constraint for the node $i_{l}$ and $i_{r}$, respectively. Considering
the unconstrained analog part, dropping the constant terms and re-
parametrizing back the transmit covariance matrices as a function of the
digital and the analog beamformers, augmenting the WSR cost function with sum-
power constraint, leads to the Lagrangian (11), given at the top of the next
page. Note that (11) does not contain the unit-modulus constraints, which will
be incorporated later.

$\displaystyle\mathbb{L}=$
$\displaystyle\sum_{i_{l}\in\mathcal{L}}\lambda_{i_{l}}p_{i_{l}}+\sum_{i_{r}\in\mathcal{R}}\lambda_{i_{r}}p_{i_{r}}+\sum_{i_{l}\in\mathcal{L}}w_{i_{l}}\mbox{ln
det}(\mathbf{I}+\mathbf{V}_{i_{r}}^{H}\mathbf{G}_{i_{r}}^{H}\mathbf{H}_{i_{l},i_{r}}^{H}\mathbf{F}_{i_{l}}^{H}\mathbf{R}_{\overline{i_{l}}}^{-1}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},i_{r}}\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}})-\mbox{Tr}(\mathbf{V}_{i_{r}}^{H}\mathbf{G}_{i_{r}}^{H}(\hat{\mathbf{A}}_{i_{r}}+\hat{\mathbf{B}}_{i_{r}}+\lambda_{i_{r}}\mathbf{I})$
(11)
$\displaystyle\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}})+\sum_{i_{r}\in\mathcal{R}}w_{i_{r}}\mbox{ln
det}(\mathbf{I}+\mathbf{V}_{i_{l}}^{H}\mathbf{G}_{i_{l}}^{H}\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}\mathbf{H}_{i_{r},i_{l}}\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}})-\mbox{Tr}(\mathbf{V}_{i_{l}}^{H}\mathbf{G}_{i_{l}}^{H}(\hat{\mathbf{A}}_{i_{l}}+\hat{\mathbf{B}}_{i_{l}}+\lambda_{i_{l}}\mathbf{I})\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}}),$

##  III Hybrid Beamforming and Combining

This sections presents a novel HYBF and combining design based on the problem
simplification in (11).

###  III-A Digital beamformers

To optimize the digital for node $i_{l}\in\mathcal{L}$ and
$i_{r}\in\mathcal{R}$, we take the derivative of the (11) with respect to
$\mathbf{V}_{i_{r}}$ and $\mathbf{V}_{i_{l}}$, respectively, which leads to
the following Karush–Kuhn–Tucker (KKT) conditions

$\displaystyle\mathbf{G}_{i_{r}}^{H}\mathbf{H}_{i_{l},i_{r}}^{H}\mathbf{F}_{i_{l}}^{H}\mathbf{R}_{\overline{i_{l}}}^{-1}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},i_{r}}\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}}\big{(}\mathbf{I}+\mathbf{V}_{i_{r}}^{H}\mathbf{G}_{i_{r}}^{H}\mathbf{H}_{i_{l},i_{r}}^{H}\mathbf{F}_{i_{l}}^{H}$
(12a)
$\displaystyle\mathbf{R}_{\overline{i_{l}}}^{-1}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},i_{r}}\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}}\big{)}^{-1}=\mathbf{G}_{i_{r}}^{H}\big{(}\hat{\mathbf{A}}_{i_{r}}+\hat{\mathbf{B}}_{i_{r}}+\lambda_{i_{r}}\mathbf{I}\big{)}\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}},$
$\displaystyle\mathbf{G}_{i_{l}}^{H}\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}\mathbf{H}_{i_{r},i_{l}}\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}}\big{(}\mathbf{I}+\mathbf{V}_{i_{l}}^{H}\mathbf{G}_{i_{l}}^{H}\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}$
(12b)
$\displaystyle\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}\mathbf{H}_{i_{r},i_{l}}\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}}\big{)}^{-1}=\mathbf{G}_{i_{l}}^{H}\big{(}\hat{\mathbf{A}}_{i_{l}}+\hat{\mathbf{B}}_{i_{l}}+\lambda_{i_{l}}\mathbf{I}\big{)}\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}},$

######  Theorem 1.

The optimal digital beamformers $\mathbf{V}_{i_{l}}$ and $\mathbf{V}_{i_{r}}$
are given by a generalized dominant eigenvectors of the pairs

$\displaystyle\mathbf{V}_{i_{r}}=\mathbf{D}_{d_{i_{l}}}($
$\displaystyle\mathbf{G}_{i_{r}}^{H}\mathbf{H}_{i_{l},i_{r}}^{H}\mathbf{F}_{i_{l}}^{H}\mathbf{R}_{\overline{i_{l}}}^{-1}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},i_{r}}\mathbf{G}_{i_{r}},$
(13)
$\displaystyle\mathbf{G}_{i_{r}}^{H}\big{(}\hat{\mathbf{A}}_{i_{r}}+\hat{\mathbf{B}}_{i_{r}}+\lambda_{r}\mathbf{I}\big{)}\mathbf{G}_{i_{r}}),$
$\displaystyle\mathbf{V}_{i_{l}}=\mathbf{D}_{d_{i_{r}}}($
$\displaystyle\mathbf{G}_{i_{l}}^{H}\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}\mathbf{H}_{i_{r},i_{l}}\mathbf{G}_{i_{l}},$
(14)
$\displaystyle\mathbf{G}_{i_{l}}^{H}\big{(}\hat{\mathbf{A}}_{i_{l}}+\hat{\mathbf{B}}_{i_{l}}+\lambda_{i_{l}}\mathbf{I}\big{)}\mathbf{G}_{i_{l}}),$

where $\mathbf{D}_{d_{i_{l}}}$ $(\mathbf{D}_{d_{i_{r}}})$ selects the
$d_{i_{l}}$ $(d_{i_{r}})$ dominant generalized eigenvectors.

###### Proof.

Given the analog part fixed, we proof the result only for
$\mathbf{V}_{i_{l}}$. The result for $\mathbf{V}_{i_{r}}$ is then
straightforward based on the proof for $\mathbf{V}_{i_{l}}$. The proof relies
on simplifying (11) with fixed analog part until the Hadamard’s inequality
applies. The Cholesky decomposition of the matrix
$\mathbf{G}_{i_{l}}^{H}\big{(}\hat{\mathbf{A}}_{i_{l}}+\hat{\mathbf{B}}_{i_{l}}+\lambda_{i_{l}}\mathbf{I}\big{)}\mathbf{G}_{i_{l}})$
is written as $\bm{L}_{i_{l}}\bm{L}_{i_{l}}^{H}$ where $\bm{L}_{i_{l}}$ is a
lower-triangular Cholesky factor. Upon defining
$\tilde{\bm{V}_{i_{l}}}=\bm{L}_{i_{l}}^{H}\bm{V}_{i_{l}}$.A similar reasoning
follows also for
$\mathbf{G}_{i_{l}}^{H}\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}\mathbf{H}_{i_{r},i_{l}}\mathbf{G}_{i_{l}}$.
By using the result provided in Proposition 1 [15], with the matrices
$\mathbf{G}_{i_{l}}^{H}\big{(}\hat{\mathbf{A}}_{i_{l}}+\hat{\mathbf{B}}_{i_{l}}+\lambda_{i_{l}}\mathbf{I}\big{)}\mathbf{G}_{i_{l}})$
and
$\mathbf{G}_{i_{l}}^{H}\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}\mathbf{H}_{i_{r},i_{l}}\mathbf{G}_{i_{l}}$,
it follows immediately that the optimal digital beamformer
$\mathbf{V}_{i_{l}}$ are given as a dominant generalized eigenvectors of these
two matrices. The same reasoning follows also for the digital beamformer
$\mathbf{V}_{i_{r}}$. ∎

Note the GEV solution provides the optimal beamforming directions but not the
optimal power allocation. Therefore, to include the optimal power allocation,
we normalize the columns of the digital beamformers to be unit norm.

###  III-B Analog beamformer

To optimize the analog beamformers for node $i_{l}\in\mathcal{L}$ and
$i_{r}\in\mathcal{R}$, we take the derivative of (11) for $\mathbf{G}_{i_{l}}$
and $\mathbf{G}_{i_{r}}$, which leads to the following KKT conditions

$\displaystyle\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}\mathbf{H}_{i_{r},i_{l}}\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}}\mathbf{V}_{i_{l}}^{H}\big{(}\mathbf{I}+\mathbf{V}_{i_{l}}^{H}\mathbf{G}_{i_{l}}^{H}\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}$
(15a)
$\displaystyle\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}\mathbf{H}_{i_{r},i_{l}}\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}}\big{)}^{-1}=\big{(}\hat{\mathbf{A}}_{i_{l}}+\hat{\mathbf{B}}_{i_{l}}+\lambda_{i_{l}}\mathbf{I}\big{)}\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}}\mathbf{V}_{i_{l}}^{H},$
$\displaystyle\mathbf{H}_{i_{l},i_{r}}^{H}\mathbf{F}_{i_{l}}^{H}\mathbf{R}_{\overline{i_{l}}}^{-1}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},i_{r}}\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}}\mathbf{V}_{i_{r}}^{H}\big{(}\mathbf{I}+\mathbf{V}_{i_{r}}^{H}\mathbf{G}_{i_{r}}^{H}\mathbf{H}_{i_{l},i_{r}}^{H}\mathbf{F}_{i_{l}}^{H}$
(15b)
$\displaystyle\mathbf{R}_{\overline{i_{l}}}^{-1}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},i_{r}}\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}}\big{)}^{-1}=\big{(}\hat{\mathbf{A}}_{i_{r}}+\hat{\mathbf{B}}_{i_{r}}+\lambda_{i_{r}}\mathbf{I}\big{)}\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}}\mathbf{V}_{i_{r}}^{H}.$

Note that to find the optimal analog beamformer, the KKT conditions are not
resolveable for $\mathbf{G}_{i_{l}}$ and $\mathbf{G}_{i_{r}}$. To do so, we
apply the identity
$\mbox{vec}(\mathbf{A}\mathbf{X}\mathbf{B})=(\mathbf{B}^{T}\otimes\mathbf{A})\mbox{vec}(\mathbf{X})$,
which shapes the the KKT conditions as

$\displaystyle[(\mathbf{V}_{i_{r}}\mathbf{V}_{i_{r}}^{H}\big{(}\mathbf{I}+\mathbf{V}_{i_{r}}^{H}\mathbf{G}_{i_{r}}^{H}\mathbf{H}_{i_{l},i_{r}}^{H}\mathbf{F}_{i_{l}}^{H}\mathbf{R}_{\overline{i_{l}}}^{-1}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},i_{r}}\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}}\big{)}^{-1})^{T}$
(16a)
$\displaystyle\otimes\mathbf{H}_{i_{l},i_{r}}^{H}\mathbf{F}_{i_{l}}^{H}\mathbf{R}_{\overline{i_{l}}}^{-1}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},i_{r}}]\mbox{vec}(\mathbf{G}_{i_{r}})$
$\displaystyle=[(\mathbf{V}_{i_{r}}\mathbf{V}_{i_{r}}^{H})^{T}\otimes\big{(}\hat{\mathbf{A}}_{i_{r}}+\hat{\mathbf{B}}_{i_{r}}+\lambda_{i_{l}}\mathbf{I}\big{)}]\mbox{vec}(\mathbf{G}_{i_{r}}),$
$\displaystyle[(\mathbf{V}_{i_{l}}\mathbf{V}_{i_{l}}^{H}\big{(}\mathbf{I}+\mathbf{V}_{i_{l}}^{H}\mathbf{G}_{i_{l}}^{H}\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}\mathbf{H}_{i_{r},i_{l}}\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}}\big{)}^{-1})^{T}$
(16b)
$\displaystyle\otimes\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}\mathbf{H}_{i_{r},i_{l}}]\mbox{vec}(\mathbf{G}_{i_{l}})$
$\displaystyle=[(\mathbf{V}_{i_{l}}\mathbf{V}_{i_{l}}^{H})^{T}\otimes\big{(}\hat{\mathbf{A}}_{i_{l}}+\hat{\mathbf{B}}_{i_{l}}+\lambda_{i_{r}}\mathbf{I}\big{)}]\mbox{vec}(\mathbf{G}_{i_{l}}).$

Given the vectorized analog beamformer, by following a similar proof as for
the digital beamformers in the Theorem (1), it can be easily seen that the
optimal unconstrained vectorized analog combiners
$\mbox{vec}(\mathbf{G}_{i_{r}})$ and $\mbox{vec}(\mathbf{G}_{i_{l}})$ are
given by

$\displaystyle\mbox{vec}(\mathbf{G}_{i_{r}})=$
$\displaystyle\mathbf{D}_{1}((\mathbf{V}_{i_{r}}\mathbf{V}_{i_{r}}^{H}\big{(}\mathbf{I}+\mathbf{V}_{i_{r}}^{H}\mathbf{G}_{i_{r}}^{H}\mathbf{H}_{i_{l},i_{r}}^{H}\mathbf{F}_{i_{l}}^{H}\mathbf{R}_{\overline{i_{l}}}^{-1}\mathbf{F}_{i_{l}}$
(17a)
$\displaystyle\mathbf{H}_{i_{l},i_{r}}\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}}\big{)}^{-1})^{T}\otimes\mathbf{H}_{i_{l},i_{r}}^{H}\mathbf{F}_{i_{l}}^{H}\mathbf{R}_{\overline{i_{l}}}^{-1}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},i_{r}},$
$\displaystyle(\mathbf{V}_{i_{r}}\mathbf{V}_{i_{r}}^{H})^{T}\otimes\big{(}\hat{\mathbf{A}}_{i_{l}}+\hat{\mathbf{B}}_{i_{l}}+\lambda_{i_{r}}\mathbf{I}\big{)}),$
$\displaystyle\mbox{vec}(\mathbf{G}_{i_{l}})=$
$\displaystyle\mathbf{D}_{1}((\mathbf{V}_{i_{l}}\mathbf{V}_{i_{l}}^{H}\big{(}\mathbf{I}+\mathbf{V}_{i_{l}}^{H}\mathbf{G}_{i_{l}}^{H}\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}$
(17b)
$\displaystyle\mathbf{H}_{i_{r},i_{l}}\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}}\big{)}^{-1})^{T}\otimes\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}\mathbf{H}_{i_{r},i_{l}},$
$\displaystyle(\mathbf{V}_{i_{r}}\mathbf{V}_{i_{r}}^{H})^{T}\otimes\big{(}\hat{\mathbf{A}}_{i_{r}}+\hat{\mathbf{B}}_{i_{r}}+\lambda_{i_{l}}\mathbf{I}\big{)}),$

which needs to be reshaped into correct dimensions with the operation
$\mbox{unvec}(\cdot)$. Furthermore, to meet the unit modulus constraint, we
apply the operation $\angle\mathbf{G}_{i_{l}}$ and $\angle\mathbf{G}_{i_{r}}$,
which preserves only the phase part.

###  III-C Analog combiner

To optimize the analog combiner, we first define the received covariance
matrices at the antenna level as $\mathbf{R}_{i_{l}}^{ant}$ and
$\mathbf{R}_{i_{r}}^{ant}$ (obtained from (2) by simply omitting the analog
combiner). After the analog combining stage, we have the following expression
for the covariance matrices
$\mathbf{R}_{i_{l}}=\mathbf{F}_{i_{l}}\mathbf{R}_{i_{l}}^{ant}\mathbf{F}_{i_{l}}^{H}$
and
$\mathbf{R}_{i_{r}}=\mathbf{F}_{i_{r}}\mathbf{R}_{i_{r}}^{ant}\mathbf{F}_{i_{r}}^{H}$.
Note that in the minorization-maximization approach, we linearize with respect
to the WRs for which each transmit covariance matrix generate interference.
However, as the combiner do not generate any interference and does not have
the sum-power constraint, we can directly solve (6), which result to be
concave for the analog combiners. Namely, we can first write

$\begin{aligned}
\underset{\begin{subarray}{c}\mathbf{F}\end{subarray}}{\text{max}}&\sum_{i_{l}\in\mathcal{L}}w_{i_{l}}[\mbox{ln
det}(\mathbf{F}_{i_{l}}\mathbf{R}_{i_{l}}^{ant}\mathbf{F}_{i_{l}}^{H})-\mbox{ln
det}(\mathbf{F}_{i_{l}}\mathbf{R}_{\overline{i_{l}}}^{ant}\bm{F}_{i_{l}}^{H})]\\\
&+\sum_{{i_{r}}\in\mathcal{R}}w_{i_{r}}[\mbox{ln
det}(\bm{F}_{i_{r}}\mathbf{R}_{i_{r}}^{ant}\mathbf{F}_{i_{r}}^{H})-\mbox{ln
det}(\mathbf{F}_{i_{l}}\mathbf{R}_{\overline{i_{l}}}^{ant}\mathbf{F}_{i_{l}}^{H})]\end{aligned}.$
(18)

in which both the terms are purely concace, in contrast to (11) in which the
trace terms was only linear. To optimize the analog combiners, we take the
derivative of (18) with respect to $\mathbf{F}_{i_{l}}$ and
$\mathbf{F}_{i_{r}}$, which leads to the following KKT conditions

$\displaystyle w_{i_{l}}$
$\displaystyle\mathbf{R}_{i_{l}}^{ant}\mathbf{F}_{i_{l}}^{H}\big{(}\bm{F}_{i_{l}}\mathbf{R}_{i_{l}}^{ant}\mathbf{F}_{i_{l}}^{H}\big{)}^{-1}$
(19a) $\displaystyle-
w_{i_{l}}\mathbf{R}_{\overline{i_{l}}}^{ant}\mathbf{F}_{i_{l}}^{H}\big{(}\mathbf{F}_{i_{l}}\mathbf{R}_{\overline{i_{l}}}|^{ant}\bm{F}_{i_{l}}^{H}\big{)}^{-1}=0,$
$\displaystyle w_{i_{r}}$
$\displaystyle\mathbf{R}_{i_{r}}^{ant}\bm{F}_{i_{r}}^{H}\big{(}\bm{F}_{i_{r}}\mathbf{R}_{i_{r}}^{ant}\bm{F}_{i_{r}}^{H}\big{)}^{-1}$
(19b) $\displaystyle-
w_{i_{r}}\mathbf{R}_{\overline{i_{r}}}^{ant}\bm{F}_{i_{r}}^{H}\big{(}\bm{F}_{i_{r}}\mathbf{R}_{\overline{i_{r}}}|^{ant}\bm{F}_{i_{r}}^{H}\big{)}^{-1}=0.$

Given the structure of the KKT conditions, it is immediate to see that the
optimal analog combiner is given as a dominant generalized eigenvectors as

$\bm{F}_{i_{l}}=\mathbf{D}_{N_{i_{l}}^{r}}(\mathbf{R}_{i_{l}}^{ant},\mathbf{R}_{\overline{i_{l}}}^{ant}),\quad\bm{F}_{i_{r}}=\mathbf{D}_{N_{i_{r}}^{r}}(\mathbf{R}_{i_{r}}^{ant},\mathbf{R}_{\overline{i_{r}}}^{ant}),$
(20)

from which we select $M_{i_{l}}^{r}$ and $M_{i_{r}}^{r}$ rows, respectively.
Given the optimal analog combiner, operation $\angle\bm{F}_{i_{l}}$ and
$\angle\bm{F}_{i_{r}}$ is required to meet the unit modulus constraint.

###  III-D Optimal power allocation

Given the optimal beamformers, in this section we present a novel power
allocation scheme for the massive MIMO FD interference channel.

Let $\mathbf{P}_{i_{l}}$ and $\mathbf{P}_{i_{r}}$ denote the power allocation
matrix for the node $i_{l}\in\mathcal{L}$ and $i_{r}\in\mathcal{R}$,
respectively. Let
$\mathbf{\Sigma}_{i_{l}}^{(1)}$,$\mathbf{\Sigma}_{i_{l}}^{(2)}$,$\mathbf{\Sigma}_{i_{r}}^{(1)}$
and $\mathbf{\Sigma}_{i_{r}}^{(2)}$ be defined as

$\mathbf{\Sigma}_{i_{l}}^{(1)}=\mathbf{V}_{i_{r}}^{H}\mathbf{G}_{i_{r}}^{H}\mathbf{H}_{i_{l},i_{r}}^{H}\mathbf{F}_{i_{l}}^{H}\mathbf{R}_{\overline{i_{l}}}^{-1}\mathbf{F}_{i_{l}}\mathbf{H}_{i_{l},i_{r}}\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}},$
(21a)
$\mathbf{\Sigma}_{i_{l}}^{(2)}=\mathbf{V}_{i_{r}}^{H}\mathbf{G}_{i_{r}}^{H}(\hat{\mathbf{A}}_{i_{l}}+\hat{\mathbf{B}}_{i_{l}}+\lambda_{i_{r}}\mathbf{I})\mathbf{G}_{i_{r}}\mathbf{V}_{i_{r}},$
(21b)
$\mathbf{\Sigma}_{i_{r}}^{(1)}=\mathbf{V}_{i_{l}}^{H}\mathbf{G}_{i_{l}}^{H}\mathbf{H}_{i_{r},i_{l}}^{H}\mathbf{F}_{i_{r}}^{H}\mathbf{R}_{\overline{i_{r}}}^{-1}\mathbf{F}_{i_{r}}\mathbf{H}_{i_{r},i_{l}}\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}},$
(21c)
$\mathbf{\Sigma}_{i_{r}}^{(2)}=\mathbf{V}_{i_{l}}^{H}\bm{G}_{i_{l}}^{H}(\hat{\mathbf{C}}_{j}+\hat{\mathbf{D}}_{j}+\lambda_{i_{l}}\mathbf{I})\mathbf{G}_{i_{l}}\mathbf{V}_{i_{l}}.$
(21d)

Given (21), the problem (11) with respect to the power matrices can be stated
as

$\displaystyle\underset{\begin{subarray}{c}\mathbf{P}_{i_{l}}\end{subarray}}{\text{max}}\sum_{i_{l}\in\mathcal{L}}w_{i_{l}}\mbox{ln
det}(\bm{I}+\mathbf{\Sigma}_{i_{l}}^{(1)}\mathbf{P}_{i_{l}})-\mbox{Tr}(\mathbf{\Sigma}_{i_{l}}^{(2)}\mathbf{P}_{i_{l}}),$
(22a)
$\displaystyle\underset{\begin{subarray}{c}\mathbf{P}_{i_{l}}\end{subarray}}{\text{max}}\sum_{i_{r}\in\mathcal{R}}w_{i_{r}}\mbox{ln
det}(\mathbf{I}+\mathbf{\Sigma}_{i_{r}}^{(1)}\mathbf{P}_{i_{r}})-\mbox{Tr}(\mathbf{\Sigma}_{i_{r}}^{(2)}\mathbf{P}_{i_{r}}).$

To include the optimal power allocation, we take the derivative of (22a) for
the power matrices $\mathbf{P}_{i_{l}}$ and $\mathbf{P}_{i_{r}}$, and solving
KKT conditions for the powers leads to the following optimal power allocation
matrices

$\mathbf{P}_{i_{l}}=(w_{i_{l}}{\mathbf{\Sigma}_{i_{l}}^{(2)}}^{-1}-{\mathbf{\Sigma}_{i_{l}}^{(1)}}^{-1})^{+},\quad\mathbf{P}_{i_{r}}=(w_{i_{r}}{\mathbf{\Sigma}_{i_{r}}^{(2)}}^{-1}-{\mathbf{\Sigma}_{i_{r}}^{(1)}}^{-1})^{+}.$
(23)

where $(x)^{+}=\mbox{max}\\{0,x\\}$. To meet the sum-power constraint, we
search the optimal multipliers satisfying the power constraints while updating
the powers with (23). The multipliers $\lambda_{i_{l}}$ and $\lambda_{i_{r}}$
should be such that the Lagrange dual function (11) is finite and the values
of the multipliers should be strictly positive. Let $\mathbf{\Lambda}$ and
$\mathbf{P}$ denote the collection of multipliers and powers. Formally, the
multipliers’ search problem can be stated as

$\displaystyle\underset{\mathbf{\Lambda}}{\text{min}}\,\,\underset{\bm{P}}{\text{max}}\quad\mathbb{L}\big{(}\mathbf{\Lambda},\bm{P}\big{)},\quad\quad\mbox{s.t.}\quad\mathbf{\Lambda}\succeq
0.$ (24)

The dual function
$\underset{\bm{P}}{\text{max}}\;\mathbb{L}(\mathbf{\Lambda},\bm{P})$ is the
pointwise supremum of a family of functions of $\mathbf{\Lambda}$, it is
convex [16] and the globally optimal values for $\mathbf{\Lambda}$ can be
found by using any of the numerous convex optimization techniques. In this
work,we adopt the bisection method. Let $\underline{\lambda_{i_{a}}}$ and
$\overline{\lambda_{i_{a}}}$, for $i_{a}\in\mathcal{L}$ or
$i_{a}\in\mathcal{R}$, denote the upper and lower bound for the Lagrange
multiplier search for node $i_{a}$. Let $[0,\overline{\lambda_{i_{a}}^{max}}]$
denote the maximum search interval for the multipliers for node $i_{a}$. The
WSR maximization hybrid beamforming design for the massive MIMO interference
channel is given in Algorithm 1

Algorithm 1 Hybrid Beamforming and Combining

Given: The CSI and rate weights.
Initialize: All the beamformers and combiners.
Repeat until convergence

for $\forall\;i_{l}\in\mathcal{L}$ and $\forall\;i_{r}\in\mathcal{R}$

Compute $\mathbf{A}_{i_{l}}$ and $\mathbf{A}_{i_{l}}$ with (10a) and (10c)

Compute $\mathbf{B}_{i_{l}}$ and $\mathbf{B}_{i_{r}}$ with (10b) and (10d)

Compute $\bm{G}_{i_{l}}$ with (17b), do $\mbox{unvec}(\bm{G}_{i_{l}})$ and
$\angle\bm{G}_{i_{l}}$

Compute $\bm{G}_{i_{r}}$ with (17a), do $\mbox{unvec}(\bm{G}_{i_{r}})$ and
$\angle\bm{G}_{i_{r}}$

Compute $\bm{F}_{i_{l}}$, $\bm{F}_{i_{r}}$ with (20) and do
$\angle\bm{F}_{i_{l}}$, $\angle\bm{F}_{i_{r}}$

Compute $\bm{V}_{i_{l}}$ and $\bm{V}_{i_{r}}$ with (14) and (13)

Set $\underline{\lambda_{i_{l}}}=0$ and
$\overline{\lambda_{i_{l}}}=\mu_{i_{max}}$ $\forall i_{l}\in\mathcal{L}_{l}$

Repeat until convergence $\forall i_{l}\in\mathcal{L}$

set
$\lambda_{i_{l}}=(\underline{\lambda_{i_{l}}}+\overline{\lambda_{i_{l}}})/2$

Compute $\bm{P}_{i_{l}}$ with (23),

If constraint for $\lambda_{i_{l}}$ is violated,

set $\underline{\lambda_{i_{l}}}=\lambda_{i_{l}}$, else
$\overline{\lambda_{i_{l}}}=\lambda_{i_{l}}$,

Repeat until convergence $\forall i_{r}\in\mathcal{R}$

set
$\lambda_{i_{r}}=(\underline{\lambda_{i_{r}}}+\overline{\lambda_{i_{r}}})/2$.

Compute $\bm{P}_{i_{l}}$ with (23),

If constraint for $\lambda_{i_{r}}$ is violated,

set $\underline{\lambda_{i_{r}}}=\lambda_{i_{r}}$, else
$\overline{\lambda_{i_{r}}}=\lambda_{i_{r}}$,

###  III-E Short Convergence Proof

For the WSR cost function (6), we construct its minorizer as in (11), which
restates the WSR maximization as a concave problem, when the remaning
variables are fixed (in the gradients. The minorizer results to be a touching
lower bound for the original WSR problem, therefore we can write

$\displaystyle\underline{WSR}$
$\displaystyle=\sum_{i_{l}\in\mathcal{L}}w_{i_{l}}\mbox{ln
det}(\bm{I}+\mathbf{\Sigma}_{i_{l}}^{(1)}\mathbf{P}_{i_{l}})-\mbox{Tr}(\mathbf{\Sigma}_{i_{l}}^{(2)}\mathbf{P}_{i_{l}}),$
(25) $\displaystyle+\sum_{i_{r}\in\mathcal{R}}w_{i_{r}}\mbox{ln
det}(\mathbf{I}+\mathbf{\Sigma}_{i_{r}}^{(1)}\mathbf{P}_{i_{r}})-\mbox{Tr}(\mathbf{\Sigma}_{i_{r}}^{(2)}\mathbf{P}_{i_{r}})$

The minorizer, which is concave in transmit covariance matrices, still has the
same gradient of the original WSR and hence the KKT conditions are not
affected. Now reparameterizing the transmit covariance matrices in terms of
beamformer with the optimal power matrices and adding the power constraints to
the minorizer, we get the Lagrangian (11). By updating all the beamformers and
combiners as a dominant generalized eigenvector, leads to an increase in the
WSR [15] at each iteration, thus assuring convergence.

##  IV SIMULATION RESULTS

This section presents the simulation results for the proposed HYBF and
combining design for WSR maximization in a mmWave FD massive MIMO interference
channel.

We consider a scenario of $2$ point-to-point communication links consisting of
$4$ nodes operating in the mmWave. We assume that all the nodes are equipped
with uniform linear arrays, with transmit antennas placed at half-wavelength.
The angle of arrival $\phi_{i_{l}}^{n_{p},n_{c}}$ and angle of departure
$\theta_{i_{l}}^{n_{p},n_{c}}$ are chosen to be uniformly distributed in the
interval $[-20^{\circ},20^{\circ}]$. The SI channel is modelled with the
Rician factor $\kappa_{i_{l}}=10^{5}$ dB and the $r_{m,n}$ is modelled as in
[4], with distance between the transmit and receive array of $20$ cm and with
a relative angle of $90^{\circ}$. The number of clusters and the number of
paths is set to be $N_{c}=3$ and $N_{p}=6$. We define the transmit SNR at node
$i_{l}$ as $\mbox{SNR}_{i_{l}}=p_{i_{l}}/\sigma_{i_{l}}^{2}$ and assume it to
be the same for all the nodes, denoted as SNR. We label the proposed hybrid
beamforming scheme as HYBF and for comparison purpose, we define the following
benchmark schemes

  * •

_Fully Digital FD_ \- with all the FD nodes having number of RF chains equal
to the number of antennas.

  * •

_Fully Digital HD_ \- with all the nodes having number of RF chains equal to
the number of antennas but operating in half-duplex (HD) mode.

Figure 2 shows the achieved WSR as a function of SNR with a different number
of RF chains in comparison with the benchmark schemes. It can be clearly seen
that the proposed hybrid beamforming and combining has the potential to
perform very close to the fully digital FD beamforming design with $32$ RF
chains. Still the proposed algorithm exhibits some gap compared to the fully
digital solution as the analog beamforming stage is constrained and must meet
the unit modulus constraint. We can see that also with $16$RF chains, the
proposed algorithm can considerably outperform the fully digital HD scheme,
with operates with $100$ transmit and receive RF chains. Figure 3 shows the
average WSR as a function of SNR with $64$ transmit and receive antennas. It
can be seen that with the same number of RF chains ($32$)less number of
transmit and receive antennas, the performance of the proposed HYB algorithm
gets strictly close to the fully digital FD scheme.

Figure 2: Average WSR as a function of SNR with transmit and receive antennas
$N_{a}^{t}=N_{a}^{r}=100$ and RF chains $M_{a}^{t}=M_{a}^{r}=32$ or $16$,
$\forall a\in\mathcal{L}$ or $\mathcal{R}$. Figure 3: Average WSR as a
function of SNR with transmit and receive antennas $N_{a}^{t}=N_{a}^{r}=64$
and RF chains $M_{a}^{t}=M_{a}^{r}=32$ or $16$, $\forall a\in\mathcal{L}$ or
$\mathcal{R}$.

##  V CONCLUSIONS

In this paper, we studied the problem of WSR maximization in a mmWave massive
MIMO FD interference channel consisting of $K$ point-to-point communication
links. The simultaneous coexistence of multiple FD nodes leads to an extremely
challenging communication scenario for which novel hybrid beamforming and
combining scheme is proposed. Simulation results show that the proposed design
can perform extremely close to the fully digital FD beamforming design
operating with $100$ antennas, with only $32$ RF chains. Moreover, the
proposed design outperforms considerably the fully digital HD system already
with only $16$ RF chains.

## References

  * [1] P. Rosson, C. K. Sheemar, N. Valecha, and D. Slock, “Towards massive mimo in-band full duplex radio,” in _ISWCS_. IEEE, 2019.
  * [2] C. K. Sheemar and D. Slock, “Receiver design and agc optimization with self interference induced saturation,” in _ICASSP_. IEEE, 2020.
  * [3] C. K. Sheemar and S. Dirk, “Beamforming for bidirectional mimo full duplex under the joint sum power and per antenna power constraints.”
  * [4] K. Satyanarayana, M. El-Hajjar, P.-H. Kuo, A. Mourad, and L. Hanzo, “Hybrid beamforming design for full-duplex millimeter wave communication,” _IEEE Transactions on Vehicular Technology_ , vol. 68, no. 2, pp. 1394–1404, 2018.
  * [5] S. Han, Y. Zhang, W. Meng, C. Li, and Z. Zhang, “Full-duplex relay-assisted macrocell with millimeter wave backhauls: Framework and prospects,” _IEEE Network_ , vol. 33, no. 5, pp. 190–197, 2019.
  * [6] S. Huang, Y. Ye, and M. Xiao, “Learning based hybrid beamforming design for full-duplex millimeter wave systems,” _arXiv preprint arXiv:2004.08285_ , 2020.
  * [7] C. K. Sheemar and D. Slock, “Hybrid beamforming for bidirectional massive mimo full duplex under practical considerations,” in _VTC Spring_ , 2021.
  * [8] C. K. Thomas, C. K. Sheemar, and D. Slock, “Multi-stage/hybrid bf under limited dynamic range for ofdm fd backhaul with mimo si nulling,” in _ISWCS_. IEEE, 2019.
  * [9] Y. Cai, K. Xu, A. Liu, M. Zhao, B. Champagne, and L. Hanzo, “Two-timescale hybrid analog-digital beamforming for mmwave full-duplex mimo multiple-relay aided systems,” _IEEE Journal on Selected Areas in Communications_ , 2020\.
  * [10] I. P. Roberts, J. G. Andrews, and S. Vishwanath, “Hybrid beamforming for millimeter wave full-duplex under limited receive dynamic range,” _arXiv preprint arXiv:2012.11647_ , 2020.
  * [11] C. K. Sheemar and C. K. T. D. Slock, “Practical hybrid beamforming for millimeter wave massive mimo full duplex with limited dynamic range,” _arXiv preprint arXiv:2104.11537_ , 2021.
  * [12] C. K. Sheemar and D. Slock, “Massive mimo mmwave full duplex relay for iab with limited dynamic range,” in _NTMS 2021_. IEEE, 2021.
  * [13] P. Stoica and Y. Selen, “Cyclic minimizers, majorization techniques, and the expectation-maximization algorithm: a refresher,” _IEEE Signal Processing Magazine_ , vol. 21, no. 1, pp. 112–114, 2004.
  * [14] S. S. Christensen, R. Agarwal, E. De Carvalho, and J. M. Cioffi, “Weighted sum-rate maximization using weighted mmse for mimo-bc beamforming design,” _IEEE Transactions on Wireless Communications_ , 2008.
  * [15] S.-J. Kim and G. B. Giannakis, “Optimal resource allocation for mimo ad hoc cognitive radio networks,” _IEEE Transactions on Information Theory_ , vol. 57, no. 5, pp. 3117–3131, 2011.
  * [16] S. Boyd, S. P. Boyd, and L. Vandenberghe, _Convex optimization_. Cambridge university press, 2004.


# Synergies between astroparticle, particle and nuclear physics

Lund University
E-mail Versions of this article have appeared on the CERN EP newsletter
(Summer 2019) and on the ECFA newsletter (September 2019).

###### Abstract:

One overarching objective of science is to further our understanding of the
universe, from its early stages to its current state and future evolution.
This depends on gaining insight on the universe’s most macroscopic components,
for example galaxies and stars, as well as describing its smallest components,
namely elementary particles and nuclei and their interactions. It is clear
that this endeavour requires combined expertise from the fields of
astroparticle physics, particle physics and nuclear physics. Pursuing common
scientific drivers also require mastering challenges related to
instrumentation (e.g. beams and detectors), data acquisition, selection and
analysis, and making data and results available to the broader science
communities. Joint work and recognition of these foundational topics will help
all communities grow towards their individual and common scientific goals. The
talk corresponding to this contribution has been presented during the special
ECFA session of EPS-HEP 2019 focused on the update of the European Strategy of
Particle Physics.

##  1 Introduction

One overarching objective of science is to further our understanding of the
universe, from its early stages to its current state and future evolution.
This depends on gaining insight on the universe’s most macroscopic components,
for example galaxies and stars, as well as describing its smallest components,
namely elementary particles and nuclei and their interactions. It is clear
that this endeavor requires combined expertise from the fields of
astroparticle physics, particle physics and nuclear physics.

A number of the contributions and discussions at the recent Granada meeting
for the update of the European Strategy of Particle Physics, as well as the
contribution at the EPS-HEP ECFA Open Session summarized in this newsletter,
highlighted a growing wish for closer collaboration between ECFA and the
astrophysics (APPEC, https://www.appec.org) and nuclear physics (NuPECC,
http://www.nupecc.org) communities.

Many physics problems where synergies between particle physics, astrophysics
and nuclear physics are required are discussed in the APPEC and NuPECC
strategy documents (see links at the bottom of this piece). Among those, this
contribution focused on the challenge of of elucidating the nature of 27% of
the matter-energy content of the universe, commonly called dark matter.
Pursuing these scientific goals also requires mastering challenges related to
instrumentation (e.g. beams and detectors), data acquisition, selection and
analysis, and making data and results available to the broader science
communities. Joint work and recognition of these foundational topics, also
covered in detail in the contributions by C. Biscari, A. Cattai and G. A.
Stewart in the ECFA newsletter [1], will help all communities grow towards
their individual and common scientific goals. This contribution presented one
of the many common challenges faced by particle physics and astrophysics: the
necessity of dealing with large, sometimes heterogeneous datasets and derive
insight from them in short periods of time.

##  2 New physics discoveries and dark matter

The Large Hadron Collider has yielded the discovery of a new particle, the
Higgs boson. Precision measurements and fits of other quantities in the
Standard Model of Particle Physics guided a search lasting decades after the
conception of the Higgs mechanism. With the European Strategy Update, the
particle physics community is currently deciding what are the best tools to
employ to test the Standard Model: how/where to look next for new physics is a
relevant question in this process, given that hints coming from the Standard
Model itself are not as telling as in the case of the Higgs boson. For this
reason, research directions for physics beyond the Standard Model can be found
in open problems in astrophysics that need a systematic exploration, for
example the determination of the nature of dark matter.

One of the many explanations for this dark matter is that it is composed by
new massive particles that interact only weakly with ordinary matter
particles, or Weakly Interacting Massive Particles (WIMPs). These new
particles can be produced at colliders, as well as detected by direct and
indirect detection astrophysics experiments in space and underground (see Ref.
[2] and links for a basic summary). By producing new particles in the lab,
colliders are well placed to understand the nature of these particle’s
interactions with ordinary matter. The necessary confirmation that these new
particles have also a cosmological origin comes from complementary
observations in direct and indirect detection experiments.

Figure 1:  Comparison of sensitivities of future collider and direct detection
experiments within a simplified model scenario of a WIMP where the interaction
between Standard Model quarks and the dark matter is mediated by a new scalar
particle. If dark matter is composed by particles within this model with a
mass between 10 GeV and 1 TeV, future colliders and direct detection
experiment can confirm each other’s discoveries in the next decades. Adapted
from Ref. [3], with the addition of ATLAS results line.

The WIMP scenario shown in Figure 1 only represents a very simple benchmark in
the landscape of theories on the nature of dark matter. Many other compelling
explanations exist: for example the WIMP can be identified with the lightest,
stable and invisible particle included in many supersymmetric models that also
answer other outstanding questions of the Standard Model. Alternatives to the
WIMP paradigm also exist, for example models where the dark matter particle is
much lighter and has a mass below the GeV, see Fig. 2 for a set of constraints
of current and planned experiments. In these cases, searches at collider and
direct detection experiments are complementary to searches at other planned
dedicated accelerator experiments (e.g. beam dump experiments such as SHIP
[5], NA64 [6], and LDMX [7] among others), as well as underground experiments
using novel sensor technologies (e.g. SENSEI [8] and DAMIC [9]).

Figure 2:  Comparison of sensitivities of future collider and non-collider
experiments for a model of light dark matter being produced via a new dark
boson, showing the complementarity between different kinds of experiments in
different ranges of dark matter particle mass and couplings with the dark
boson. Taken from Ref. [3, 4].

Axions (and axion-like particles) [10, 11] also may be connected to solutions
of the dark matter problem, being the DM particle candidate themselves or the
mediators of the SM-DM interaction. Synergies between many different
experiments and theoretical frameworks is evident in the case of those
particles. Depending on the mass range and coupling of those particles, a
discovery of such particles may occur at the high-luminosity LHC, at lepton
colliders (e.g. Belle II [12]) or at high-precision experiments that search
for these particles directly (e.g. IAXO [13] and ADMX [14]), or measure
fundamental constants sensitive to fifth forces. Interaction between the
members of these different experimental communities, as well as with theory
and astrophysics, are needed to shape the future search program of these
complementary experiments.

In general, connecting results and potential discoveries from different
experiments within a coherent framework requires both particle and
astroparticle physics theory involvement. This effort has been started by the
LHC Dark Matter Working Group [15], where the Astroparticle community wishes
to be further involved (see Ref. [16]). A parallel effort for non-WIMP, non-
collider dark matter and dark sector searches is ongoing within the Physics
Beyond Colliders Working Group.

A connection to nuclear physics in these and other (e.g. beam dump) dark
matter experiments is needed to fully understand instrumental and beam
backgrounds, as well as simulation. Particle and astroparticle experiments
searching for dark matter also benefit from cross-talk in terms of
instrumentation (e.g. sensors and cryogenics) and interpretation of results.

##  3 Data firehoses and shared solutions in high energy physics and
astrophysics

Another example of a common challenge for different fields is the ever-
increasing volume of data available to different fields of research. Examples
of current data firehoses are the LHC, especially in light of the planned
high-luminosity upgrade, and upcoming astrophysics surveys such as LSST [20]
and SKA [21] to name but two. Similar challenges in data acquisition and
recording are present in neutrino physics experiments, in the case of their
supernova detection data streams. In all these cases, a fast and close-to-
real-time analysis of the data is necessary, so that events of interest can be
recorded or investigated further in a timely and cost-effective way, and
common real-time analysis solutions are being investigated and deployed by
multiple experiments.

Another point of contact is when software solutions are shared across fields,
for example in the case of gravitational waves and high energy physics with
the CernVM software appliance [17] and the RUCIO distributed data management
system that is in use by LHC experiments and will be adopted by neutrino
experiments as well [18].

##  4 Collaborative efforts

A number of platforms and fora exist at CERN and in Europe to facilitate
cross-talk among different communities. In addition to the already-mentioned
Dark Matter Working Group, there are (just to name a few) the recently
inaugurated European Center for Astroparticle Physics currently hosted by CERN
[19], the European Science Cluster for Astronomy and Particle physics ESFRI
research infrastructures project, the HEP Software Foundation to facilitate
cooperation and common effort in software and computing, as well as the very
successful CERN neutrino platform.

##  5 Conclusions

The examples brought forward in this EPS-HEP contribution are only a very
limited subset of how the particle, astroparticle and nuclear physics
communities can be answering challenging scientific questions together. Other
topics where synergies exist mentioned during the session were axion-like
particles, the theory and experimental efforts bridging the gap between
nuclear and high energy physics, and the opportunities offered by astrophysics
experiments (e.g. Auger) spanning a much higher and complementary energy
regime with respect to nuclear and particle physics experiments.

Since detector technologies are often common to different communities, the
CERN expertise stemming from the current world-leading collider program can be
reused. Moreover, data collection and analysis benefit from becoming faster,
more efficient and more open: using versatile computing strategies and tools
to solve diverse problems encourages common expertise that lasts beyond a
single experiment.

In conclusion, there is the common wish that the European Strategy process
will facilitate closer collaboration between the particle, astroparticle and
nuclear physics communities, in a context where the design of detectors, data
acquisition systems and computing are an integral part in our quest to
understand the universe.

## References

  * [1] ECFA Newsletter #3, Summer 2019. Online, accessed 09/2019.
  * [2] Searching for Dark Matter with the ATLAS detector, ATLAS experiment feature article, published online on March 2019, accessed 09/2019. Zenodo: http://doi.org/10.5281/zenodo.3594499
  * [3] Richard Keith Ellis et al., “Physics Briefing Book : Input for the European Strategy for Particle Physics Update 2020,” [arXiv:1910.11775 [hep-ex]].
  * [4] J. Beacham et al., J. Phys. G 47, no. 1, 010501 (2020) doi:10.1088/1361-6471/ab4cd2 [arXiv:1901.09966 [hep-ex]].
  * [5] M. Anelli et al. [SHiP Collaboration], “A facility to Search for Hidden Particles (SHiP) at the CERN SPS,” [arXiv:1504.04956 [physics.ins-det]].
  * [6] D. Banerjee et al. [NA64 Collaboration], “Search for invisible decays of sub-GeV dark photons in missing-energy events at the CERN SPS,” Phys. Rev. Lett. 118, no. 1, 011802 (2017) doi:10.1103/PhysRevLett.118.011802 [arXiv:1610.02988 [hep-ex]].
  * [7] T. Åkesson et al. [LDMX Collaboration], “Light Dark Matter eXperiment (LDMX),” [arXiv:1808.05219 [hep-ex]].
  * [8] O. Abramoff et al. [SENSEI Collaboration], “SENSEI: Direct-Detection Constraints on Sub-GeV Dark Matter from a Shallow Underground Run Using a Prototype Skipper-CCD,” Phys. Rev. Lett. 122, no. 16, 161801 (2019) doi:10.1103/PhysRevLett.122.161801 [arXiv:1901.10478 [hep-ex]].
  * [9] A. Aguilar-Arevalo et al. [DAMIC Collaboration], “Constraints on Light Dark Matter Particles Interacting with Electrons from DAMIC at SNOLAB,” Phys. Rev. Lett. 123, no. 18, 181802 (2019) doi:10.1103/PhysRevLett.123.181802 [arXiv:1907.12628 [astro-ph.CO]].
  * [10] M. Tanabashi et al. [Particle Data Group], Phys. Rev. D 98, 030001 (2018) and 2019 update. Review on axions: http://pdg.lbl.gov/2019/reviews/rpp2018-rev-axions.pdf
  * [11] CERN EP Newsletter, “Axion-like particle searches at the LHC,” September 2018. Online, accessed 09/2019.
  * [12] T. Abe et al. [Belle-II Collaboration], “Belle II Technical Design Report,” arXiv:1011.0352 [physics.ins-det].
  * [13] E. Armengaud et al. [IAXO Collaboration], “Physics potential of the International Axion Observatory (IAXO),” JCAP 1906, 047 (2019) doi:10.1088/1475-7516/2019/06/047 [arXiv:1904.09155 [hep-ph]].
  * [14] T. Braine et al. [ADMX Collaboration], “Extended Search for the Invisible Axion with the Axion Dark Matter Experiment,” [arXiv:1910.08638 [hep-ex]].
  * [15] LHC Dark Matter Working Group website, Online, accessed 09/2019.
  * [16] APPEC News, T. Montaruli, “EPPSU and the synergy of Astroparticle and Particle Physics,” September 2018. Online, accessed 09/2019.
  * [17] CERN EP Newsletter, “CernVM 2019 Workshop,” June 2019. Online, accessed 09/2019.
  * [18] CERN News, “Managing scientific data at the exascale with Rucio,” March 2018. Online, accessed 09/2019.
  * [19] CERN News, “A new centre for astroparticle physics theory,” July 2019. Online, accessed 09/2019.
  * [20] Z. Ivezic et al. [LSST Collaboration], “LSST: from Science Drivers to Reference Design and Anticipated Data Products,” Astrophys. J. 873, no. 2, 111 (2019) doi:10.3847/1538-4357/ab042c [arXiv:0805.2366 [astro-ph]].
  * [21] A. Weltman et al., “Fundamental Physics with the Square Kilometre Array,” arXiv:1810.02680 [astro-ph.CO].