Electromagnetically Consistent Optimization Algorithms for the Global Design of RIS

Abstract

The reconfigurable intelligent surface is an emerging technology for wireless communications. We model it as an inhomogeneous boundary of surface impedance, and consider various optimization problems that offer different tradeoffs in terms of performance and implementation complexity. The considered non-convex optimization problems are reformulated as a sequence of approximating linear quadratically constrained or semidefinite programs, which are proved to have a polynomial complexity and to converge monotonically in the objective value.

I Introduction

The reconfigurable intelligent surface (RIS) is a physical layer technology that allows information and communication providers to optimize the propagation of electromagnetic waves, hence sculpting favorable communication channels for the efficient transmission and processing of information. Several optimization algorithms for RIS-aided systems are available, but most of them rely on simple (often simplistic) communication models [1]. The development of electromagnetically consistent and tractable communication models for evaluating the performance and optimizing RIS-aided wireless networks, from a signal-level and system-level perspective, is, on the other hand, an important subject open to research [2].

In communication engineering, in addition, current optimization criteria for RIS-aided channels are usually based on the so-called local design, i.e., a specified constraint is imposed to each reconfigurable element of the RIS [3]. In the electromagnetic community, however, the local design criterion is known to be sub-optimal in terms of scattered power [4], and to possibly result in scattered electromagnetic waves towards directions different from the intended one [1], [5], [6]. To overcome these limitations, an RIS needs to be optimized based on the so-called global design [1], [4], which ensures that the total power reflected towards the direction of interest is as close as possible (ideally equal) to the total incident power. This is ensured by imposing a single reflection constraint that encompasses all the reconfigurable elements of the RIS simultaneously. Existing optimization algorithms for the global design of an RIS are, however, based on communication models that are either not electromagnetically consistent [3] or rely on general-purpose optimization functions with no performance guarantee in terms of computational complexity, convergence properties, and optimality of the solution [1], [7].

In this context, we embrace the electromagnetically consistent communication model introduced in [1], which models an RIS as an inhomogeneous boundary of surface impedance. Accordingly, we formulate several optimization problems with the following distinguishing features: (i) electromagnetic consistency of the solution; (ii) specified power efficiency towards the intended direction of reflection; (iii) specified maximum power towards unwanted directions of reflection; and (iv) specified physical implementation constraints, which are all imposed by design in the problem formulation. The considered optimization problems are shown not to be convex. To tackle them efficiently, we reformulate them as a sequence of approximating linear quadratically constrained or semidefinite programs, which are proved to have a polynomial complexity and to converge monotonically in the objective value [8].

The significance of the proposed optimization algorithms is that the global design solution for anomalous reflectors is known only under ideal assumptions, i.e., unitary power efficiency, no parasitic scattering, and no physical constraints on the surface impedance, whose real part needs to be positive and negative [4]. To the best of our knowledge, there exist no optimization algorithms to the design of RISs for which the power efficiency, undesired reradiations, and implementation constrains (e.g., the real part of the surface impedance shall not be negative) are specified as optimization constraints.

Notation: Matrices and column vectors are denoted by bold uppercase and lowercase fonts. $|\cdot|$, ${(\cdot)}^{*}$, $\operatorname{Re}(\cdot)$ denote the absolute value, conjugate, real part. ${(\cdot)}^{H}$, ${(\cdot)}^{T}$, $\operatorname{tr}(\cdot)$ denote the hermitian, transpose, trace operators. $\lVert\cdot\rVert$, $\lVert\cdot\rVert_{*}$, $\lVert\cdot\rVert_{\rm{F}}$ denote the spectral, nuclear, Frobenius norms. $\succeq$ denotes positive semidefinite. ${O}(\cdot)$ stands for the big-O notation. $j$ is the imaginary unit. ${\bf{1}}$ is the all ones vector. $\nabla\left({f\left({\bf{x}}\right)}\right)$ denotes the gradient of ${f\left({\bf{x}}\right)}$ with respect to ${{{\bf{x}}^{*}}}$. $\nabla_{n}\left({f\left({\bf{x}}\right)}\right)$ is the $n$th entry of $\nabla\left({f\left({\bf{x}}\right)}\right)$. $f\left({\bf{x}}\right)\approx f\left({{\bar{\bf{x}}}}\right)+2\operatorname{Re}\left(\nabla^{T}\left({f\left(\bar{\bf{x}}\right)}\right)\left({{\bf{x}}-{\bar{\bf{x}}}}\right)^{*}\right)$ is the first-order Taylor approximation of $f\left({\bf{x}}\right)$ at the point $\bar{\bf{x}}$.

II Electromagnetic Model

II-A RIS Model

We consider the same system model as in [1, Sec. III-A], which encompasses a single-antenna transmitter, a single-antenna receiver, and an RIS (a flat surface $\mathcal{S}$) that is modeled as an inhomogeneous boundary of surface impedance with negligible thickness with respect to the considered wavelength. The RIS is modeled as a rectangle that lies in the $xy$-plane (i.e., $z=0$) with its center located at the origin. Specifically, $\mathcal{S}$ is defined as $\mathcal{S}=\left\{(x,y):|x|\leq L_{x},|y|\leq L_{y}\right\}$, with $2L_{x}$ and $2L_{y}$ being the lengths of $\mathcal{S}$ along the $x$-axis and $y$-axis, respectively. We consider a reflecting RIS, i.e., the transmitter and receiver are located on the same side of $\mathcal{S}$.

The transmitter and receiver are located in the Fraunhofer far-field region of each other and of the RIS. Thus, the incident and reflected signals are modeled as plane waves, whose angles of incidence and reflection, with respect to the normal (i.e., the $z$-axis) to $\mathcal{S}$, are denoted by $\theta_{i}$ and $\theta_{r}$, respectively. As in [1], the incident and reflected signals propagate in the $yz$-plane, so that the dependence on the azimuth angle is ignored. We consider only the signal reflected by the RIS and ignore the transmitter-receiver direct link due to the presence of blocking objects. A free space propagation environment is considered.

Since the RIS is modeled as an inhomogeneous boundary of surface impedance, it is characterized by a surface impedance $Z(x,y)$ or, equivalently, by a surface reflection coefficient $\Gamma(x,y)$ for $(x,y)\in\mathcal{S}$, according to the definitions in [1, Eq. (37)]. Because of the independence from the azimuth angle, $Z(x,y)$ and $\Gamma(x,y)$ are constant functions along the $x$-axis. For system optimization, we can hence consider $Z(x,y)=Z(y)$ and $\Gamma(x,y)=\Gamma(y)$. For ease of writing and numerical implementation when solving the optimization problems, the surface impedance and reflection coefficient are discretized with spatial sampling $\Delta_{x}$ and $\Delta_{y}$ along the $x$-axis and $y$-axis, respectively, as detailed in [1, Sec. III-C]. Accordingly, $Z(y)$ and $\Gamma(y)$ are represented by two column vectors $\mathbf{z}=[z_{1},z_{2},\dots,z_{N}]^{T}$ and $\boldsymbol{\gamma}=[\gamma_{1},\gamma_{2},\dots,\gamma_{N}]^{T}$, respectively. Specifically, $z_{n}=Z(y_{n})$ with $y_{n}=-L_{y}-\Delta_{y}/2+n\Delta_{y}$, $n=1,2,\ldots,N$, and $N=2L_{y}/\Delta_{y}$. The same applies to $\boldsymbol{\gamma}$.

