SCA-Based Beamforming Optimization for IRS-Enabled Secure
Integrated Sensing and Communication
Vaibhav Kumar1, Marwa Chafii23,
A. Lee Swindlehurst4, Le-Nam Tran1,
and Mark F. Flanagan1

1School of Electrical and Electronic Engineering,
University College Dublin, Belfield, Dublin 4, Ireland

2Engineering Division, New York University (NYU),
Abu Dhabi, UAE

3NYU WIRELESS, NYU Tandon School of Engineering,
New York, USA

3Center for Pervasive Communications and Computing,
University of California, Irvine, CA 92697 USA

Email<EMAIL_ADDRESS><EMAIL_ADDRESS><EMAIL_ADDRESS>
<EMAIL_ADDRESS><EMAIL_ADDRESS>This work was supported by the Irish Research Council under Grant IRCLA/2017/209, and also in part by Science Foundation Ireland under Grant 17/CDA/4786.

Integrated sensing and communication (ISAC) is expected to be offered
as a fundamental service in the upcoming sixth-generation (6G) communications
standard. However, due to the exposure of information-bearing signals
to the sensing targets, ISAC poses unique security challenges. In
recent years, intelligent reflecting surfaces (IRSs) have emerged
as a novel hardware technology capable of enhancing the physical layer
security of wireless communication systems. Therefore, in this paper,
we consider the problem of transmit and reflective beamforming design
in a secure IRS-enabled ISAC system to maximize the beampattern gain
at the target. The formulated non-convex optimization problem is challenging
to solve due to the intricate coupling between the design variables.
Moreover, alternating optimization (AO) based methods are inefficient
in finding a solution in such scenarios, and convergence to a stationary
point is not theoretically guaranteed. Therefore, we propose a novel
successive convex approximation (SCA)-based second-order cone programming
(SOCP) scheme in which all of the design variables are updated simultaneously
in each iteration. The proposed SCA-based method significantly outperforms
a penalty-based benchmark scheme previously proposed in this context.
Moreover, we also present a detailed complexity analysis of the proposed
scheme, and show that despite having slightly higher per-iteration
complexity than the benchmark approach the average problem-solving
time of the proposed method is notably lower than that of the benchmark

Intelligent reflecting surface (IRS), integrated sensing and communication
(ISAC), physical layer security, successive convex approximation (SCA),
second-order cone programming (SOCP).

§ INTRODUCTION

The sixth-generation (6G) wireless standard is being developed not
only to improve the quality of user experience compared to that offered
by the fifth-generation networks, but also to support a range of new
wireless communication services, such as autonomous vehicles, drone
monitoring, human activity recognition, environmental monitoring,
enhanced localization and tracking, and many more. Supporting such
new services will require the integration of communication, sensing
and localization capabilities as fundamental services in a single
network architecture rather than as auxiliary functionalities [1].
Integrated sensing and communication (ISAC) has recently emerged as
a potential enabler in this direction, combining the communication
and sensing capabilities in a single hardware platform using a common
waveform [2, 3]. Preliminary results
have confirmed that ISAC can improve the spectral efficiency of a
network by virtue of exploiting a common hardware, signal processing
and spectral framework, thereby offering a low-cost solution to the
spectrum scarcity problem. Furthermore, by exploiting the possibility
of communication-centric and sensing-centric designs, it can also
enjoy significant coordination gains compared to conventional networks.
However, due to the broadcast nature of wireless channels and the
inclusion of information-bearing signaling in the sensing waveform,
susceptibility to eavesdropping targets poses unique security challenges
in ISAC.

Intelligent reflecting surfaces (IRSs) have recently emerged as a
groundbreaking hardware technology to robustify wireless communication
systems against eavesdroppers via passive beamforming [4].
The benefits of IRS in a secure communication-only multiple-input
multiple-output (MIMO) system have been well established in the literature [5, 6].
Hence, it is worth exploring the advantages of IRSs in an ISAC system
in terms of physical layer security [7]. However,
it is interesting to note that in contrast to the somewhat rich literature
on IRS-aided ISAC systems [8, 9, 10, 11, 12],
there is a dearth of literature on secure IRS-aided ISAC
system design [13, 14]. As one
of the few examples, the authors in [13]
considered an active IRS-aided multiuser multiple-input single-output
(MU-MISO) ISAC system, where the aim was to obtain an optimal beamforming
design that maximizes the achievable secrecy rate of the communication
users while guaranteeing a minimum radar signal-to-interference-plus-noise
ratio (SINR). In [14], the authors considered
the problem of beampattern optimization for an eavesdropping target
in an IRS-enabled MU-MISO secure ISAC system, subject to SINR constraints
at the communication users and information leakage constraints at
the target. Two different scenarios were considered in [14];
in the first scenario, full channel state information (CSI) and the
target location were assumed to be known at the base station (BS),
while imperfect CSI and uncertain target location were assumed in
the second scenario. In this paper, we will focus on the first scenario
only, where the CSI and target location are known at the BS. The more
practical setting in which these quantities are imprecisely known
will be considered in future work. To optimize the transmit and IRS
beamforming in the first scenario, the authors in [14]
proposed a penalty-based alternating optimization (AO) algorithm to
obtain a semi-closed-form solution using Lagrange duality and a majorization-minimization
(MM) algorithm.

System model for IRS-enabled secure ISAC system.

Even though the use of AO in [14] makes the optimization
problem much easier to solve, it may not produce a high-quality solution
because of the complicated interdependence between the design variables [15].
Moreover, as we will show in Section <ref>, the use
of a penalty-based method requires a large number of iterations to
achieve convergence and therefore has a very high problem-solving
time. Note also that a feasible solution is not guaranteed if the
algorithm terminates prematurely. To tackle these issues, in this
paper we propose a successive convex approximation (SCA) based beampattern
optimization scheme which results in a high-performance solution and
also requires a much shorter convergence time. The main contributions
of the paper are listed as follows:

* We propose a provably convergent SCA-based algorithm to maximize the
beampattern gain at the eavesdropping target in the secure ISAC system,
subject to the SINR requirements at the communication users and information
leakage constraints at the target. In contrast to the AO-based scheme
of [14] where all design variables are updated
in an alternating fashion, we derive a second-order cone program (SOCP)
where all of the optimization variables are updated simultaneously
in each iteration.

* We present a complexity analysis of the proposed scheme which demonstrates
that the per-iteration complexity grows as $\mathcal{O}\big(N^{3.5}\big)$,
while that of the benchmark solution is given by $\mathcal{O}\big(N^{3}\big)$ <cit.>,
where $N$ is the number of reflecting elements in the IRS. Although
the per-iteration complexity of the proposed approach is slightly
higher than that of the benchmark scheme, we show that it requires
significantly fewer iterations to converge, resulting in a much shorter
problem-solving time.

* We present extensive numerical results to confirm that the proposed
SCA-based SOCP approach results in a high-performance solution, and
significantly outperforms the penalty-based AO algorithm in [14].

Bold uppercase and lowercase letters are used to denote matrices and
vectors, respectively. By $\mathbb{C}^{M\times N}$, we denote the
vector space of all $M\times N$ complex-valued matrices. By $\mathbf{X}\trans$,
$\mathbf{X}\herm$, $\|\mathbf{X}\|$, $\Re\{\mathbf{X}\}$ and $\Im\{\mathbf{X}\}$,
we respectively denote the transpose, conjugate transpose, Frobenius
norm, real and imaginary components of a matrix $\mathbf{X}$. $|x|$
denotes the absolute value of a complex number $x$, and $\diag(x)$
denotes the diagonal matrix whose main diagonal comprises the elements
of $\mathbf{x}$. $\mathcal{O}(\cdot)$ denotes the Bachmann–Landau

§ SYSTEM MODEL AND PROBLEM FORMULATION

Consider the ISAC system shown in Fig. <ref> consisting
of a multi-antenna dual-function base station (BS), one IRS, $K$
single-antenna communication users (denoted by $\text{U}_{k},k\in\mathcal{K}\triangleq\{1,2,\ldots,K\}$),
and one single-antenna eavesdropping target.[Although we consider a single eavesdropping target in this paper,
it is straightforward to use the proposed algorithm for a system with
multiple eavesdropping targets.] Let $L$ and $N$ denote the number of antennas at the BS and the
number of elements in the IRS, respectively. We assume that the BS
transmits a linear superposition of radar and information signals
for the purpose of joint sensing and communication. The signal vector
transmitted from the BS is given by

\begin{equation}
\mathbf{s}=\sum\nolimits _{k\in\mathcal{K}}\mathbf{x}_{k}w_{k}+\sum\nolimits _{l\in\mathcal{L}}\hat{\mathbf{x}}_{l}\hat{w}_{l},\label{eq:txSignal}
\end{equation}

where $w_{k}$ is the communication signal intended for $\text{U}_{k}$
and $\hat{w}_{l}$ is the $l^{\text{th}}$ radar signal with $l\in\mathcal{L}\triangleq\{1,2,\ldots,L\}$.
Moreover, $\mathbf{x}_{k}\in\mathbb{C}^{L\times1}$ and $\hat{\mathbf{x}}_{l}\in\mathbb{C}^{L\times1}$
are the beamforming vectors corresponding to $w_{k}$ and $\hat{w}_{l}$,
respectively. It is assumed that $\mathbb{E}\big\{ w_{k}\big\}=0$,
$\mathbb{E}\big\{|w_{k}|^{2}\big\}=1\ \forall k\in\mathcal{K}$, $\mathbb{E}\big\{\hat{w}_{l}\big\}=0$,
$\mathbb{E}\big\{|\hat{w}_{l}|^{2}\big\}=1\ l\in\mathcal{L}$, and
$\mathbb{E}\big\{ w_{k}\hat{w}_{l}\herm\big\}=0\ \forall k\in\mathcal{K},l\in\mathcal{L}$,
i.e., the communication and radar signals are mutually independent
and uncorrelated.

Why we have to assume $\mathbb{E}\big\{|w_{k}|^{2}\big\}=0\ \forall k\in\mathcal{K}$
and $\mathbb{E}\big\{|\hat{w}_{l}|^{2}\big\}=0\ l\in\mathcal{L}$.
This is only true if signals are zero.

Denoting the BS-IRS, BS-$\text{U}_{k}$ and IRS-$\text{U}_{k}$ links
by $\mathbf{G}\in\mathbb{C}^{N\times L}$, $\mathbf{h}_{\text{D}k}\in\mathbb{C}^{1\times L}$,
and $\mathbf{h}_{\text{R}k}\in\mathbb{C}^{1\times N}$, respectively,
the signal received at $\text{U}_{k}$ is given by

\begin{equation}
\end{equation}