The entries of $\mathbf{z}$ and $\boldsymbol{\gamma}$ are related to one another as follows:

$$z_{n}=\eta_{0}\frac{1+{\gamma_{n}}}{\cos\theta_{i}-{\gamma_{n}}\cos\theta_{r}},\quad\gamma_{n}=\frac{{z_{n}}\cos\theta_{i}-\eta_{0}}{{z_{n}}\cos\theta_{r}+\eta_{0}}\vspace{-0.1cm}$$ (1)

where $\eta_{0}$ is the free space impedance and $n=1,2,\ldots,N$.

II-B Electromagnetic Consistency

The entries of $\mathbf{z}$ and $\boldsymbol{\gamma}$ are arbitrary complex values, with the only requirement that they need to produce reflected electric and magnetic fields, given the incident electric and magnetic fields, that are electromagnetically consistent, i.e., that fulfill Maxwell’s equations [1, p. 1183]. This implies that the entries of $\boldsymbol{\gamma}$ need to satisfy the constraint [1, Eq. (44), Table 5]

$$H_{n}\left(\boldsymbol{\gamma}\right)=\frac{\left|f_{n}^{\prime\prime}-2j\kappa f_{n}^{\prime}\sin\theta_{r}\right|}{\kappa^{2}\left|\gamma_{n}\right|}=0,\quad n=1,\ldots,N-2\vspace{-0.1cm}$$ (2)

where $f_{n}=\gamma_{n}e^{j\kappa\left(\sin\theta_{r}-\sin\theta_{i}\right)y_{n}}$, $f_{n}^{\prime}=({f_{n+1}-f_{n}})/{\Delta_{y}}$, $f_{n}^{\prime\prime}=({f_{n+1}^{\prime}-f_{n}^{\prime}})/{\Delta_{y}}$, $\kappa=2\pi/\lambda$, and $\lambda$ is the considered wavelength. The corresponding constraint that the entries of $\mathbf{z}$ need to fulfill can be found by inserting $\gamma_{n}$ in (1) into (2).

As detailed in [1], the condition $H_{n}=0$ results in the optimal solution, at the highest implementation complexity, as the surface impedance is characterized by large variations, and by positive and negative values of its real part, which make it difficult to implement the resulting RIS in practice [4], [7]. For this reason, the constraint $H_{n}=0$ that ensures the electromagnetic consistency of the solution is usually relaxed, by replacing it with the constraint $\varepsilon_{\rm L}\leq H_{n}\leq\varepsilon_{\rm U}$, where $0\leq\varepsilon_{\rm L}\leq\varepsilon_{\rm U}$ are small positive constants that control the tradeoff between the optimality and electromagnetically consistency (i.e., $\varepsilon_{\rm U}\to 0$) of the obtained design against the implementation complexity (i.e., $\varepsilon_{\rm L}>0$) of the RIS.

II-C Performance Metrics

As detailed in [1, Sec. III-B], the performance of an RIS is completely characterized by two performance metrics.

Surface Net Power Flow – The surface net power flow is defined as the difference between the power reradiated by the whole RIS towards the intended direction of reflection and the total incident power. If the power reradiated towards the intended direction of reflection is hence equal to the total incident power, the surface net power flow is zero. An RIS for which the surface net power flow is equal to zero is defined as globally optimum. To optimize an RIS according to the global design criterion, the surface net power flow needs then to be zero. Considering typical performance versus implementation tradeoffs, the global design criterion is hence tantamount to either minimizing the surface net power flow or ensuring that it is as close as possible to zero within a specified tolerance.

Based on the considered electromagnetic model, the surface net power flow can be formulated, in terms of $\boldsymbol{\gamma}$, as follows:

$$P_{\mathcal{S}}(\boldsymbol{\gamma})={a_{X}}{\Delta_{y}}\left({{c_{i}}+{\alpha_{r}}{{\boldsymbol{\gamma}}^{H}}{\boldsymbol{\gamma}}+0.5{\alpha_{ir}}\left({{{\bf{1}}^{T}}{\boldsymbol{\gamma}}+{{\boldsymbol{\gamma}}^{H}}{\bf{1}}}\right)}\right)\vspace{-0.1cm}$$ (3)

where $a_{X}={\left|E_{0}\right|^{2}L_{x}}/{\eta_{0}}\geq 0$, $\alpha_{i}=\cos\theta_{i}\geq 0$, $\alpha_{r}=\cos\theta_{r}\geq 0$, $c_{i}=-2L_{y}\alpha_{i}/{\Delta_{y}}\leq 0$, $\alpha_{ir}=\alpha_{r}-\alpha_{i}\in[-1,+1]$, $E_{0}$ is the amplitude of the incident electric field (a plane wave).

Power Flux – The power flux characterizes the amount of power reradiated by an RIS at a specified point of observation. We focus on observation points located in the far-field region of the RIS. The power flux provides information on the angular response of the RIS, i.e., how the total incident power is reradiated towards different directions. In the far-field, the power flux is proportional to the radiation pattern of the RIS. The power flux is an essential performance indicator in wireless communications since it determines the amount of received power and hence the signal-to-interference ratio.

Based on the considered electromagnetic model, the power flux evaluated towards a generic direction of reradiation $\theta_{k}$ can be formulated, as a function of $\boldsymbol{\gamma}$, as follows:

$${P_{{\theta_{k}}}}\left({\boldsymbol{\gamma}}\right)={a_{k}}\Delta_{y}^{2}{\chi_{ik}}{\left|{{{\boldsymbol{\gamma}}^{T}}{\bf{u}}_{ik}}\right|^{2}}\vspace{-0.1cm}$$ (4)

where the $n$th entry of vector ${\bf{u}}_{ik}$ is $u_{ik,n}=e^{j\kappa\left(\sin\theta_{k}-\sin\theta_{i}\right)y_{n}}$ for $n=1,2,\ldots,N$, ${a}_{k}=\frac{\kappa^{2}}{\eta_{0}}\frac{\left|E_{0}\right|^{2}L_{x}^{2}}{8\pi^{2}R_{k}^{2}}$, $R_{k}$ is the distance from the RIS towards the direction $\theta_{k}$, and ${\chi}_{ik}=\cos^{2}\theta_{r}+\cos^{2}{\theta}_{k}+2\cos\theta_{r}\cos{\theta}_{k}\geq 0$. If the direction of observation coincides with the intended direction of reflection, i.e., ${\theta}_{k}={\theta}_{r}$, then $R_{k}=R_{r}$, ${\chi}_{ik}={\chi}_{ir}=4\cos^{2}\theta_{r}$, and $P_{{\theta}_{k}}\left(\boldsymbol{\gamma}\right)=P_{{\theta}_{r}}\left(\boldsymbol{\gamma}\right)$.