where $\mathbf{h}_{k}\triangleq\mathbf{h}_{\text{D}k}+\mathbf{h}_{\text{R}k}\boldsymbol{\Theta}\mathbf{G}$,
$\boldsymbol{\Theta}\triangleq\diag(\boldsymbol{\theta})$, $\boldsymbol{\theta}=[\theta_{1},\theta_{2},\ldots,\theta_{N}]\trans$,
$\theta_{n}\triangleq\exp(j2\pi\phi_{n})$ with $\phi_{n}\in[0,2\pi)$
denoting the phase shift induced by the $n^{\text{th}}$ IRS element,
and $z_{k}$ denotes the zero-mean complex additive white Gaussian
noise (AWGN) at $\text{U}_{k}$ with variance $\sigma_{k}^{2}$ .

Defining $\mathbf{X}\triangleq[\mathbf{x}_{1},\mathbf{x}_{2},\ldots,\mathbf{x}_{K},\hat{\mathbf{x}}_{1},\hat{\mathbf{x}}_{2},\ldots\hat{\mathbf{x}}_{L}]\in\mathbb{C}^{L\times(K+L)}$,
and $\tilde{\mathbf{x}}_{m}$ as the $m^{\text{th}}$ column of $\mathbf{X}$,
the signal-to-interference-plus-noise ratio (SINR) at $\text{U}_{k}$
to decode the intended message is given by

\begin{equation}
\gamma_{k}=\frac{|\mathbf{h}_{k}\mathbf{x}_{k}|^{2}}{\sigma_{k}^{2}+\sum_{\ell\in\mathcal{M}\setminus\{k\}}|\mathbf{h}_{k}\tilde{\mathbf{x}}_{\ell}|^{2}},\label{eq:defSINR-Uk}
\end{equation}

where $\mathcal{M}\triangleq\{1,2,\ldots,K,K+1,\ldots,K+L\}.$ Similarly,
by denoting the IRS-target channel by $\mathbf{g}_{\text{R}}\in\mathbb{C}^{1\times N}$
and assuming that the BS-target (direct) link is blocked due to obstacles,
the signal received at the target is given by

\begin{equation}
\end{equation}

where $\mathbf{g}\triangleq\mathbf{g}_{\text{R}}\boldsymbol{\Theta}\mathbf{G}\in\mathbb{C}^{1\times L}$,
$\mathbf{g}_{\mathrm{R}}$ is the steering vector
from the IRS in the direction of the target, and $z_{\text{T}}$
is the zero-mean complex AWGN at the target with variance $\sigma_{\text{T}}^{2}$.
Therefore, the SINR at the target to wiretap the signal intended for
$\text{U}_{k}$ is given by

\begin{equation}
\hat{\gamma}_{k}=\frac{|\mathbf{g}\mathbf{x}_{k}|^{2}}{\sigma_{\text{T}}^{2}+\sum\nolimits _{\ell\in\mathcal{M}\setminus\{k\}}|\mathbf{g}\tilde{\mathbf{x}}_{\ell}|^{2}}.\label{eq:defSINR-Tk}
\end{equation}

We assume that all of the channels and the target location are perfectly
known to the BS. The beampattern gain toward the target is then given
by (c.f. [14])

\begin{equation}
\mathscr{G}(\mathbf{X},\boldsymbol{\theta})=\mathbb{E}\{|\mathbf{gs}|^{2}\}=\sum_{m\in\mathcal{M}}|\mathbf{g}\tilde{\mathbf{x}}_{m}|^{2}.\label{eq:defBPG}
\end{equation}

Therefore, the problem of maximizing the beampattern gain toward the
target is given by

\begin{align}
\underset{\mathbf{X},\boldsymbol{\theta}}{\maximize}\  & \mathscr{G}(\mathbf{X},\boldsymbol{\theta}),\label{eq:mainObj}\\
\st\  & \gamma_{k}\geq\Gamma_{k}\ \forall k\in\mathcal{K},\label{eq:SINRC-Uk}\\
& \hat{\gamma}_{k}\leq\hat{\Gamma}_{k}\ \forall k\in\mathcal{K},\label{eq:SINRC-Tk}\\
& \|\mathbf{X}\|\leq\sqrt{P}\label{eq:TPC}\\
& |\theta_{n}|^{2}=1\ \forall n\in\mathcal{N}\triangleq\{1,2,\ldots,N\},\label{eq:UMC}
\end{align}

where (<ref>) ensures that the SINR at $\text{U}_{k}$
is greater than or equal to the predefined threshold $\Gamma_{k}$, (<ref>)
enforces the constraint that the maximum leakage of the information
intended for $\text{U}_{k}$ at the eavesdropping target is below
the tolerance level $\hat{\Gamma}_{k}$ and $P$ is the transmit power
budget at the BS. Note that in a system with heterogeneous secrecy
requirements, considering information leakage constraints results
in a more flexible resource allocation compared to that offered via
imposing constraints on the achievable secrecy rate [16].
It is easy to note that due to coupling between the design variables
$\mathbf{X}$ and $\boldsymbol{\theta}$ in (<ref>)–(<ref>)
and the non-convex constraints in (<ref>), the problem in (<ref>)
is non-convex and challenging to solve.

Hua et al. [14] proposed a penalty-based
dual-loop AO algorithm to obtain a solution to (<ref>).
More specifically, in the inner loop, auxiliary variables were updated
by solving a quadratically-constrained quadratic program (QCQP), the
beamformers (i.e., $\mathbf{X}$) were updated using a bisection search,
and the IRS reflection coefficients ($\boldsymbol{\theta}$) were
updated via MM; the outer loop was used to update the penalty parameter
only. Although the use of auxiliary variables and the penalty method
in [14] resulted in a reformulated optimization
problem where the design variables were decoupled in the constraints,
obtaining a high quality solution is not guaranteed via AO. Additionally,
although the per-iteration complexity of the penalty-based solution
in [14] was $\mathcal{O}\big(N^{3}\big)$, because
of the use of bisection search and the dual-loop structure, the number
of iterations required for convergence is large. This in turn results
in a high problem-solving time because if the iterations are terminated
prematurely (i.e., before the penalty terms becomes nearly zero),
the obtained solution may not be feasible.

§ PROPOSED SOLUTION

In this section, we apply a series of convex approximations to tackle
the non-convexity of (<ref>) and to obtain a high-performance
solution. In this regard, for two arbitrary complex-valued vectors
$\mathbf{u}$ and $\mathbf{v}$, we recall the following (in)equalities
(c.f. <cit.>)

\begin{align}
\|\mathbf{u}\|^{2} & \ \geq2\Re\{\mathbf{v}\herm\mathbf{u}\}-\|\mathbf{v}\|^{2},\label{eq:Inequality}\\
\Re\{\mathbf{u}\herm\mathbf{v}\} & \ =\frac{1}{4}\big(\|\mathbf{u}+\mathbf{v}\|^{2}-\|\mathbf{u}-\mathbf{v}\|^{2}\big),\label{eq:EQ1}\\
\Im\{\mathbf{u}\herm\mathbf{v}\} & \ =\frac{1}{4}\big(\|\mathbf{u}-j\mathbf{v}\|^{2}-\|\mathbf{u}+j\mathbf{v}\|^{2}\big).\label{eq:EQ2}
\end{align}

Next, we note that the term in (<ref>) is neither convex
nor concave. Since we want to maximize the function in (<ref>),
we obtain a corresponding concave lower bound as follows:

\begin{align}
& \mathscr{G}(\mathbf{X},\boldsymbol{\theta})=\sum_{m\in\mathcal{M}}|\mathbf{g}\tilde{\mathbf{x}}_{m}|^{2}\nonumber \\
\overset{(\texttt{a})}{\geq} & \ \sum_{m\in\mathcal{M}}\big[2\Re\{a_{m}^{(i)}\!\!\ \herm\mathbf{g}\tilde{\mathbf{x}}_{m}\}-|a_{m}^{(i)}|^{2}\big]\nonumber \\
\overset{(\texttt{b})}{=} & \ \sum_{m\in\mathcal{M}}\Big[\frac{1}{2}\big\{\|a_{m}^{(i)}\mathbf{g}\herm+\tilde{\mathbf{x}}_{m}\|^{2}\!-\!\|a_{m}^{(i)}\mathbf{g}\herm\!-\tilde{\mathbf{x}}_{m}\|^{2}\big\}\!-\!|a_{m}^{(i)}|^{2}\Big]\nonumber \\
\overset{(\texttt{c})}{\geq} & \ \sum_{m\in\mathcal{M}}\Big[\Re\big\{\mathbf{b}_{m}^{(i)}\!\!\ \herm\big[a_{m}^{(i)}\mathbf{g}\herm+\tilde{\mathbf{x}}_{m}\big]\big\}-\frac{1}{2}\|\mathbf{b}_{m}^{(i)}\|^{2}\nonumber \\
& \qquad\qquad-\frac{1}{2}\|a_{m}^{(i)}\mathbf{g}\herm-\tilde{\mathbf{x}}_{m}\|^{2}-|a_{m}^{(i)}|^{2}\Big]\nonumber \\
\triangleq & \ \sum_{m\in\mathcal{M}}f_{m}(\tilde{\mathbf{x}}_{m},\boldsymbol{\theta};\tilde{\mathbf{x}}_{m}^{(i)},\boldsymbol{\theta}^{(i)}),\label{eq:objLowerBound}
\end{align}

where $\tilde{\mathbf{x}}_{m}^{(i)}$ and $\boldsymbol{\theta}^{(i)}$
denote the value of $\tilde{\mathbf{x}}_{m}$ and $\boldsymbol{\theta}$
in the $i^{\text{th}}$ iteration of the SCA process, respectively.
Moreover, $(\texttt{a})$ and $(\texttt{c})$ follow from (<ref>),
and $(\texttt{b})$ follows from (<ref>). Additionally in (<ref>),
$\mathbf{b}_{m}^{(i)}\triangleq a_{m}^{(i)}\mathbf{g}^{(i)}\!\!\ \herm+\tilde{\mathbf{x}}_{m}^{(i)}$,
and $\mathbf{g}^{(i)}\triangleq\mathbf{g}_{\text{R}}\boldsymbol{\Theta}^{(i)}\mathbf{G}$.
Note that $f_{m}(\tilde{\mathbf{x}}_{m},\boldsymbol{\theta};\tilde{\mathbf{x}}_{m}^{(i)},\boldsymbol{\theta}^{(i)})$
is jointly concave with respect to (w.r.t.) $\tilde{\mathbf{x}}_{m}$
and $\boldsymbol{\theta}$.

Next, we turn our attention to the non-convex constraints in (<ref>).
Using (<ref>), for any $k\in\mathcal{K}$, we can
equivalently represent (<ref>) as follows:

\begin{align}
\frac{1}{\Gamma_{k}}\big|\mathbf{h}_{k}\mathbf{x}_{k}\big|^{2} & \geq\sigma_{k}^{2}+\sum_{\ell\in\mathcal{M}\setminus\{k\}}\big(\wp_{k\ell}^{2}+\bar{\wp}_{k\ell}^{2}\big),\label{eq:cSINR-1-1}\\
\wp_{k\ell} & \geq\big|\Re\big\{\mathbf{h}_{k}\tilde{\mathbf{x}}_{\ell}\big\}\big|\ \forall\ell\in\mathcal{M}\setminus\{k\},\label{eq:cSINR-1-2}\\
\bar{\wp}_{k\ell} & \geq\big|\Im\big\{\mathbf{h}_{k}\tilde{\mathbf{x}}_{\ell}\big\}\big|\ \forall\ell\in\mathcal{M}\setminus\{k\}.\label{eq:cSINR-1-3}
\end{align}

It is easy to see that if (<ref>) is feasible, then
so is (<ref>) and vice versa. Note that the right-hand
side (RHS) of (<ref>) is convex, and we only need to
obtain a concave lower bound on the left-hand side (LHS) of (<ref>).
Following a similar line of argument to (<ref>),
this can be done as follows:

\begin{align}
& \frac{1}{\Gamma_{k}}\big|\mathbf{h}_{k}\mathbf{x}_{k}\big|^{2} \geq\frac{1}{\Gamma_{k}}\Re\big\{\big(\mathbf{d}_{k}^{(i)}\!\!\ \herm\big)\big[c_{k}^{(i)}\mathbf{h}_{k}\herm+\mathbf{x}_{k}\big]\big\} -\frac{1}{2}\|\mathbf{d}_{k}^{(i)}\|^{2}\nonumber \\
& \!\!\!-\frac{1}{2}\|c_{k}^{(i)}\mathbf{h}_{k}\herm-\mathbf{x}_{k}\|^{2}-\big|c_{k}^{(i)}\big|^{2} \triangleq\frac{1}{\Gamma_{k}}\bar{f}_{k}\big(\mathbf{x}_{k},\boldsymbol{\theta};\mathbf{x}_{k}^{(i)},\boldsymbol{\theta}^{(i)}\big), \!\label{eq:fBarDef}
\end{align}

where $c_{k}^{(i)}\triangleq\mathbf{h}_{k}^{(i)}\mathbf{x}_{k}^{(i)}$
and $\mathbf{d}_{k}^{(i)}\triangleq c_{k}^{(i)}\mathbf{h}_{k}^{(i)}\!\!\ \herm+\mathbf{x}_{k}^{(i)}$.

Using the fact that $u\geq|v|$ iff $u\geq v$ or $u\geq|-v|$, and
following (<ref>), $\wp_{k\ell}$ in (<ref>)
can be equivalently written as

\begin{align}
\wp_{k\ell} & \geq\Re\big\{\mathbf{h}_{k}\tilde{\mathbf{x}}_{\ell}\big\}=\frac{1}{4}\big(\|\mathbf{h}_{k}\herm+\tilde{\mathbf{x}}_{\ell}\|^{2}-\|\mathbf{h}_{k}\herm-\tilde{\mathbf{x}}_{\ell}\|^{2}\big),\label{eq:tkl-1-1}\\
\wp_{k\ell} & \geq-\Re\big\{\mathbf{h}_{k}\tilde{\mathbf{x}}_{\ell}\big\}=\frac{1}{4}\big(\|\mathbf{h}_{k}\herm-\tilde{\mathbf{x}}_{\ell}\|^{2}-\|\mathbf{h}_{k}\herm+\tilde{\mathbf{x}}_{\ell}\|^{2}\big).\label{eq:tkl-1-2}
\end{align}

Since the negative quadratic term in the RHS of (<ref>)
results in its non-convexity, we use the inequality in (<ref>)
to convexify (<ref>) as follows:

\begin{align}
\wp_{k\ell} & \geq\frac{1}{4}\big[\|\mathbf{h}_{k}\herm+\tilde{\mathbf{x}}_{\ell}\|^{2}-2\Re\big\{\big(\mathbf{h}_{k}^{(i)}-\tilde{\mathbf{x}}_{\ell}^{(i)}\!\!\ \herm\big)\big(\mathbf{h}_{k}\herm-\tilde{\mathbf{x}}_{\ell}\big)\big\}\nonumber \\
& \quad+\|\mathbf{h}_{k}^{(i)}\!\!\ \herm-\tilde{\mathbf{x}}_{\ell}^{(i)}\|^{2}\big]\triangleq\mu_{k\ell}\big(\tilde{\mathbf{x}}_{\ell},\boldsymbol{\theta};\tilde{\mathbf{x}}_{\ell}^{(i)},\boldsymbol{\theta}^{(i)}\big).\label{eq:muklDef}
\end{align}

Following a similar argument, (<ref>) yields

\begin{align}
\wp_{k\ell} & \geq\frac{1}{4}\big[\|\mathbf{h}_{k}\herm-\tilde{\mathbf{x}}_{\ell}\|^{2}-2\Re\big\{\big(\mathbf{h}_{k}^{(i)}+\tilde{\mathbf{x}}_{\ell}^{(i)}\!\!\ \herm\big)\big(\mathbf{h}_{k}\herm+\tilde{\mathbf{x}}_{\ell}\big)\big\}\nonumber \\
& \quad+\|\mathbf{h}_{k}^{(i)}\!\!\ \herm+\tilde{\mathbf{x}}_{\ell}^{(i)}\|^{2}\big]\triangleq\bar{\mu}_{k\ell}\big(\tilde{\mathbf{x}}_{\ell},\boldsymbol{\theta};\tilde{\mathbf{x}}_{\ell}^{(i)},\boldsymbol{\theta}^{(n)}\big).\label{eq:muklHatDef}
\end{align}

Analogously, (<ref>) yields the following inequalities:

\begin{align}
\bar{\wp}_{k\ell} & \geq\frac{1}{4}\big[\|\mathbf{h}_{k}\herm-j\tilde{\mathbf{x}}_{\ell}\|^{2}-2\Re\big\{\big(\mathbf{h}_{k}^{(i)}-j\tilde{\mathbf{x}}_{\ell}^{(i)}\!\!\ \herm\big)\big(\mathbf{h}_{k}\herm+j\tilde{\mathbf{x}}_{\ell}\big)\big\}\nonumber \\
& \quad+\|\mathbf{h}_{k}^{(i)}\!\!\ \herm+j\tilde{\mathbf{x}}_{\ell}^{(i)}\|^{2}\big]\triangleq\upsilon_{k\ell}\big(\tilde{\mathbf{x}}_{\ell},\boldsymbol{\theta};\tilde{\mathbf{x}}_{\ell}^{(i)},\boldsymbol{\theta}^{(i)}\big),\label{eq:uklDef}\\
\bar{\wp}_{k\ell} & \geq\frac{1}{4}\big[\|\mathbf{h}_{k}\herm+j\tilde{\mathbf{x}}_{\ell}\|^{2}-2\Re\big\{\big(\mathbf{h}_{k}^{(i)}+j\tilde{\mathbf{x}}_{\ell}^{(i)}\!\!\ \herm\big)\big(\mathbf{h}_{k}\herm-j\tilde{\mathbf{x}}_{\ell}\big)\big\}\nonumber \\
& \quad+\|\mathbf{h}_{k}^{(i)}\!\!\ \herm-j\tilde{\mathbf{x}}_{\ell}^{(i)}\|^{2}\big]\triangleq\bar{\upsilon}_{k\ell}\big(\tilde{\mathbf{x}}_{\ell},\boldsymbol{\theta};\tilde{\mathbf{x}}_{\ell}^{(i)},\boldsymbol{\theta}^{(i)}\big).\label{eq:uklBarDef}
\end{align}

We now focus on the non-convex constraint in (<ref>),
which for any $k\in\mathcal{K}$, can be written as

\begin{equation}
\hat{\gamma}_{k}\leq\hat{\Gamma}_{k}\Rightarrow\sigma_{\text{T}}^{2}+\sum_{\ell\in\mathcal{M}\setminus\{k\}}|\mathbf{g}\tilde{\mathbf{x}}_{\ell}|^{2}\geq\frac{1}{\hat{\Gamma}}|\mathbf{g}\mathbf{x}_{k}|^{2}.\label{eq:SINRC-Tk-1}
\end{equation}

Note that we need a concave lower bound on the LHS of (<ref>),
and a convex upper bound on the RHS. Similar to (<ref>),
the former can be obtained by linearizing the quadratic term in the
LHS as follows:

\begin{align}
& \sigma_{\text{T}}^{2}+\sum_{\ell\in\mathcal{M}\setminus\{k\}}|\mathbf{g}\tilde{\mathbf{x}}_{\ell}|^{2} \geq \sigma_{\text{T}}^{2}+\sum_{\ell\in\mathcal{M}\setminus\{k\}}f_{\ell}(\tilde{\mathbf{x}}_{\ell},\boldsymbol{\theta};\tilde{\mathbf{x}}_{\ell}^{(i)},\boldsymbol{\theta}^{(i)}).\label{eq:SINRC-Tk-LHS}
\end{align}

On the other hand, a convex upper bound on $|\mathbf{g}\mathbf{x}_{k}|^{2}/\hat{\Gamma}$
is given by $(\tau_{k}^{2}+\bar{\tau}_{k}^{2})/\hat{\Gamma}$, where
$\tau_{k}\geq|\Re\{\mathbf{g}\mathbf{x}_{k}\}|$ and $\bar{\tau}_{k}\geq|\Im\{\mathbf{g}\mathbf{x}_{k}\}|$.
Therefore, using (<ref>), (<ref>)
and the preceding arguments, for a given $k\in\mathcal{K}$, the constraint
in (<ref>) can be equivalently written as

\begin{align}
\sigma_{\mathrm{T}}^{2}+\sum_{\ell\in\mathcal{M}\setminus\{k\}}f_{\ell}\big(\tilde{\mathbf{x}}_{\ell},\boldsymbol{\theta};\tilde{\mathbf{x}}_{\ell}^{(i)},\boldsymbol{\theta}^{(i)}\big) & \geq\frac{1}{\hat{\Gamma}_{k}}\big(\tau_{k}^{2}+\bar{\tau}_{k}^{2}\big),\label{eq:secureSINRC}\\
\tau_{k} & \geq\big|\Re\big\{\mathbf{g}\mathbf{x}_{k}\big\}\big|,\label{eq:secureSINRC-real}\\
\bar{\tau}_{k} & \geq\big|\Im\big\{\mathbf{g}\mathbf{x}_{k}\big\}\big|.\label{eq:secureSINRC-imag}
\end{align}

Again, it can be noted that if (<ref>) is feasible,
then so is (<ref>) and vice versa. Moreover,
following a similar set of arguments to those in (<ref>)–(<ref>),
lower bounds on $\tau_{k}$ and $\bar{\tau}_{k}$ in (<ref>)
and (<ref>), respectively, are given by