Next, we formulate and solve optimization problems that aim to minimize $P_{\mathcal{S}}(\boldsymbol{\gamma})$ in (3) (Sec. III) and to maximize $P_{{\theta}_{r}}\left(\boldsymbol{\gamma}\right)$ in (4) (Sec. IV) subject to specified design constraints.

III Optimization: Surface Net Power Flow

In [4], it is shown that minimizing $P_{\mathcal{S}}(\boldsymbol{\gamma})$ in (3) results in an optimal anomalous reflector that steers the total incident power towards the intended direction of reflection with no undesired scattering towards other directions. The optimal surface impedance $\bf{z}$ is known in a closed-form expression only if the constraint in (2) is strictly fulfilled, leading to a high implementation complexity [4], [7]. To the best of our knowledge, there exists no general and efficient optimization framework that aims to minimize $P_{\mathcal{S}}(\boldsymbol{\gamma})$ in (3) by imposing specified design constraints. This is tackled in this section.

III-A Global Design – Helmholtz Constraint

We commence generalizing [4], by relaxing the constraint in (2) within the specified tolerances $\varepsilon_{\rm{HC_{L}}}\geq 0$ and $\varepsilon_{\rm{HC_{U}}}\geq 0$. This keeps under control the variations of the surface impedance $\bf{z}$ along $\mathcal{S}$, facilitating the implementation of the resulting RIS, while ensuring a quasi-Maxwellian solution [1].

The considered problem can be stated as follows:

$$(\textbf{S-HC})\quad \underset{\boldsymbol{\gamma}}{\min}\quad\left|P_{\mathcal{S}}(\boldsymbol{\gamma})\right|$$ (5)

s.t. $$\mathcal{H}_{n}(\boldsymbol{\gamma})\leq\varepsilon_{\rm{HC_{U}}},\quad n=1,2,\ldots,N-2\hskip 36.98866pt\text{(a)}$$

$$\mathcal{H}_{n}(\boldsymbol{\gamma})\geq\varepsilon_{\rm{HC_{L}}},\quad n=1,2,\ldots,N-2\hskip 36.98866pt\text{(b)}\vspace{-0.15cm}$$

The problem S-HC is not convex. To tackle it efficiently, the objective function in (5) is rewritten in epigraph form [9, Eq. (4.11)]. Also, the function $\mathcal{H}_{n}(\boldsymbol{\gamma})$ in (5) is formulated explicitly in terms of the optimization variable $\boldsymbol{\gamma}$, as follows:

$$\mathcal{H}_{n}(\boldsymbol{\gamma})=\frac{{\left|{g_{n}}\left({{\gamma_{n}},{\gamma_{n+1}},{\gamma_{n+2}}}\right)\right|}}{{{\kappa^{2}}\Delta_{y}^{2}\left|{{\gamma_{n}}}\right|}}=\frac{{\left|{{g_{n}}\left({\boldsymbol{\gamma}}_{n}\right)}\right|}}{{{\kappa^{2}}\Delta_{y}^{2}\left|{{\gamma_{n}}}\right|}}\vspace{-0.15cm}$$ (6)

where $\beta_{1}=1+2j\kappa\sin\theta_{r}\Delta_{y}$, $\beta_{2}=-\left(2+2j\kappa\sin\theta_{r}\Delta_{y}\right)$, and ${{g_{n}}\left({\boldsymbol{\gamma}}_{n}\right)}={g_{n}}\left({{\gamma_{n}},{\gamma_{n+1}},{\gamma_{n+2}}}\right)={u_{n+2}}{\gamma_{n+2}}+{\beta_{2}}{u_{n+1}}{\gamma_{n+1}}+{\beta_{1}}{u_{n}}{\gamma_{n}}$ with ${{\boldsymbol{\gamma}}_{n}}={\left({{\gamma_{n}},{\gamma_{n+1}},{\gamma_{n+2}}}\right)^{T}}$.

The problem S-HC can then be rewritten as follows:

$$(\textbf{S-HC-a})\quad\underset{\mathbf{\boldsymbol{\gamma}},t}{\min}\quad t\quad\text{s.t.}\quad-t\leq 0\;\text{(a)}$$ (7)

$$P_{\mathcal{S}}(\boldsymbol{\gamma})-t\leq 0\;\text{(b)},\quad\left|{{g_{n}}\left({\boldsymbol{\gamma}}_{n}\right)}\right|-{{\tilde{\varepsilon}}_{{\mathop{\rm{HC_{U}}}\nolimits}}}\left|{{\gamma_{n}}}\right|\leq 0\;\text{(c)}$$

$$-P_{\mathcal{S}}(\boldsymbol{\gamma})-t\leq 0\;\text{(d)},\quad-\left|{{g_{n}}\left({\boldsymbol{\gamma}}_{n}\right)}\right|+{{\tilde{\varepsilon}}_{{\mathop{\rm{HC_{L}}}\nolimits}}}\left|{{\gamma_{n}}}\right|\leq 0\;\text{(e)}\vspace{-0.15cm}$$

where ${{\tilde{\varepsilon}}_{{\mathop{\rm{HC_{L,U}}}\nolimits}}}={\varepsilon_{{\mathop{\rm{HC_{L,U}}}\nolimits}}}{\kappa^{2}}\Delta_{y}^{2}$ and $n=1,2,\ldots,N-2$.

Problem S-HC-a is still not convex due to the constraints (7c)-(7e). To tackle it, we embrace the iterative inner approximation framework [8]. Specifically, we consider convex upper bounds for the concave functions in (7c)-(7e), replacing them with their first-order Taylor approximation. As for the function ${{{\boldsymbol{\gamma}}^{H}}{\boldsymbol{\gamma}}}$, we note that $\nabla\left({{{\boldsymbol{\gamma}}^{H}}{\boldsymbol{\gamma}}}\right)={\boldsymbol{\gamma}}$ and introduce the following lower bound ${a_{X}}{\Delta_{y}}{\alpha_{r}}{{{{\boldsymbol{\gamma}}}}^{H}}{{\boldsymbol{\gamma}}}\geq{p_{\mathcal{S}}}\left({\boldsymbol{\gamma}},\bar{\boldsymbol{\gamma}}\right)$ at the point ${\bar{\boldsymbol{\gamma}}}$:

$${p_{\mathcal{S}}}\left({\boldsymbol{\gamma}},\bar{\boldsymbol{\gamma}}\right)={a_{X}}{\Delta_{y}}{\alpha_{r}}\left({{{\bar{\boldsymbol{\gamma}}}}^{H}}{\bar{\boldsymbol{\gamma}}}+2\operatorname{Re}\left({\nabla^{T}}\left({{\bar{\boldsymbol{\gamma}}^{H}}{\bar{\boldsymbol{\gamma}}}}\right)\left({{\boldsymbol{\gamma}}-{\bar{\boldsymbol{\gamma}}}}\right)^{*}\right)\right)\vspace{-0.15cm}$$

As for the function $\left|{{g_{n}}\left({\boldsymbol{\gamma}}_{n}\right)}\right|$, we note that ${\nabla_{n}}\left({\left|{{g_{n}}\left({{{\boldsymbol{\gamma}}_{n}}}\right)}\right|}\right)={{0.5\beta_{1}^{*}u_{n}^{*}{g_{n}}\left({{{\boldsymbol{\gamma}}_{n}}}\right)}\mathord{\left/{\vphantom{{0.5\beta_{1}^{*}u_{n}^{*}{g_{n}}\left({{{\boldsymbol{\gamma}}_{n}}}\right)}{\left|{{g_{n}}\left({{{\boldsymbol{\gamma}}_{n}}}\right)}\right|}}}\right.\kern-1.2pt}{\left|{{g_{n}}\left({{{\boldsymbol{\gamma}}_{n}}}\right)}\right|}}$ and introduce the following lower bound $\left|{{g_{n}}\left({{{\boldsymbol{\gamma}}_{n}}}\right)}\right|\geq\left|{g_{n}^{\rm{L}}\left({{{\boldsymbol{\gamma}}_{n}}}\right)}\right|$ evaluated at the point ${{{{\boldsymbol{\bar{\gamma}}}}_{n}}}$:

$$\left|{g_{n}^{\rm{L}}\left({{{\boldsymbol{\gamma}}_{n}}}\right)}\right|=\left|{{g_{n}}\left({{{{\bar{\boldsymbol{\gamma}}}}_{n}}}\right)}\right|+2{\mathop{\rm Re}\nolimits}\left({{\nabla_{n}}\left({\left|{{g_{n}}\left({{{{\bar{\boldsymbol{\gamma}}}}_{n}}}\right)}\right|}\right){{\left({{\gamma_{n}}-{{\bar{\gamma}}_{n}}}\right)}^{*}}}\right)\vspace{-0.15cm}$$

Therefore, the following reformulation for S-HC-a is obtained at the generic iteration of the inner algorithm:

$$(\textbf{S-HC-b})\quad\underset{\mathbf{\boldsymbol{\gamma}},t}{\min}\quad t\quad\text{s.t.}\quad(\ref{Eq:S-HC-a}\text{a}),\;(\ref{Eq:S-HC-a}\text{b})$$ (8)

$$-{a_{X}}{\Delta_{y}}\left({{c_{i}}+0.5{\alpha_{ir}}\left({{{\bf{1}}^{T}}{\boldsymbol{\gamma}}+{{\boldsymbol{\gamma}}^{H}}{\bf{1}}}\right)}\right)-{p_{\mathcal{S}}}\left({\boldsymbol{\gamma}},\bar{\boldsymbol{\gamma}}\right)-t\leq 0\hskip 3.69899pt\text{(a)}$$

$$\left|{g_{n}}\left({{{{\boldsymbol{\gamma}}}_{n}}}\right)\right|-{{\tilde{\varepsilon}}_{{\mathop{\rm{HC_{U}}}\nolimits}}}\left(\left|{{{\bar{\gamma}}_{n}}}\right|+2\operatorname{Re}\left(\nabla_{n}\left({\left|{{{\bar{\gamma}}_{n}}}\right|}\right)\left({{\gamma_{n}}-{{\bar{\gamma}}_{n}}}\right)^{*}\right)\right)\leq 0\;\text{(b)}$$

$$-\left|{{g_{n}^{\rm L}}\left({\boldsymbol{\gamma}}_{n}\right)}\right|+{{\tilde{\varepsilon}}_{{\mathop{\rm{HC_{L}}}\nolimits}}}\left|{{\gamma_{n}}}\right|\leq 0\;\text{(c)},\quad\left\|{{\boldsymbol{\gamma}}-{\bar{\boldsymbol{\gamma}}}}\right\|\leq{\varepsilon_{{\rm{TR}}}}\;\text{(d)}\vspace{-0.15cm}$$

where ${\bar{\boldsymbol{\gamma}}}$ is the point at which Taylor’s approximation is made (the solution of the preceding iteration of the inner algorithm), ${\bar{\boldsymbol{\gamma}}_{n}}={\left({{{\bar{\gamma}}_{n}},{{\bar{\gamma}}_{n+1}},{{\bar{\gamma}}_{n+2}}}\right)^{T}}$, and ${\nabla_{n}}\left({\left|{{\gamma_{n}}}\right|}\right)=0.5{{{\gamma_{n}}}\mathord{\left/{\vphantom{{{\gamma_{n}}}{\left|{{\gamma_{n}}}\right|}}}\right.\kern-1.2pt}{\left|{{\gamma_{n}}}\right|}}$. Eq. (8d) is the trust region, which ensures that, at each iteration, the set of feasible solutions is limited to the points for which the approximation is sufficiently accurate [10]. The radius of the trust region is controlled via the small positive constant ${\varepsilon_{{\rm{TR}}}}$, which is updated at each iteration as detailed in Sec. V.

S-HC-b is convex and is solved as detailed in Sec. III-D.

III-B Global Design – Reradiation Mask

Fulfilling the constraint in (2) ensures the electromagnetic consistency of the obtained solution. The optimization problem has, however, a number of constraints that is equal to the number of optimization variables $N$. Motivated by recent results [1], [5], [6], [7], we show that the constraint in (2) can be replaced by keeping the power scattered towards specified (undesired) directions of reradiation below a given maximum level. The advantage of this formulation is that the unwanted directions of reradiation may be known apriori [5], [11], and their number may be much less than the number of optimization variables $N$. The resulting optimization constraint is referred to as reradiation mask constraint [1], [6].

The considered problem can be stated as follows:

$$(\textbf{S-RM})\quad \underset{\boldsymbol{\gamma}}{\min}\quad\left|P_{\mathcal{S}}(\boldsymbol{\gamma})\right|$$ (9)

s.t. $${P_{{\theta_{k}}}}\left({\boldsymbol{\gamma}}\right)\leq{\varepsilon_{{\rm{RM}}}},\quad{\theta_{k}}\in\mathcal{K}\vspace{-0.1cm}$$

where $\mathcal{K}$ is the set of undesired directions of radiation and ${\varepsilon_{{\rm{RM}}}}$ is the maximum amount of radiated power towards them.

The problem S-RM can be tackled efficiently by rewriting the objective function in epigraph form [9, Eq. (4.11)], and by making explicit the constraint in (9) by using (4), as follows:

$$(\textbf{S-RM-a})\quad \underset{\mathbf{\boldsymbol{\gamma}},t}{\min}\quad t\quad\text{s.t.}\quad(\ref{Eq:S-HC-a}\text{a}),\;(\ref{Eq:S-HC-a}\text{b}),\;(\ref{Eq:S-HC-b}\text{a}),\;(\ref{Eq:S-HC-b}\text{d})$$ (10)

$$\left|{{\boldsymbol{\gamma}}^{T}}{{\bf{u}}_{ik}}\right|-\sqrt{{{\tilde{\varepsilon}}_{{\rm{RM}}}}}\leq 0,\quad{\theta_{k}}\in\mathcal{K}\hskip 48.36958pt(\text{a})\vspace{-0.1cm}$$

where ${{\tilde{\varepsilon}}_{{\rm{RM}}}}=\frac{{{\varepsilon_{{\rm{RM}}}}}}{{{a_{k}}\Delta_{y}^{2}{\chi_{ik}}}}$. S-RM-a is convex and can be solved efficiently as detailed in Sec. III-D.

III-C Approximated Global Design – Reactive Impedance

The optimization problems S-HC and S-RM aim to maximize the radiation efficiency of the RIS, by simultaneously maximizing and minimizing the power reradiated towards the intended and undesired directions of reradiation, respectively. S-HC and S-RM impose, however, mild or no implementation constraint to the feasible set of solutions, respectively. As illustrated in [1], [4], the surface impedance $\bf{z}$ obtained by solving P-HC and P-RM has a negative real part, which makes the obtained design difficult to be implemented [4]. Thus, we consider an optimization problem that imposes specified implementation constraints by design. A sought after implementation requirement is that the real part of $\bf{z}$ is not negative and as small as possible, ideally equal to zero, in order to minimize the power losses. Since the real part of $\bf{z}$ cannot take negative values, the obtained design cannot be deemed globally optimum but only approximately globally optimum.

The considered problem can be stated as follows:

$$(\textbf{S-RI})\quad \underset{\boldsymbol{\gamma}}{\min}\quad\left|P_{\mathcal{S}}(\boldsymbol{\gamma})\right|$$ (11)

s.t. $${P_{{\theta_{k}}}}\left({\boldsymbol{\gamma}}\right)\leq{\varepsilon_{{\rm{RM}}}},\quad{\theta_{k}}\in\mathcal{K}\hskip 83.93553pt\text{(a)}$$

$$\operatorname{Re}\left({{z}}_{n}\right)\geq 0,\quad\;\;\;n=1,2,\ldots,N\hskip 51.21504pt\text{(b)}$$

$$\operatorname{Re}\left({{z}}_{n}\right)\leq{\varepsilon_{{\rm{RI}}}},\quad n=1,2,\ldots,N\hskip 50.64586pt\text{(c)}\vspace{-0.1cm}$$

where ${\varepsilon_{{\rm{RI}}}}\geq 0$ is used to adjust the tradeoff between the power losses and implementation complexity. If ${\varepsilon_{{\rm{RI}}}}=0$, the losses are zero but the implementation complexity is the highest.

The objective function in (11) and the constraint in (11a) can be tackled as in S-RM-a. The constraints in (11b) and (11c) can be reformulated in terms of ${\boldsymbol{\gamma}}$ from (1). For each pair $(z_{n},\gamma_{n})$, (11b) and (11c) are equivalent to the following:

$$\text{(\ref{Eq:S-RI}\text{b})}:\quad{\alpha_{r}}{\left|{{\gamma_{n}}}\right|^{2}}-\left({{\alpha_{i}}-{\alpha_{r}}}\right)\operatorname{Re}\left({{\gamma_{n}}}\right)-{\alpha_{i}}\leq 0$$ (12)

$$\text{(\ref{Eq:S-RI}\text{c})}:\quad-\left({1+{{\tilde{\varepsilon}}_{{\rm{RI}}}}{\alpha_{r}}}\right){\alpha_{r}}{\left|{{\gamma_{n}}}\right|^{2}}+{\alpha_{i}}\left({1-{{\tilde{\varepsilon}}_{{\rm{RI}}}}{\alpha_{i}}}\right)\vspace{-0.1cm}$$

$$\hskip 42.67912pt+\left({{\alpha_{i}}-{\alpha_{r}}+2{{\tilde{\varepsilon}}_{{\rm{RI}}}}{\alpha_{i}}{\alpha_{r}}}\right){\mathop{\rm Re}\nolimits}\left({{\gamma_{n}}}\right)\leq 0\vspace{-0.1cm}$$ (13)

where ${{\tilde{\varepsilon}}_{{\rm{RI}}}}={{{\varepsilon_{{\rm{RI}}}}}\mathord{\left/{\vphantom{{{\varepsilon_{{\rm{RI}}}}}{{\eta_{0}}}}}\right.\kern-1.2pt}{{\eta_{0}}}}$. The constraint in (12) is convex. The constraint in (13) is concave, and it is tackled through an upper bound obtained from Taylor’s approximation applied to ${\left|{{\gamma_{n}}}\right|^{2}}$.

By introducing the function $\psi\left({{\gamma_{n}}}\right)={\alpha_{i}}\left({1-{{\tilde{\varepsilon}}_{{\rm{RI}}}}{\alpha_{i}}}\right)+\left({{\alpha_{i}}-{\alpha_{r}}+2{{\tilde{\varepsilon}}_{{\rm{RI}}}}{\alpha_{i}}{\alpha_{r}}}\right){\mathop{\rm Re}\nolimits}\left({{\gamma_{n}}}\right)$ and considering $\nabla_{n}\left({{{\boldsymbol{\gamma}}^{H}}{\boldsymbol{\gamma}}}\right)={{\gamma}_{n}}$, S-RI can be reformulated as follows (for $n=1,2,\ldots,N$):

$$(\textbf{S-RI-a})\quad\underset{\mathbf{\boldsymbol{\gamma}},t}{\min}\quad t\quad\text{s.t.}\;(\ref{Eq:S-HC-a}\text{a}),\;(\ref{Eq:S-HC-a}\text{b}),\;(\ref{Eq:S-HC-b}\text{a}),\;(\ref{Eq:S-HC-b}\text{d}),\;(\ref{Eq:S-RM-a}\text{a}),\;\eqref{Eq:ReZ_1}$$

$$\psi\left({{\gamma_{n}}}\right)-{\tilde{\alpha}_{r}}\left({\left|{{{\bar{\gamma}}_{n}}}\right|^{2}}+2\operatorname{Re}\left(\nabla_{n}\left({{\bar{\boldsymbol{\gamma}}^{H}}{\bar{\boldsymbol{\gamma}}}}\right)\left({{\gamma_{n}}-{{\bar{\gamma}}_{n}}}\right)^{*}\right)\right)\leq 0\vspace{-0.1cm}$$

where ${\tilde{\alpha}_{r}}={\left({1+{{\tilde{\varepsilon}}_{{\rm{RI}}}}{\alpha_{r}}}\right){\alpha_{r}}}$. The problem S-RI-a is convex and can be solved efficiently as detailed in Sec. III-D.

III-D Convergence and Complexity

The complete algorithm to solve S-HC-b, S-RM-a, and S-RI-a is an instance of sequential programming. The trust region is, in fact, applied in the neighborhood of $\bar{\boldsymbol{\gamma}}$, which is, by design, a feasible point at any iteration. Thus, the proposed algorithm generates a monotonically decreasing sequence of objective values that converges in the objective [8]. The convex problem solved at each iteration has a linear objective and quadratic constraints, which can be tackled using the interior-point method. Thus, the total arithmetic cost per iteration is ${\mathcal{O}}\left({c^{1/2}}\left({c{v^{2}}+{v^{3}}}\right)\right)$ [12, p. 232], with $v$ and $c$ being the numbers of optimization variables and constraints, respectively.

IV Optimization: Power Flux

In electromagnetic theory, the surface net power flow is the optimality criterion usually used to design perfect anomalous reflectors [4]. In communication theory, the optimization problem is usually formulated by maximizing the power scattered towards the intended direction of reflection (known as the power flux), while imposing specified constraints to the power efficiency of the RIS [3]. The two problem formulations are naturally related to one another, but the corresponding solutions and algorithms are different. Often, the problem formulation considered in communication theory is more general, as (i) it can be extended to different objective functions, besides the reflected power, and scenarios, and (ii) it can account for the surface net power flow as an optimization constraint, resulting in specified performance versus implementation complexity tradeoffs by design [3]. However, the computational complexity and memory requirements of the power flux optimization problem are usually higher.

In this section, therefore, we analyze the globally optimum design of an RIS under the lenses of communication theory, by considering alternative formulations for S-RM and S-RI. S-HC is not considered for brevity and because it is de facto equivalent to S-RM from a communication perspective.

IV-A Global Design – Reradiation Mask

By considering the power scattered towards the intended direction of reflection as the objective function and the surface net power flow as an optimization constraint, a problem de facto equivalent to S-RM can be formulated as follows:

$$(\textbf{P-RM})\quad\underset{\boldsymbol{\gamma}}{\max}\quad{P_{{\theta_{r}}}}\left({\boldsymbol{\gamma}}\right)$$ (14)

$$\text{s.t.}\quad{P_{{\theta_{k}}}}\left({\boldsymbol{\gamma}}\right)\leq{\varepsilon_{{\rm{RM}}}},\;{\theta_{k}}\in\mathcal{K}\;\text{(a)},\quad\left|P_{\mathcal{S}}(\boldsymbol{\gamma})\right|\leq\varepsilon_{\rm{SP}}\;\text{(b)}\vspace{-0.1cm}$$

where $\varepsilon_{\rm{SP}}\geq 0$ is used to adjust the tradeoff between the amount of power scattered towards the intended direction of reflection and the implementation complexity of the RIS.

P-RM is not convex. To tackle it, we introduce ${\boldsymbol{\Gamma}}={\boldsymbol{\gamma}}{{\boldsymbol{\gamma}}^{H}}$, ${\bf{U}}_{ik}={\bf{u}}_{ik}{{\bf{u}}_{ik}^{H}}$. Since ${{{\boldsymbol{\gamma}}^{H}}{\boldsymbol{\gamma}}}=\operatorname{tr}\left({\boldsymbol{\Gamma}}\right)$ in (3), we define

$${{\hat{P}}_{\mathcal{S}}}\left({{\boldsymbol{\Gamma}},{\boldsymbol{\gamma}}}\right)={a_{X}}{\Delta_{y}}\left({{c_{i}}+{\alpha_{r}}\operatorname{tr}\left({\boldsymbol{\Gamma}}\right)+0.5{\alpha_{ir}}\left({{{\bf{1}}^{T}}{\boldsymbol{\gamma}}+{{\boldsymbol{\gamma}}^{H}}{\bf{1}}}\right)}\right)\vspace{-0.1cm}$$

Then, P-RM can be equivalently formulated as follows:

$$(\textbf{P-RM-a})\quad\underset{\boldsymbol{\Gamma},\boldsymbol{\gamma}}{\max}\quad{a_{r}}\Delta_{y}^{2}{\chi_{ir}}\operatorname{tr}\left({{\boldsymbol{\Gamma}}{\bf{U}}_{ir}^{*}}\right)$$ (15)

$$\text{s.t.}\quad{a_{k}}\Delta_{y}^{2}{\chi_{ik}}\operatorname{tr}\left({{\boldsymbol{\Gamma}}{\bf{U}}_{ik}^{*}}\right)-{\varepsilon_{{\rm{RM}}}}\leq 0,\quad{\theta_{k}}\in\mathcal{K}\hskip 12.80365pt\text{(a)}$$

$$\left|{\hat{P}}_{\mathcal{S}}(\boldsymbol{\Gamma},\boldsymbol{\gamma})\right|-\varepsilon_{\rm{SP}}\leq 0\quad\text{(b)},\quad{\boldsymbol{\Gamma}}={\boldsymbol{\gamma}}{{\boldsymbol{\gamma}}^{H}}\quad\text{(c)}\vspace{-0.15cm}$$

The only non-convex constraint in P-RM-a is (15c). To tackle it, we reformulate it equivalently as follows:

$$\text{(\ref{Eq:P-RM-a}\text{c})}:\;{\left\|{\boldsymbol{\Gamma}}\right\|_{*}}-{\left\|{\boldsymbol{\Gamma}}\right\|_{\rm{F}}}\leq 0\;\text{(a)}$$ (16)

$$\text{(\ref{Eq:P-RM-a}\text{c})}:\;\left[{\begin{array}[]{*{20}{c}}{\boldsymbol{\Gamma}}&{\boldsymbol{\gamma}}\\ {{{\boldsymbol{\gamma}}^{H}}}&1\end{array}}\right]\succeq 0\;\text{(b)},\quad{\left\|{\boldsymbol{\Gamma}}\right\|}-{\left\|{\boldsymbol{\gamma}}\right\|^{2}}\leq 0\;\text{(c)}\vspace{-0.1cm}$$ (19)