\begin{align}
\tau_{k} & \geq\frac{1}{4}\big[\|\mathbf{g}\herm+\mathbf{x}_{k}\|^{2}-2\Re\big\{\big(\mathbf{g}^{(i)}-\mathbf{x}_{k}^{(i)}\!\!\ \herm\big)\big(\mathbf{g}\herm-\mathbf{x}_{k}\big)\big\}\nonumber \\
& \quad+\|\mathbf{g}^{(i)}\!\!\ \herm-\mathbf{x}_{k}^{(i)}\|^{2}\big]\triangleq\eta_{k}\big(\mathbf{x}_{k},\boldsymbol{\theta};\mathbf{x}_{k}^{(i)},\boldsymbol{\theta}^{(i)}\big),\label{eq:etaDef}\\
\tau_{k} & \geq\frac{1}{4}\big[\|\mathbf{g}\herm-\mathbf{x}_{k}\|^{2}-2\Re\big\{\big(\mathbf{g}^{(i)}+\mathbf{x}_{k}^{(i)}\!\!\ \herm\big)\big(\mathbf{g}\herm+\mathbf{x}_{k}\big)\big\}\nonumber \\
& \quad+\|\mathbf{g}^{(i)}\!\!\ \herm+\mathbf{x}_{k}^{(i)}\|^{2}\big]\triangleq\bar{\eta}_{k}\big(\mathbf{x}_{k},\boldsymbol{\theta};\mathbf{x}_{k}^{(i)},\boldsymbol{\theta}^{(i)}\big),\label{eq:etaBarDef}\\
\bar{\tau}_{k} & \geq\frac{1}{4}\big[\|\mathbf{g}\herm-j\mathbf{x}_{k}\|^{2}-2\Re\big\{\big(\mathbf{g}^{(i)}-j\mathbf{x}_{k}^{(i)}\!\!\ \herm\big)\big(\mathbf{g}\herm+j\mathbf{x}_{k}\big)\big\}\nonumber \\
& \quad+\|\mathbf{g}^{(i)}\!\!\ \herm+j\mathbf{x}_{k}^{(i)}\|^{2}\big]\triangleq\chi_{k}\big(\mathbf{x}_{k},\boldsymbol{\theta};\mathbf{x}_{k}^{(i)},\boldsymbol{\theta}^{(i)}\big),\label{eq:chiDef}\\
\bar{\tau}_{k} & \geq\frac{1}{4}\big[\|\mathbf{g}\herm+j\mathbf{x}_{k}\|^{2}-2\Re\big\{\big(\mathbf{g}^{(i)}+j\mathbf{x}_{k}^{(i)}\!\!\ \herm\big)\big(\mathbf{g}\herm-j\mathbf{x}_{k}\big)\big\}\nonumber \\
& \quad+\|\mathbf{g}^{(i)}\!\!\ \herm-j\mathbf{x}_{k}^{(i)}\|^{2}\big]\triangleq\bar{\chi}_{k}\big(\mathbf{x}_{k},\boldsymbol{\theta};\mathbf{x}_{k}^{(i)},\boldsymbol{\theta}^{(i)}\big).\label{eq:chiBarDef}
\end{align}

Next, since the constraint in (<ref>) is already convex,
we are left only with the non-convexity of (<ref>). To tackle
this, we first relax the equality constraint in (<ref>) by
a (convex) inequality constraint. In order to ensure that the inequality
constraint is satisfied with equality (i.e., the constraint is binding
at convergence), we add a regularization term in the objective and
handle the resulting non-convex objective by the first-order approximation
of the regularization term around $\boldsymbol{\theta}^{(i)}$. Therefore,
an equivalent reformulation of the problem in (<ref>)
can be given by

\begin{align}
\underset{\mathbf{X},\boldsymbol{\theta},\boldsymbol{\wp},\bar{\boldsymbol{\wp}},\boldsymbol{\tau},\bar{\boldsymbol{\tau}}}{\maximize}\  & \sum_{m\in\mathcal{M}}f_{m}\big(\tilde{\mathbf{x}}_{m},\boldsymbol{\theta};\tilde{\mathbf{x}}_{m}^{(i)},\boldsymbol{\theta}^{(i)}\big)\nonumber \\
& \qquad+\zeta\big[2\Re\big\{\boldsymbol{\theta}^{(n)}\ \herm\boldsymbol{\theta}\big\}-\|\boldsymbol{\theta}^{(n)}\|^{2}\big],\label{eq:transObj-1}\\
\st\  & \frac{1}{\Gamma_{k}}\bar{f}_{k}\big(\mathbf{x}_{k},\boldsymbol{\theta};\mathbf{x}_{k}^{(i)},\boldsymbol{\theta}^{(i)}\big)\nonumber \\
& \quad\geq\sigma_{k}^{2}+\!\!\!\!\sum_{\ell\in\mathcal{M}\setminus\{k\}}\!\!\!\!\big(\wp_{k\ell}^{2}+\bar{\wp}_{k\ell}^{2}\big)\ \forall k\in\mathcal{K},\label{eq:transSINRC-1}\\
& \eqref{eq:muklDef}-\eqref{eq:uklBarDef}\ \forall k\in\mathcal{K},\forall\ell\in\mathcal{M}\setminus\{k\},\label{eq:transSINRC-2}\\
& \eqref{eq:secureSINRC},\eqref{eq:etaChi}\ \forall k\in\mathcal{K},\nonumber \\
& \eqref{eq:TPC},\nonumber \\
& |\theta_{n}|\leq1\ \forall n\in\mathcal{N},\label{eq:relaxedUMC}
\end{align}

where $\boldsymbol{\wp}\triangleq[\wp_{11},\wp_{12},\ldots,\wp_{KL}]\trans$,
and $\zeta>0$ is the regularization parameter. It is straightforward
to show that all of the constraints in (<ref>) can
be represented by quadratic cones, and therefore (<ref>)
is an SOCP problem which can be solved efficiently using off-the-shelf
solvers, e.g., MOSEK [18]. The proposed SCA-based SOCP method
is outlined in Algorithm <ref>.
Proposed SCA-based Method to Solve (<ref>).

$\mathbf{X}^{(0)}$, $\boldsymbol{\theta}^{(0)}$, $\xi>0$

Solve (<ref>) and denote the solution as $\mathbf{X}^{\star}$,

Update: $\mathbf{X}^{(i+1)}\leftarrow\mathbf{X}^{\star}$, $\boldsymbol{\theta}^{(i+1)}\leftarrow\boldsymbol{\theta}^{\star}$

$i\leftarrow i+1$

$\mathbf{X}^{\star}$, $\boldsymbol{\theta}^{\star}$

One needs to find feasible starting points $\mathbf{X}^{(0)}$ and
$\boldsymbol{\theta}^{(0)}$ to run  Algorithm <ref>,
which is not straightforward. Therefore, below we describe a practical
way to obtain a set of initial points. Consider the following optimization

\begin{align}
\underset{\mathbf{X},\boldsymbol{\theta},\boldsymbol{\delta},\bar{\boldsymbol{\delta}}}{\minimize}\  & \sum_{k\in\mathcal{K}}(\delta_{k}+\bar{\delta}_{k}),\label{eq:initialPointObj}\\
\st\  & \delta_{k}+\frac{1}{\Gamma_{k}}\bar{f}_{k}\big(\mathbf{x}_{k},\boldsymbol{\theta};\mathbf{x}_{k}^{(i)},\boldsymbol{\theta}^{(i)}\big)\nonumber \\
& {\small \quad\geq\sigma_{k}^{2}+ \sum \nolimits_{\ell\in\mathcal{M}\setminus\{k\}}\!\!\!\!\big(\wp_{k\ell}^{2}+\bar{\wp}_{k\ell}^{2}\big)\ \forall k\in\mathcal{K},\label{eq:intialC1}}\\
\  & \bar{\delta}_{k}+\sigma_{\mathrm{T}}^{2}+\sum_{\ell\in\mathcal{M}\setminus\{k\}}f_{\ell}\big(\tilde{\mathbf{x}}_{\ell},\boldsymbol{\theta};\tilde{\mathbf{x}}_{\ell}^{(i)},\boldsymbol{\theta}^{(i)}\big)\nonumber \\
& \qquad\geq\frac{1}{\hat{\Gamma}_{k}}\big(\tau_{k}^{2}+\bar{\tau}_{k}^{2}\big)\ \forall k\in\mathcal{K},\label{eq:initialC2}\\
& \eqref{eq:TPC},\eqref{eq:UMC},\eqref{eq:etaChi},\eqref{eq:transSINRC-2},\nonumber \\
& \delta_{k}\geq0,\bar{\delta}_{k}\geq0\ \forall k\in\mathcal{K}.\label{eq:initialC3}
\end{align}

Note that the problem in (<ref>) is always feasible
for sufficiently large $\boldsymbol{\delta}$ and $\bar{\boldsymbol{\delta}}$.
We solve the problem in (<ref>) by following
a similar procedure to that of Algorithm <ref>, with
random $\mathbf{X}$ and $\boldsymbol{\theta}$ as initial points.
The minimization in <ref> forces $\delta_{k}$
and $\bar{\delta}_{k}$ to approach $0$. At convergence, if $\delta_{k}=\bar{\delta}_{k}=0\ \forall k\in\mathcal{K}$,
the problem in (<ref>) is obviously feasible. Thus
we can choose the final values of $\mathbf{X}$ and $\boldsymbol{\theta}$
in (<ref>) as initial points for Algorithm <ref>.
However, if the objective $\sum_{k\in\mathcal{K}}(\delta_{k}+\bar{\delta}_{k})$
is not zero at convergence, then we simply declare that the considered
problem is infeasible and will not run Algorithm <ref>[We note that the considered problem may be feasible even though $\sum_{k\in\mathcal{K}}(\delta_{k}+\bar{\delta}_{k})>0$.
The reason is that the SCA-based method applied to solve (<ref>)
can only guarantee a stationary solution. In general, checking (<ref>)
is feasible or not is an NP hard problem since the feasible set is
non-convex. For practical purposes, if the SCA-based method cannot
find a feasible solution, we can simply say that the problem is infeasible
and ignore this realization.]

The convergence of the proposed SCA-based method in Algorithm <ref>
can be readily proved following the set of arguments in <cit.>.

§.§ Complexity Analysis

It is straightforward to show that the total number of (real-valued)
optimization variables in (<ref>) is $2\big(L^{2}+K^{2}+2KL+N\big)+1$,
and the total number of (second-order) conic constraints is $4K^{2}+4KL+2K+N+2$.
Therefore, following the arguments in <cit.>,
the overall per-iteration complexity of the proposed SOCP-based method
is given by

\begin{align}
& \mathcal{O}\big[\big(4K^{2}+4KL+N\big)^{0.5}\big(2K^{2}+4KL+2L^{2}+2N)\nonumber \\
& \big(4K^{5}+8K^{4}L+4K^{3}L^{2}+48K^{3}L+60K^{2}L^{2}+24KL^{3}\nonumber \\
& +52KL^{2}+4L^{4}+\big(2K^{2}+4KL+2L^{2}+2N\big)^{2}\big)\big].\label{eq:exactComplexity}
\end{align}

However, in a practical setup, the number of elements in the IRS is
expected to be much larger than the number of BS antennas and the
number of users, i.e., $N\gg\max\{L,K\}$. Hence, the complexity of
the proposed SCA-based method can be well-approximated by $\mathcal{O}\big(N^{3.5}\big)$.
On the other hand, the per-iteration computational complexity of <cit.>
can be approximated by $\mathcal{O}(N^{3})$ (see <cit.>).
Although the order of complexity of the proposed SCA-based method
is slightly higher than that of the penalty-based benchmark, we will
show in the simulation section that the proposed method requires fewer
iterations, resulting in a significantly reduced problem-solving time.

Convergence results for $L=4$, $K=3$, $N=100$, $P=40$ dBm, $\Gamma_{k}=10$ dB $\forall k\in\mathcal{K}$ and $\hat{\Gamma}_{k}=0$ dB $\forall k\in\mathcal{K}$.

Average beampattern gain for $L=4$, $K=3$,$\Gamma_{k}=10$ dB $\forall k\in\mathcal{K}$ and $\hat{\Gamma}_{k}=0$ dB $\forall k\in\mathcal{K}$.

Average problem solving time for $L=4$, $P=40$ dBm$,\Gamma_{k}=10$ dB $\forall k\in\mathcal{K}$ and $\hat{\Gamma}_{k}=0$ dB $\forall k\in\mathcal{K}$.

§ NUMERICAL RESULTS AND DISCUSSION

In this section, we present a detailed performance comparison between
the proposed SCA-based method and the penalty-based benchmark approach
of <cit.>. The location of the nodes
and the channel model assumed here are the same as those in [14].
The simulations are performed on a high-performance computing cluster
with a Intel Xeon Gold 6152 processor, using Python v3.9.7 and MOSEK
Fusion API for Python Rel.-10.0.40 [18]. In Figs. <ref>
and <ref>, the results are obtained by averaging over
100 independent channel realizations.

In Fig. <ref>, we show the convergence behavior of both
the proposed and penalty-based benchmark methods. For the given set
of channels, the proposed SCA-based method converges in less than
30 iterations, whereas the penalty-based benchmark requires around
270 iterations. Nevertheless, even with significantly fewer iterations,
the proposed method results in nearly a 30% higher beampattern gain
as compared to that offered by the penalty-based benchmark. More interestingly,
each iteration of the proposed SCA-based method returns a set of feasible
points and therefore the iterations of the proposed method can be
terminated even before convergence has been attained, if this is required.
On the other hand, the penalty-based benchmark returns a feasible
solution only in the final outer-loop iteration, and therefore the
algorithm cannot be stopped earlier to achieve a feasible solution.
Therefore, the benchmark is not suitable in rapidly changing environments
with very small coherence times where at least a suboptimal solution
is required within a certain fraction of the channel's coherence time.

The impact of the number of IRS elements on the average beampattern
gain for the two algorithms- is shown in Fig. <ref>. An
increase in the number of IRS elements increases the degrees-of-freedom
at the IRS, allowing the IRS to perform highly focused beamforming.
This in turn results in increasing beampattern gain with increasing
$N$. On the other hand, since a fixed amount of transmit power is
required to achieve the SINR constraints at the communication users,
a higher transmit power budget results in a higher surplus power at
the BS, which is then used to attain a larger beampattern gain toward
the target. Therefore, increasing the value of $P$ increases the
average beampattern gain, which is also clearly evident from the figure.
The performance gap between the SCA-based and penalty-based methods
increases with an increase in the number of IRS elements. As the value
of $N$ increases, the impact of coupling between $\mathbf{X}$ and
$\boldsymbol{\theta}$ becomes more intricate. Therefore, the solution
obtained via the AO-based approach of [14] returns
a highly suboptimal beampattern gain. On the other hand, as clearly
observed in the figure, the simultaneous update of all variables in
the proposed algorithm outperforms the penalty-based benchmark.

In Fig. <ref>, we plot the average problem-solving time
versus the number of IRS elements for different numbers of communication
users $K$. As the value of $N$ and/or $K$ increases, the size of
the optimization problem to be solved also increases for both methods.
This in turn increases the average problem solving time for both approaches.
Although the per-iteration complexity of the proposed SCA-based method
is slightly higher than that of the penalty-based benchmark, the proposed
approach requires a much smaller time to find the solution due to
its convergence in fewer iterations.

§ CONCLUSION

In this paper, we have considered the problem of optimal transmit
and reflective beamforming design in a secure IRS-enabled ISAC system.
More specifically, we aim to maximize the beampattern gain toward
the eavesdropping target subject to the SINR constraints at the communication
users and information leakage constraints at the target. In contrast
to the conventional AO-based approach, we proposed a novel SCA-based
optimization in which all variables are updated simultaneously. The
superiority of the proposed method was clearly established with the
help of numerical experiments in terms of both achieving a high-performance
solution and low problem-solving time compared to that of the penalty-based
benchmark. Moreover, the performance gap between the proposed SCA-based
approach and penalty-based benchmark was shown to be increasing with
the number of IRS elements or the transmit power budget.

[1]
Z. Wei, F. Liu, C. Masouros, N. Su, and A. P. Petropulu, “Toward
multi-functional 6G wireless networks: Integrating sensing,
communication, and security,” IEEE Commun. Mag., vol. 60, no. 4, pp.
65–71, 2022.

[2]
D. Ma, N. Shlezinger, T. Huang, Y. Liu, and Y. C. Eldar, “Joint
radar-communication strategies for autonomous vehicles: Combining two key
automotive technologies,” IEEE Signal Process. Mag., vol. 37, no. 4,
pp. 85–97, 2020.

[3]
X. Fang, W. Feng, Y. Chen, N. Ge, and Y. Zhang, “Joint communication and
sensing toward 6G: Models and potential of using MIMO,” IEEE
Internet Things J., vol. 10, no. 5, pp. 4093–4116, 2023.

[4]
Z. Chen et al., “Reconfigurable-intelligent-surface-assisted B5G/6G
wireless communications: Challenges, solution, and future opportunities,”
IEEE Commun. Mag., vol. 61, no. 1, pp. 16–22, 2023.

[5]
A. Mukherjee, V. Kumar, and L.-N. Tran, “Secrecy rate maximization for
intelligent reflecting surface assisted MIMOME wiretap channels,” in
IEEE MILCOM, 2021, pp. 261–266.

[6]
A. Mukherjee, V. Kumar, D. W. K. Ng, and L.-N. Tran, “On the energy-efficiency
maximization for IRS-assisted MIMOME wiretap channels,” in IEEE
VTC (Fall), 2022, pp. 1–6.

[7]
H. Du, J. Kang, D. Niyato, J. Zhang, and D. I. Kim, “Reconfigurable
intelligent surface-aided joint radar and covert communications:
Fundamentals, optimization, and challenges,” IEEE Veh. Technol.
Mag., vol. 17, no. 3, pp. 54–64, 2022.

[8]
R. Liu et al., “Integrated sensing and communication with
reconfigurable intelligent surfaces: Opportunities, applications, and
future directions,” IEEE Wireless Commun., vol. 30, no. 1, pp.
50–57, 2023.

[9]
X. Song et al., “Joint transmit and reflective beamforming for
IRS-assisted integrated sensing and communication,” in IEEE WCNC,
2022, pp. 189–194.

[10]
H. Luo, R. Liu, M. Li, Y. Liu, and Q. Liu, “Joint beamforming design for
RIS-assisted integrated sensing and communication systems,” IEEE
Trans. Veh. Technol., vol. 71, no. 12, pp. 13 393–13 397, 2022.

[11]
R. Liu, M. Li, and A. L. Swindlehurst, “Joint beamforming and reflection
design for RIS-assisted ISAC systems,” in EUSIPCO, 2022, pp.

[12]
R. P. Sankar and S. P. Chepuri, “Beamforming in hybrid RIS assisted
integrated sensing and communication systems,” in EUSIPCO, 2022, pp.

[13]
A. A. Salem, M. H. Ismail, and A. S. Ibrahim, “Active reconfigurable
intelligent surface-assisted MISO integrated sensing and communication
systems for secure operation,” IEEE Trans. Veh. Technol., to appear.

[14]
M. Hua, Q. Wu, W. Chen, O. A. Dobre, and A. L. Swindlehurst, “Secure
intelligent reflecting surface aided integrated sensing and communication,”
arXiv preprint arXiv:2207.09095, 2022.

[15]
B. Feng et al., “Optimization techniques in reconfigurable intelligent
surface aided networks,” IEEE Wireless Commun., vol. 28, no. 6, pp.
87–93, 2021.

[16]
X. Yu et al., “Robust and secure wireless communications via
intelligent reflecting surfaces,” IEEE J. Sel. Areas Commun.,
vol. 38, no. 11, pp. 2637–2652, 2020.

[17]
V. Kumar, R. Zhang, M. D. Renzo, and L.-N. Tran, “A novel SCA-based method
for beamforming optimization in IRS/RIS-assisted MU-MISO downlink,”
IEEE Wireless Commun. Lett., vol. 12, no. 2, pp. 297–301, 2023.

[18]
MOSEK ApS, MOSEK Fusion API for Python. Release 10.0.40, 2023.
[Online]. Available: <http://docs.mosek.com/10.0/pythonfusion.pdf>

[19]
A. Ben-Tal and A. Nemirovski, Lectures on modern convex
optimization.1em plus 0.5em minus 0.4emPhiladelphia: SIAM,


[a]Pedro M. F. Pereira

# Non-Unitary Mixing Matrices in Neutrino and Vector-like Quark Models

###### Abstract

This work analyses a class of extensions of the Standard Model through the
lens of a novel exact parameterisation. This parameterisation is specially
suitable to the analysis of extensions of the Standard Model with non-unitary
mixing matrices, like models with vector-like fermions or right-handed
neutrinos. The usefulness of this parameterisation is motivated with two
example models: A) Standard Model with the addition of $n_{u}$ up and $n_{d}$
down singlet vector-like quarks and B) Standard Model with the addition of
$n_{R}$ right-handed neutrinos.

##  1 Introduction