where (16a) ensures that ${\boldsymbol{\Gamma}}$ has rank one and (16b) ensures that ${\boldsymbol{\Gamma}}$ is positive semidefinite. Also, (16b) and (16c) ensure that ${\left\|{\boldsymbol{\Gamma}}\right\|}-{\left\|{\boldsymbol{\gamma}}\right\|^{2}}=0$, i.e., the only singular value of ${\left\|{\boldsymbol{\Gamma}}\right\|}$ is equal to ${\left\|{\boldsymbol{\gamma}}\right\|^{2}}$. This is because the spectral norm is monotone, i.e., ${\boldsymbol{\Gamma}}\succeq{\boldsymbol{\gamma}}{{\boldsymbol{\gamma}}^{H}}$ in (16b) implies $\left\|{\boldsymbol{\Gamma}}\right\|\geq\left\|{{\boldsymbol{\gamma}}{{\boldsymbol{\gamma}}^{H}}}\right\|={\left\|{\boldsymbol{\gamma}}\right\|^{2}}$.

The constraints in (16a) and (16d) are not convex, since they are given by the difference of two convex functions. We tackle them by applying the iterative inner approximation framework [8]. Specifically, we consider the following convex lower bounds for the functions ${\left\|{\boldsymbol{\gamma}}\right\|^{2}}$ and ${\left\|{\boldsymbol{\Gamma}}\right\|_{\rm{F}}}$:

$${f_{\rm{V}}}\left({\boldsymbol{\gamma}},\bar{\boldsymbol{\gamma}}\right)={\left\|{{\bar{\boldsymbol{\gamma}}}}\right\|^{2}}+2{\mathop{\rm Re}\nolimits}\left({{\nabla^{T}}\left({{{\left\|\bar{\boldsymbol{\gamma}}\right\|}^{2}}}\right)\left({{\boldsymbol{\gamma}}-{\bar{\boldsymbol{\gamma}}}}\right)^{*}}\right)$$ (20)

$${f_{\rm{F}}}\left({\boldsymbol{\Gamma}},\bar{\boldsymbol{\Gamma}}\right)={\left\|{{\bar{\boldsymbol{\Gamma}}}}\right\|_{\rm{F}}}+2{\mathop{\rm Re}\nolimits}\left({\sum\nolimits_{n,m}{{\delta_{n,m}}\left(\bar{\boldsymbol{\Gamma}}\right){{\left({{\Gamma_{n,m}}-{{\bar{\Gamma}}_{n,m}}}\right)}^{*}}}}\right)\vspace{-0.1cm}$$

with ${\left\|{\boldsymbol{\gamma}}\right\|^{2}}\geq{f_{\rm{V}}}\left({\boldsymbol{\gamma}},\bar{\boldsymbol{\gamma}}\right)$, ${\left\|{\boldsymbol{\Gamma}}\right\|_{\rm{F}}}\geq{f_{\rm{F}}}\left({\boldsymbol{\Gamma}},\bar{\boldsymbol{\Gamma}}\right)$, ${{\Gamma_{n,m}}}$ is the $(n,m)$th entry of ${\boldsymbol{\Gamma}}$, $\nabla\left({{{\left\|{\boldsymbol{\gamma}}\right\|}^{2}}}\right)={\boldsymbol{\gamma}}$, ${\delta_{n,m}}\left({\boldsymbol{\Gamma}}\right)=0.5{{{\Gamma_{n,m}}}\mathord{\left/{\vphantom{{{\Gamma_{n,m}}}{{{\left\|{\boldsymbol{\Gamma}}\right\|}_{\rm{F}}}}}}\right.\kern-1.2pt}{{{\left\|{\boldsymbol{\Gamma}}\right\|}_{\rm{F}}}}}$, and ${\bar{\boldsymbol{\Gamma}}}$ and ${\bar{\boldsymbol{\gamma}}}$ are the points at which Taylor’s approximation is made (the solution of the preceding iteration).

At the generic iteration of the inner algorithm, therefore, P-RM-a is tackled by solving the following problem:

$$(\textbf{P-RM-b})\quad\underset{\boldsymbol{\Gamma},\boldsymbol{\gamma}}{\max}\quad{a_{r}}\Delta_{y}^{2}{\chi_{ir}}\operatorname{tr}\left({{\boldsymbol{\Gamma}}{\bf{U}}_{ir}^{*}}\right)$$ (21)

$$\text{s.t.}\quad(\ref{Eq:P-RM-a}\text{a}),\;(\ref{Eq:P-RM-a}\text{b}),(\ref{Eq:SDP_Rank1Constraint}\text{b})\quad\text{(a)}$$

$${\left\|{\boldsymbol{\Gamma}}\right\|_{*}}-{f_{\rm{F}}}\left({{\boldsymbol{\Gamma}},{\bar{\boldsymbol{\Gamma}}}}\right)\leq 0\quad\text{(b)},\quad{\left\|{\boldsymbol{\Gamma}}\right\|}-{f_{\rm{V}}}\left({{\boldsymbol{\gamma}},{\bar{\boldsymbol{\gamma}}}}\right)\leq 0\quad\text{(c)}$$

$$\left\|{{\boldsymbol{\gamma}}-{\bar{\boldsymbol{\gamma}}}}\right\|\leq{\varepsilon_{{\rm{TR,}}{\boldsymbol{\gamma}}}}\quad\text{(d)},\quad\left\|{{\boldsymbol{\Gamma}}-{\bar{\boldsymbol{\Gamma}}}}\right\|\leq{\varepsilon_{{\rm{TR,}}{\boldsymbol{\Gamma}}}}\quad\text{(e)}\vspace{-0.1cm}$$

where ${\varepsilon_{{\rm{TR,}}{\boldsymbol{\gamma}}}}\geq 0$, ${\varepsilon_{{\rm{TR,}}{\boldsymbol{\Gamma}}}}\geq 0$, and (21e) and (21d) are trust region constraints, similar to the constraint in (8d).

P-RM-b is convex and is solved as detailed in Sec. IV-C.

IV-B Approximated Global Design – Reactive Impedance

By considering the power scattered towards the intended direction of reflection as the objective function and the surface net power flow as an optimization constraint, a problem de facto equivalent to S-RI can be formulated by adding the constraints in (11b) and (11c), or, equivalently, the constraints in (12) and (13), to the problem P-RM-b in (21).

By setting ${\left|{{\gamma_{n}}}\right|^{2}}={\Gamma_{n,n}}$, we then obtain ($n=1,2,\ldots,N$):

$$(\textbf{P-RI})\quad \underset{\boldsymbol{\Gamma},\boldsymbol{\gamma}}{\max}\quad{a_{r}}\Delta_{y}^{2}{\chi_{ir}}\operatorname{tr}\left({{\boldsymbol{\Gamma}}{\bf{U}}_{ir}^{*}}\right)$$ (22)

$$\text{s.t.}\quad(\ref{Eq:P-RM-a}\text{a}),\;(\ref{Eq:P-RM-a}\text{b}),\;(\ref{Eq:SDP_Rank1Constraint}\text{b}),\;(\ref{Eq:P-RM-b}\text{b})\text{-}(\ref{Eq:P-RM-b}\text{e})$$

$${\alpha_{r}}{\Gamma_{n,n}}-\left({{\alpha_{i}}-{\alpha_{r}}}\right){\mathop{\rm Re}\nolimits}\left({{\gamma_{n}}}\right)-{\alpha_{i}}\leq 0$$