This paper is based on refs. [3, 1, 2] as well as work in progress done in
collaboration with J.M. Alves, G.C. Branco, A.L. Cherchiglia, C.C. Nishi, J.T.
Penedo, M.N. Rebelo, and J.I. Silva-Marcos, to whom I am deeply grateful.
Studying a class of extensions of the Standard Model with non-unitary mixing
matrices can easily become a very demanding task when the region of interest
of the parameter space is very close to regions where common approximations,
widely used in the literature, fail. From this problem originated the need to
develop an exact parameterisation, that would work in all cases.
For instance, the classic neutrino seesaw mechanism papers [4, 5] work under
the assumption that the deviations from unitarity of the leptonic mixing
matrix are negligible, which imply that the heavy neutrino masses are around
the GUT scale, in order to have order one Yukawa couplings in the neutrino
Dirac mass matrix. Furthermore, the famous Casas-Ibarra parameterisation [6]
is also an approximation, if used naively, as explained in section 2.2 of ref.
[1]. Therefore, if one wants to study a model where some of the heavy
neutrinos have masses around the GeV scale or below, the aforementioned
assumptions start to crumble. Recently, a parameterisation where the mixing
matrix is parameterised via a power series was introduced [7]. The small
parameter of the power series is the deviations from unitarity of the mixing
matrix and a first order truncation on this parameter yields an accurate
result for models with GeV neutrinos. However, if one wants to study a model
where the spectra of heavy neutrinos includes eV, keV or MeV masses, none of
the parameterisations available in the literature was suitable for it.
In the case of vector-like quarks, many works in the literature use the
assumption that vector-like quarks only couple to the third generation of SM
quarks [8]. This does not need to be the case, as experimental evidence
signals that if vector-like quarks are the explanation for the CKM unitarity
problem, then the vector like quark might couple more to the 1st generation,
as discussed in ref. [2].
Hence, the following exact parameterisation is suitable when one wants to
study a region of the parameter space where common approximations fail or when
one wants to perform a general scan of the parameter space, without biases.

##  2 Mass Matrices and their Diagonalisation

The usefulness of the parameterisation will be motivated with two example
models:

  1. (A)

SM with the addition of $n_{u}$ up and $n_{d}$ down singlet vector-like quarks

  2. (B)

SM with the addition of $n_{R}$ right-handed neutrinos.

The procedure to derive the parameterisation will be performed in the next
sections. Keep in mind that no approximations will be performed in any step,
the results are exact. The first step is to define the mass matrices and then
diagonalise them. Considering the mass Lagrangian of (A)

$\mathcal{L}_{\text{m}}^{q}\,=\,-\begin{pmatrix}\overline{d}^{0}_{L}&\overline{D}^{0}_{L}\end{pmatrix}\,\mathcal{M}_{d}\,\begin{pmatrix}d^{0}_{R}\\\\[2.84526pt]
D^{0}_{R}\end{pmatrix}-\begin{pmatrix}\overline{u}^{0}_{L}&\overline{U}^{0}_{L}\end{pmatrix}\,\mathcal{M}_{u}\,\begin{pmatrix}u^{0}_{R}\\\\[2.84526pt]
U^{0}_{R}\end{pmatrix}+\text{h.c.}\,,$ (1)

and of (B), in the weak basis (WB) where charged leptons are diagonal,

$\displaystyle\mathcal{L}_{\text{m}}^{\nu}$
$\displaystyle=\,-\left[\frac{1}{2}n_{L}^{T}C^{*}\mathcal{M}^{*}n_{L}+\overline{l_{L}}d_{l}l_{R}\right]+\text{h.c.}\,,$
(2)

where $n_{L}^{T}=(\nu_{L}^{T}\,\,\,\overline{\nu_{R}}\,C^{T})$,
$n_{L}=(\nu_{L}\,,\,\,C\,\overline{\nu_{R}}^{T})$ is a column vector and $C$
is the charge conjugation matrix. The up or down ($q=u,d$) quark
$(3+n_{q})\times(3+n_{q})$ mass matrix is given by

$\displaystyle\mathcal{M}_{q}\,=\,\left(\begin{array}[]{c;{2pt/2pt}c}\,m_{q}&\,\overline{m}\hskip
0.86108pt_q\\\\[4.2679pt]
\hdashline[2pt/2pt]\\!{\,\,\,\overline{\\!M}\mkern-2.0mu_{q}\,\,\,}\\!&\\!{\,\,\,M_q\,\,\,}\\!\\\\[2.84526pt]
\end{array}\right)\,,$ (5)

and the $(3+n_{R})\times(3+n_{R})$ neutrino mass matrix is given by

$\displaystyle\mathcal{M}\,=\,\left(\begin{array}[]{c;{2pt/2pt}c}\,\textbf{0}&\,m\\\\[4.2679pt]
\hdashline[2pt/2pt]{\,\,\,m^{T}\,\,\,}\\!&{\,\,\,M\,\,\,}\\!\\\\[2.84526pt]
\end{array}\right)\,,$ (8)

the diagonalisation equations are

${\mathcal{V}^{q}_{L}}^{\dagger}\,\,\mathcal{M}_{q}\,\,\mathcal{V}_{R}^{q}\,=\,\mathcal{D}_{q}\,\,~{}~{},~{}~{}\mathcal{V}^{T}\mathcal{M}^{*}\mathcal{V}=\mathcal{D}~{}.$
(9)

The matrices that diagonalise the mass matrices, $\mathcal{V}$, are
$(3+n)\times(3+n)$ 111$n=n_{q}$ for quarks and $n=n_{R}$ for neutrinos unitary
matrices and can be decomposed into two other matrices $A$ and $B$

$\mathcal{V}^{q}_{\chi}\,=\,\left(\begin{array}[]{c}{}\hfil\\\\[-11.38109pt]
\qquad A^{q}_{\chi}\\\\[5.69054pt] \hdashline[2pt/2pt]{}\\\\[-11.38109pt]
\qquad B^{q}_{\chi}\\\\[5.69054pt] \end{array}\right)\,,$ (10)

$\mathcal{V}\,=\,\left(\begin{array}[]{c}{}\hfil\\\\[-11.38109pt] \qquad
A\\\\[5.69054pt] \hdashline[2pt/2pt]{}\\\\[-11.38109pt] \qquad
B\\\\[5.69054pt] \end{array}\right)\,,$ (11)

where $\chi=L,R$, $A$ is a $3\times(3+n)$ and $B$ a $n\times(3+n)$ matrix. $A$
and $B$ are non-unitary matrices, in general. After this decomposition into
$A$ and $B$, one can rewrite the equations in eq. (9) as

$\begin{split}&m_{q}=A^{q}_{L}\,\mathcal{D}_{q}\,{A_{R}^{q}}^{\dagger}\,,\quad\\\
&\overline{m}\hskip
0.86108pt_{q}=A^{q}_{L}\,\mathcal{D}_{q}\,{B_{R}^{q}}^{\dagger}\,,\quad\\\
&\overline{\\!M}\mkern-2.0mu_{q}=B^{q}_{L}\,\mathcal{D}_{q}\,{A_{R}^{q}}^{\dagger}\,,\quad\\\
&M_{q}=B^{q}_{L}\,\mathcal{D}_{q}\,{B_{R}^{q}}^{\dagger}\,,\end{split}$ (12)

for quarks, and

$\begin{split}&\textbf{0}=A\,\mathcal{D}\,A^{T}\,,\quad\\\
&m=A\,\mathcal{D}\,B^{T}\,,\quad\\\ &M=B\,\mathcal{D}\,B^{T}\,,\end{split}$
(13)

for neutrinos.
The last thing to do is a proper parameterisation of these $A$ and $B$
matrices. Before proceeding, the next section contains, for completeness, the
charged currents and neutral interactions Lagrangians for quarks and
neutrinos, in the physical basis. Naturally, these will only depend on the $A$
part of the unitary matrices $\mathcal{V}$, as this contains all the physical
parameters of the theory, as it will become clear in the Parameterisation
section.

##  3 Interactions

The charged current Lagrangian for the quarks is given by

$\mathcal{L}_{W}^{q}=-\frac{g}{\sqrt{2}}\begin{pmatrix}\overline{u}_{L}&\overline{U}_{L}\end{pmatrix}\,V\,\gamma^{\mu}\begin{pmatrix}d_{L}\\\\[5.69054pt]
D_{L}\end{pmatrix}W_{\mu}^{+}\,+\,\text{h.c.}~{},$ (14)

where

$V=A_{L}^{u^{\dagger}}A_{L}^{d}~{},$ (15)

and for the neutrinos

$\mathcal{L}_{W}^{\nu}=-\frac{g}{\sqrt{2}}\overline{l}_{L}\,V\,\gamma^{\mu}\begin{pmatrix}n_{L}\\\\[5.69054pt]
N_{L}\end{pmatrix}W_{\mu}^{+}\,+\,\text{h.c.}~{},$ (16)

where

$V=A~{}.$ (17)

The part of the Lagrangian that contains Z interactions for the quarks is
given by

$\begin{split}\mathcal{L}_{Z}^{q}&=-\frac{g}{2\cos{\theta_{W}}}Z_{\mu}\left[\begin{pmatrix}\overline{u}_{L}&\overline{U}_{L}\end{pmatrix}\,F^{u}\,\gamma^{\mu}\begin{pmatrix}u_{L}\\\\[5.69054pt]
U_{L}\end{pmatrix}-\begin{pmatrix}\overline{d}_{L}&\overline{D}_{L}\end{pmatrix}\,F^{d}\,\gamma^{\mu}\begin{pmatrix}d_{L}\\\\[5.69054pt]
D_{L}\end{pmatrix}-2\sin^{2}{\theta_{W}}J^{\mu~{}q}_{\text{em}}\right]\,,\end{split}$
(18)

where

$\displaystyle J^{\mu~{}q}_{\text{em}}$
$\displaystyle=\,\frac{2}{3}\Big{(}\overline{u_{i}}\gamma^{\mu}{u_{i}}+\overline{U_{r}}\gamma^{\mu}{U_{r}}\Big{)}-\frac{1}{3}\Big{(}\overline{d_{i}}\gamma^{\mu}{d_{i}}+\overline{D_{s}}\gamma^{\mu}{D_{s}}\Big{)}\,,$
(19)

and $\psi\equiv\psi_{L}+\psi_{R}$ for $\psi\in\\{u_{i},d_{i},U_{r},D_{s}\\}$,
the index $i$ runs from 1 to 3 as in the SM, while the indices $r,s$ run from
1 to $n_{u,d}$. While the part that contains Higgs interactions is given by

$\mathcal{L}_{H}^{q}\,=\,-\frac{h}{v}\bigg{[}\begin{pmatrix}\overline{u}_{L}&\overline{U}_{L}\end{pmatrix}\,F^{u}\,\mathcal{D}_{u}\begin{pmatrix}u_{R}\\\\[5.69054pt]
U_{R}\end{pmatrix}+\begin{pmatrix}\overline{d}_{L}&\overline{D}_{L}\end{pmatrix}\,F^{d}\,\mathcal{D}_{d}\begin{pmatrix}d_{R}\\\\[5.69054pt]
D_{R}\end{pmatrix}\bigg{]}+\,\text{h.c.}\,~{},$ (20)

where

$F^{q}=A_{L}^{q^{\dagger}}A_{L}^{q}~{}.$ (21)

Note that $F^{q}$ is present in interactions with the Higgs and the Z boson
and is the matrix that controls Flavour Changing Neutral Couplings (FCNCs).
For neutrinos,

$\begin{split}\mathcal{L}_{Z}^{\nu}&=-\frac{g}{2\cos{\theta_{W}}}Z_{\mu}[\left(\begin{array}[]{cc}\overline{n_{L}}&\overline{N_{L}}\end{array}\right)~{}F~{}\gamma^{\mu}\left(\begin{array}[]{cc}n_{L}\\\
N_{L}\end{array}\right)]+\,\text{h.c.}~{},\\\ \end{split}$ (22)

$\mathcal{L}_{H}^{\nu}\,=\,-\frac{h}{v}\bigg{[}\begin{pmatrix}\overline{n}_{L}&\overline{N}_{L}\end{pmatrix}\,F~{}\mathcal{D}\,\begin{pmatrix}n_{L}^{c}\\\\[5.69054pt]
N_{L}^{c}\end{pmatrix}\,\bigg{]}\,+\,\text{h.c.}\,,$ (23)

where

$F=A^{\dagger}A~{}.$ (24)

##  4 Parameterisation

In ref. [3], where the neutrino case was studied, it was shown that a unitary
matrix $\mathcal{V}$ can be parameterised as

$\displaystyle\mathcal{V}\,=\,\left(\begin{array}[]{c;{2pt/2pt}c}\,\,\,K&\,\,\,K\,X^\textdagger\\\\[4.2679pt]
\hdashline[2pt/2pt]{-\overline{\\!K}\,\,{X}\,\,}\\!&{\quad\,\,\overline{\\!K}\,\quad\,\,}\\!\\\\[2.84526pt]
\end{array}\right)~{},$ (27)

for a non-singular $3\times 3$ general complex matrix, $K$, a non-singular
$n_{R}\times n_{R}$ general complex matrix $\overline{\\!K}\,$ and a
$n_{R}\times 3$ matrix, $X$, that will be defined next. From eq. (11) one can
identify

$A=(K~{}~{}KX^{\dagger})~{},~{}B=(-\overline{\\!K}\,X~{}~{}\overline{\\!K}\,)~{},$
(28)

where the unitarity of $\mathcal{V}$ yields

$\displaystyle\mathcal{V}\mathcal{V}^{\dagger}=\begin{pmatrix}\,A\,\\\\[5.69054pt]
B\end{pmatrix}\begin{pmatrix}\,{A}^{\dagger}&{B}^{\dagger}\,\end{pmatrix}\,$
$\displaystyle=\,\begin{pmatrix}\,A\,{A}^{\dagger}&A\,{B}^{\dagger}\,\\\\[5.69054pt]
\,B\,{A}^{\dagger}&B\,{B}^{\dagger}\,\end{pmatrix}\,=\,\begin{pmatrix}\,\mathds{1}_{3}\,&\,0\,\\\\[5.69054pt]
\,0\,&\,\mathds{1}_{n_{R}}\,\end{pmatrix}\,,$ (29)
$\displaystyle\mathcal{V}^{\dagger}\mathcal{V}=\begin{pmatrix}\,{A}^{\dagger}&{B}^{\dagger}\,\end{pmatrix}\begin{pmatrix}\,A\,\\\\[5.69054pt]
B\end{pmatrix}\,$
$\displaystyle=\,\,{A}^{\dagger}\,A+{B}^{\dagger}\,B\,\,=\,\,\mathds{1}_{3+n_{R}}\,.$
(30)

Using eq. (28) and the unitarity relations in eqs. (29) and (30), eq. (13)
becomes

$\displaystyle\textbf{0}=d+X^{\dagger}DX^{*}~{},$ (31) $\displaystyle
m\,=\,K\,X^{\dagger}D\left(\overline{\\!K}\,^{-1}\right)^{*}~{},$
$\displaystyle
M\,=\,\overline{\\!K}\,\,(D+X\,d\,X^{T})\,\overline{\\!K}\,^{T}\,.$

The last two equations are exact equations for the $3\times n_{R}$ Dirac mass
matrix, $m$, and the $n_{R}\times n_{R}$ Majorana mass matrix, $M$. The
equation involving the null matrix has the solution

$\displaystyle X\,=\pm\,i\,\sqrt{D^{-1}}\,O_{c}\,\sqrt{d}\,,$ (32)

where $O_{c}$ is an orthogonal complex matrix and $d$ ($D$) is a diagonal
matrix with the masses of the light (heavy) neutrinos in the diagonal. This
unique solution defines the matrix $X$. From the unitarity relations of
$\mathcal{V}$ one can also obtain the following definitions

$\begin{split}&K=U_{K}(\mathds{1}_{3}+X^{\dagger}X)^{-1/2}~{},\\\
&\overline{\\!K}\,=U_{\overline{\\!K}\,}(\mathds{1}_{n_{R}}+XX^{\dagger})^{-1/2}~{},\\\
&K_{PMNS}=K~{},\end{split}$ (33)

where $U_{K}$ and $U_{\overline{\\!K}\,}$ are unitary matrices, since the
unitarity relations only define $K$ and $\overline{\\!K}\,$ up to a unitary
matrix on the left. $K$ will play the role of the PMNS mixing matrix and is
only unitary when $X\xrightarrow[]{}\textbf{0}$.

For the quarks, one can perform the same steps and write

$\displaystyle\mathcal{V}_{\chi}^{q}\,=\,\left(\begin{array}[]{c;{2pt/2pt}c}\,\,\,K_{\chi}^{q}&\,\,\,K_\chi^q\,{X_\chi^q}^\textdagger\\\\[4.2679pt]
\hdashline[2pt/2pt]{-\overline{\\!K}\,_{\\!\chi}^{\\!q}\,X_{\chi}^{q}\,\,}\\!&{\quad\,\,\overline{\\!K}\,_{\\!\chi}^{\\!q}\quad\,\,}\\!\\\\[2.84526pt]
\end{array}\right)~{}.$ (36)

for a non-singular $3\times 3$ general complex matrix, $K^{q}_{\chi}$, a non-
singular $n_{q}\times n_{q}$ general complex matrix
$\overline{\\!K}\,^{q}_{\chi}$ and a $n_{q}\times 3$ matrix, $X^{q}_{\chi}$.
Again, from eq. (10) one obtains

$A^{q}_{\chi}=(K^{q}_{\chi}~{}~{}K^{q}_{\chi}{X^{q}_{\chi}}^{\dagger})~{},~{}B^{q}_{\chi}=(-\overline{\\!K}\,^{q}_{\chi}X^{q}_{\chi}~{}~{}~{}\overline{\\!K}\,^{q}_{\chi})~{},$
(37)

where the unitarity of $\mathcal{V}^{q}_{\chi}$ yields

$\displaystyle\mathcal{V}^{q}_{\chi}{\mathcal{V}^{q}_{\chi}}^{\dagger}=\begin{pmatrix}\,A^{q}_{\chi}\,\\\\[5.69054pt]
B^{q}_{\chi}\end{pmatrix}\begin{pmatrix}\,{A^{q}_{\chi}}^{\dagger}&{B^{q}_{\chi}}^{\dagger}\,\end{pmatrix}\,$
$\displaystyle=\,\begin{pmatrix}\,A^{q}_{\chi}\,{A^{q}_{\chi}}^{\dagger}&A^{q}_{\chi}\,{B^{q}_{\chi}}^{\dagger}\,\\\\[5.69054pt]
\,B^{q}_{\chi}\,{A^{q}_{\chi}}^{\dagger}&B^{q}_{\chi}\,{B^{q}_{\chi}}^{\dagger}\,\end{pmatrix}\,=\,\begin{pmatrix}\,\mathds{1}_{3}\,&\,0\,\\\\[5.69054pt]
\,0\,&\,\mathds{1}_{n_{q}}\,\end{pmatrix}\,,$ (38)
$\displaystyle{\mathcal{V}^{q}_{\chi}}^{\dagger}\mathcal{V}^{q}_{\chi}=\begin{pmatrix}\,{A^{q}_{\chi}}^{\dagger}&{B^{q}_{\chi}}^{\dagger}\,\end{pmatrix}\begin{pmatrix}\,A^{q}_{\chi}\,\\\\[5.69054pt]
B^{q}_{\chi}\end{pmatrix}\,$
$\displaystyle=\,\,{A^{q}_{\chi}}^{\dagger}\,A^{q}_{\chi}+{B^{q}_{\chi}}^{\dagger}\,B^{q}_{\chi}\,\,=\,\,\mathds{1}_{3+n_{q}}\,,$
(39)

for each $q=u,d$ and $\chi=L,R$.

Using eq. (37) and the unitarity relations in eqs. (38) and (39), eq. (12)
becomes

$\begin{split}&m_{q}\,=\,K_{L}^{q}\left(d_{q}+{X_{L}^{q}}^{\dagger}\,D_{q}\,X_{R}^{q}\right){K_{R}^{q}}^{\dagger},\\\
&\overline{m}\hskip
0.86108pt_{q}\,=\,K_{L}^{q}\left({X_{L}^{q}}^{\dagger}\,D_{q}-d_{q}\,{X_{R}^{q}}^{\dagger}\right){\overline{\\!K}\,_{R}^{q}}^{\dagger},\\\
&\overline{\\!M}\mkern-2.0mu_{q}\,=\,\overline{\\!K}\,_{L}^{q}\left(D_{q}\,X_{R}^{q}-X_{L}^{q}\,d_{q}\right){K_{R}^{q}}^{\dagger},\\\
&M_{q}\,=\,\overline{\\!K}\,_{L}^{q}\left(D_{q}+X_{L}^{q}\,d_{q}\,{X_{R}^{q}}^{\dagger}\right){\overline{\\!K}\,_{R}^{q}}^{\dagger},\end{split}$
(40)

where $m_{q}$ is the $3\times 3$ Dirac mass matrix for the quarks and the
entries of the $3\times n_{q}$ matrix, $\overline{m}\hskip 0.86108pt_{q}$, are
proportional to the Higgs vacuum expectation value, as well. The $n_{q}\times
3$ mass matrix, $\overline{\\!M}\mkern-2.0mu_{q}$, and the $n_{q}\times n_{q}$
mass matrix, $M_{q}$, are bare mass terms involving only singlet quark fields.
It is always possible to go to a WB where $\overline{m}\hskip 0.86108pt_{q}$
is 0 222This is always possible for $\overline{\\!M}\mkern-2.0mu_{q}$, as
well.. In that WB, one can proceed like in the neutrino case to obtain a
formula for $X_{\chi}^{q}$,

$\displaystyle\begin{aligned}
&X_{L}^{q}=\sqrt{D_{q}^{-1}}P^{q}\sqrt{d_{q}}~{},\\\
&X_{R}^{q}=\sqrt{D_{q}}P^{q}\sqrt{d_{q}^{-1}}~{},\end{aligned}$ (41)