$$-\left({1+{{\tilde{\varepsilon}}_{{\rm{RI}}}}}\right){\alpha_{r}}{\Gamma_{n,n}}+{\alpha_{i}}\left({1-{{\tilde{\varepsilon}}_{{\rm{RI}}}}{\alpha_{i}}}\right)$$

$$\hskip 28.45274pt+\left({{\alpha_{i}}-{\alpha_{r}}+2{{\tilde{\varepsilon}}_{{\rm{RI}}}}{\alpha_{i}}{\alpha_{r}}}\right){\mathop{\rm Re}\nolimits}\left({{\gamma_{n}}}\right)\leq 0\vspace{-0.1cm}$$

P-RI is convex and is solved as detailed in Sec. IV-C.

IV-C Convergence and Complexity

Problems P-RM-b and P-RI fulfill the same convergence properties as the problems analyzed in Sec. III-C. The only difference is that the convex problem solved at each iteration is a semidefinite program [13], which can be tackled by using the interior-point method. By utilizing the same notation as that in Sec. III-C, the total arithmetic cost per iteration is ${\mathcal{O}}\left({c^{1/2}}\left({{c^{3}}v+{c^{2}}{v^{2}}+{v^{3}}}\right)\right)$ [12, p. 247].

V Numerical Results

In this section, we focus on the problems S-RI-a and P-RI, as they encompass all the others, are the most challenging to solve, and analytical solutions are not known. The algorithms are implemented in CVX. The simulation parameters are $f=28$ GHz, $\theta_{i}=0^{\circ}$, $\theta_{r}=60^{\circ}$, $R_{r}=R_{k}=100$ m, $\eta_{0}=377$ $\Omega$, $E_{0}=1$ Watt/m², $L_{x}=0.5$ m, $L_{y}=4.9652\lambda$ m, $\Delta_{y}=\lambda/6.0420$, $N=60$. Also, $\varepsilon_{\rm{RM}}=2\cdot 10^{-8}$, $\varepsilon_{\rm{RI}}=10^{-2}$, $\varepsilon_{\rm{SP}}=10^{-9}$, and the reradiation mask comprises the set of angles $[-2,2]$ and $[-58,-62]$ with angular resolution 0.1^∘.

As for the trust region, the radii in (8d), (21d), (21e) are set to large values at the first iteration and are progressively reduced at each iteration by 1.2 for S-RI-a and by 1.1 for P-RI. If the objective does not decrease due to numerical inaccuracies, the following approach is used: the algorithm steps back to the solution and setup attained three iterations earlier, while the radii of the trust region are reduced by 1.2 for S-RI-a and by 1.1 for P-RI. At convergence, we attained ${\varepsilon_{\rm{TR}}}=8.79\cdot 10^{-13}$, ${\varepsilon_{{\rm{TR,}}{\boldsymbol{\gamma}}}}=8.53\cdot 10^{-6}$, ${\varepsilon_{{\rm{TR,}}{\boldsymbol{\Gamma}}}}=1.02\cdot 10^{-4}$.

To ease the convergence of P-RI, (21b) and (21c) are rewritten as ${\left\|{\boldsymbol{\Gamma}}\right\|_{*}}-{f_{\rm{F}}}\left({{\boldsymbol{\Gamma}},{\bar{\boldsymbol{\Gamma}}}}\right)\leq{\varepsilon_{{\rm{rk1-b}}}}$ and $\left\|{\boldsymbol{\Gamma}}\right\|-{f_{\rm{V}}}\left({{\boldsymbol{\gamma}},{\bar{\boldsymbol{\gamma}}}}\right)\leq{\varepsilon_{{\rm{rk1-c}}}}$, with ${\varepsilon_{{\rm{rk1-b}}}}\geq 0$ and ${\varepsilon_{{\rm{rk1-c}}}}\geq 0$ small positive values. The impact is minor, as it only implies that the difference of the single singular value of ${\left\|{\boldsymbol{\Gamma}}\right\|}$ and ${\left\|{\boldsymbol{\gamma}}\right\|^{2}}$ is less than $\varepsilon_{{\rm{rk1-c}}}\geq 0$. Again, ${\varepsilon_{{\rm{rk1-b}}}}$ and ${\varepsilon_{{\rm{rk1-c}}}}$ are set to large values and are progressively reduced by 5 at each iteration, attaining ${\varepsilon_{{\rm{rk1-b}}}}=1.25\cdot 10^{-9}$, ${\varepsilon_{{\rm{rk1-c}}}}=1.27\cdot 10^{-7}$ at convergence.

As performance metric, we consider the reradiation pattern of the RIS, i.e., we plot the power flux $P_{{\theta}}\left(\boldsymbol{\gamma}\right)$ in (4) as a function of the angle of observation $\theta$, with $\boldsymbol{\gamma}$ obtained from the proposed optimization algorithms. Also, three benchmark schemes are considered: (i) GO is the geometric optics solution in [1, Eq. (102)], with non-negative values of the real part of the surface impedance in (1); (ii) GD is the global design solution in [1, Eq. (81)], with positive and negative values of the real part of the surface impedance in (1); and (iii) GO-RI is the GO solution in [1, Eq. (102)], by setting the real part of the surface impedance equal to zero. The initial values of the Taylor series approximations in S-RI-a and P-RI are set to GO-RI for ensuring the feasibility of the initial point.

Figure 1: Fig. 1: Reradiation pattern of the RIS: (left) S-RI-a and (right) P-RI

The results are illustrated in Fig. 1. According to theory [1], GO and GD provide reradiation patterns with no beams towards undesired directions of radiation, with GD offering the best beamforming gain at the highest implementation complexity (local amplifications due to the negative values of $\operatorname{Re}(\bf{z})$ in (1)). It is worth nothing that the difference of beamforming gain (at $\theta_{r}=60^{\circ}$) between GO and GD is 3 dB. GD-RI, which is often utilized as a simple solution [7], offers a good beamforming gain towards the direction of interest, but strong beams towards two undesired directions [5]. The proposed algorithms ensure a good, close to optimal, beamforming gain, no beams towards unwanted directions, and a real part of the surface impedance almost equal to zero, i.e., positive and less than $\varepsilon_{\rm{RI}}$. The problem P-RI usually needs higher memory requirements due to the SDP formulation (in matrix form), but it is more stable than the problem S-RI-a from the numerical point of view. Therefore, the proposed algorithms provide an efficient approach for optimization, ensuring close to optimal performance while fulfilling specified design constraints.

VI Conclusion

We have introduced a suite of algorithms for optimizing reconfigurable anomalous reflectors based on the global design criterion. Several extensions of this work can be envisioned, including multi-beam anomalous reflecting and refracting surfaces, multiple antenna systems, near-field RIS-aided channels.