where $P^{q}$ is a general complex matrix and $d_{q}$ ($D_{q}$) is a diagonal
matrix with the masses of the light (heavy) quarks in the diagonal. The
following definitions can also be obtained from the unitarity relations in
eqs. (38) and (39),

$\begin{split}&K^{q}_{\chi}=U^{q}_{K\chi}(\mathds{1}_{3}+{X_{\chi}^{q}}^{\dagger}X_{\chi}^{q})^{-1/2}\,,\\\
&\overline{\\!K}\,^{q}_{\chi}=U^{q}_{\overline{\\!K}\,\chi}(\mathds{1}_{n_{q}}+X_{\chi}^{q}{X_{\chi}^{q}}^{\dagger})^{-1/2}~{},\\\
&K_{CKM}={K_{L}^{u}}^{\dagger}K_{L}^{d}~{},\end{split}$ (42)

where $U^{q}_{K\chi}$ and $U^{q}_{\overline{\\!K}\,\chi}$ are unitary
matrices, since the unitarity relations only define $K$ and
$\overline{\\!K}\,$ up to a unitary matrix on the left. The combination
${K_{L}^{u}}^{\dagger}K_{L}^{d}$ will play the role of the CKM mixing matrix
and is only unitary when $X_{L}^{u},X_{L}^{d}\xrightarrow[]{}\textbf{0}$.

The matrices relevant for FCNCs, involving the $Z$ and the Higgs, as defined
in the section Interactions, are given by

$F^{q}=\begin{pmatrix}(\mathds{1}_{3}+{X_{L}^{q}}^{\dagger}X_{L}^{q})^{-1}&(\mathds{1}_{3}+{X_{L}^{q}}^{\dagger}{X_{L}^{q}})^{-1}{X_{L}^{q}}^{\dagger}\\\
{X_{L}^{q}}(\mathds{1}_{3}+{X_{L}^{q}}^{\dagger}{X_{L}^{q}})^{-1}&{X_{L}^{q}}(\mathds{1}_{3}+{X_{L}^{q}}^{\dagger}{X_{L}^{q}})^{-1}{X_{L}^{q}}^{\dagger}\\\
\end{pmatrix}~{},$ (43)

for quarks, and

$F=\begin{pmatrix}(\mathds{1}_{3}+X^{\dagger}X)^{-1}&(\mathds{1}_{3}+X^{\dagger}X)^{-1}X^{\dagger}\\\
X(\mathds{1}_{3}+X^{\dagger}X)^{-1}&X(\mathds{1}_{3}+X^{\dagger}X)^{-1}X^{\dagger}\\\
\end{pmatrix}~{},$ (44)

for neutrinos.
Finally, the exact formula for the $(3+n)\times(3+n)$ unitary matrix assumes
the form

$\mathcal{V}^{q}_{\chi}=\begin{pmatrix}U^{q}_{K\chi}(\mathds{1}_{3}+{X_{\chi}^{q}}^{\dagger}X_{\chi}^{q})^{-1/2}&U^{q}_{K\chi}(\mathds{1}_{3}+{X_{\chi}^{q}}^{\dagger}X_{\chi}^{q})^{-1/2}X_{\chi}^{q{\dagger}}\\\
-U^{q}_{\overline{\\!K}\,\chi}(\mathds{1}_{n_{q}}+X_{\chi}^{q}{X_{\chi}^{q}}^{\dagger})^{-1/2}X^{q}_{\chi}&U^{q}_{\overline{\\!K}\,\chi}(\mathds{1}_{n_{q}}+X_{\chi}^{q}{X_{\chi}^{q}}^{\dagger})^{-1/2}\\\
\end{pmatrix}~{},$ (45)

for quarks, and

$\mathcal{V}=\begin{pmatrix}U_{K}(\mathds{1}_{3}+X^{\dagger}X)^{-1/2}&U_{K}(\mathds{1}_{3}+X^{\dagger}X)^{-1/2}X^{\dagger}\\\
-U_{\overline{\\!K}\,}(\mathds{1}_{n_{R}}+XX^{\dagger})^{-1/2}X&U_{\overline{\\!K}\,}(\mathds{1}_{n_{R}}+XX^{\dagger})^{-1/2}\\\
\end{pmatrix}~{},$ (46)

for neutrinos.
Note that parameterisations with a similar structure, used in the leptonic
sector, existed in the literature prior to this work [9][10], but are either
approximations or a special case of this one.

##  5 Procedure and Usefulness

Now that all the necessary equations were derived, one needs to obtain
observables out of them. One procedure is the following

  * •

Start with $d$, $D$, $U_{K}$ and $O_{c}$/$P^{q}$,

  * •

Calculate $X$,

  * •

Calculate mass matrices, charged currents matrix $V$ and the FCNCs matrix $F$.

With these, one should be able to obtain any observable, exact at tree level.
Another valid approach is to proceed in the reverse order, but starting from
the mass matrices would be inherently more difficult due to the redundant
extra parameters.
It is important to emphasize that all physical parameters are contained in
$d$, $D$, $U_{K}$ and $O_{c}$/$P^{q}$. For the neutrino case, a counting for
$n_{R}=3$ was performed in table 1 of ref. [1].
Performing the depicted prescription yields results at tree level that are
exact. However, there are some caveats.
In the case of neutrinos, one should mind radiative corrections on the light
neutrino masses, as stated in section 4 of ref. [1]. The upshot from that
discussion is that one needs some kind of lepton-number like softly broken
symmetry to protect these radiative corrections of becoming too large, when
the mass of the heavy neutrinos is below the GeV scale. This is known in the
literature as “symmetry protected seesaw models” [11].
A concern that exists in every model is perturbativity. This is easy to
understand when one analyses the Dirac mass matrix, $m$, equation in (31) or
(40): If the Yukawa couplings333$m/v$, where $v$ is the Higgs vacuum
expectation value. are $O(1)$ then $X^{\dagger}D$ must also be $O(1)$, which
implies that deviations from unitarity need to decrease when the masses of the
heavy neutrinos increase, or the inverse. This fact was used to set an upper
bound on the mass of the hypothetical up vector-like quark singlet introduced
in ref. [2] to solve the CKM unitarity problem.
The perturbativity issue was discussed in some detail, for quarks, on section
2.5 of ref. [2] and, for neutrinos, on section 2.3 of ref. [1].

##  6 Conclusions

In this article, two models with non-unitary mixing matrices were analysed
through the lens of an exact parameterisation that proved to be extremely
useful. The usefulness of this parameterisation resides in the exact formulas
at tree level, which are easy to implement numerically. Moreover, deviations
from unitarity of the mixing matrix are controlled by the matrix $X$ and no
approximations regarding them or the mass scale of heavy fields is needed to
obtain exact results. Furthermore, it can be used in any model with non-
unitary mixing matrices: Inverse Seesaw, Linear Seesaw, type-II and type-III
seesaw or models with vector like fermions and scalars.
More details regarding this analysis can be found in the following works [3,
1, 2].

## Acknowledgements

This work was partially supported by Fundação para a Ciência e a Tecnologia
(FCT, Portugal) through the projects CFTP-FCT Unit 777 (UIDB/00777/2020 and
UIDP/00777/2020), PTDC/FISPAR/29436/2017, CERN/FIS-PAR/0004/2019 and CERN/FIS-
PAR/0008/2019, which are partially funded through POCTI (FEDER), COMPETE, QREN
and EU. P.M.F.P. acknowledges support from FCT through the PhD grant
SFRH/BD/145399/2019

## References

  * [1] G. C. Branco, J. T. Penedo, P. M. F. Pereira, M. N. Rebelo and J. I. Silva-Marcos, “Type-I Seesaw with eV-Scale Neutrinos,” JHEP 07 (2020), 164 doi:10.1007/JHEP07(2020)164 [arXiv:1912.05875 [hep-ph]].
  * [2] G. C. Branco, J. T. Penedo, P. M. F. Pereira, M. N. Rebelo and J. I. Silva-Marcos, “Addressing the CKM unitarity problem with a vector-like up quark,” JHEP 07 (2021), 099 doi:10.1007/JHEP07(2021)099 [arXiv:2103.13409 [hep-ph]].
  * [3] N. R. Agostinho, G. C. Branco, P. M. F. Pereira, M. N. Rebelo and J. I. Silva-Marcos, “Can one have significant deviations from leptonic 3 $\times$ 3 unitarity in the framework of type I seesaw mechanism?,” Eur. Phys. J. C 78 (2018) no.11, 895 doi:10.1140/epjc/s10052-018-6347-2 [arXiv:1711.06229 [hep-ph]].
  * [4] T. Yanagida, “Horizontal gauge symmetry and masses of neutrinos,” Conf. Proc. C 7902131 (1979), 95-99 KEK-79-18-95.
  * [5] P. Minkowski, “$\mu\to e\gamma$ at a Rate of One Out of $10^{9}$ Muon Decays?,” Phys. Lett. B 67 (1977), 421-428 doi:10.1016/0370-2693(77)90435-X
  * [6] J. A. Casas and A. Ibarra, “Oscillating neutrinos and $\mu\to e,\gamma$,” Nucl. Phys. B 618 (2001), 171-204 doi:10.1016/S0550-3213(01)00475-8 [arXiv:hep-ph/0103065 [hep-ph]].
  * [7] E. Fernandez-Martinez, J. Hernandez-Garcia and J. Lopez-Pavon, “Global constraints on heavy neutrino mixing,” JHEP 08 (2016), 033 doi:10.1007/JHEP08(2016)033 [arXiv:1605.08774 [hep-ph]].
  * [8] J. A. Aguilar-Saavedra, R. Benbrik, S. Heinemeyer and M. Pérez-Victoria, “Handbook of vectorlike quarks: Mixing and single production,” Phys. Rev. D 88 (2013) no.9, 094010 doi:10.1103/PhysRevD.88.094010 [arXiv:1306.0572 [hep-ph]].
  * [9] J. G. Korner, A. Pilaftsis and K. Schilcher, “Leptonic CP asymmetries in flavor changing H0 decays,” Phys. Rev. D 47 (1993), 1080-1086 doi:10.1103/PhysRevD.47.1080 [arXiv:hep-ph/9301289 [hep-ph]].
  * [10] W. Grimus and L. Lavoura, “The Seesaw mechanism at arbitrary order: Disentangling the small scale from the large scale,” JHEP 11 (2000), 042 doi:10.1088/1126-6708/2000/11/042 [arXiv:hep-ph/0008179 [hep-ph]].
  * [11] S. Antusch and O. Fischer, “Testing sterile neutrino extensions of the Standard Model at future lepton colliders,” JHEP 05 (2015), 053 doi:10.1007/JHEP05(2015)053 [arXiv:1502.05915 [hep-ph]].