Macroscopic Modeling of Anomalously Reflecting Metasurfaces: Angular Response and Far-Field Scattering

In view of extremely challenging requirements on the design and optimization of future mobile communication systems, researchers are considering the possibilities of creating intelligent radio environments by using reconfigurable and smart metasurfaces integrated into walls, ceilings, or facades. In this novel communication paradigm, tunable metasurfaces redirect incident waves into the desired directions. To design and characterize such smart radio environments in any realistic scenario, it is necessary to know how these metasurfaces behave when illuminated from other directions and how scattering from finite-sized anomalous reflectors can be estimated. In this work, we analyze the far-field scattering of reflective metasurfaces and study the angular response of anomalous reflectors for arbitrary illumination angles. Using the surface-impedance model, we explain the dependence of the reflection coefficients of phase-gradient metasurfaces on the illumination angle and present numerical examples for typical structures. We also consider scattering from finite-sized metasurfaces and define a route toward including the full-angle response of anomalous reflections into the ray-tracing models of the propagation channel. The developed models apply to other diffraction gratings of finite size.


I. INTRODUCTION
W ITH the development of mobile communications, there is a continuous challenge to improve the communications systems making them more efficient and versatile.New research directions point toward the use of intelligent radio environments supported by reconfigurable metasurfaces; see, e.g., [1]- [9].In contrast to conventional wireless communication systems where only transmitters and receivers can be optimized, these novel approaches will potentially enable the controllability, programmability, and optimization of the propagation channel.Metasurfaces can control the properties of reflected and scattered fields and, by embedding tunable and active elements in the constitutive meta-atoms, allow changing their electromagnetic response, optimizing the functionality.Among different types of metasurfaces, reflective metasurfaces, typically backed by a metallic plane that blocks transmitted waves, are good candidates to be integrated into flat or curved surfaces of walls, ceilings, building facades, etc.The most common approach to design metasurfaces for manipulating the direction of reflected waves is to make the reflection (Corresponding author: A. Díaz-Rubio) This work was supported in part by the European Commission through the H2020 ARIADNE project under grant 871464.
phase nonuniform along the metasurface.As it is known from the phased-array antenna theory, in the particular case of anomalous reflection (a plane wave is reflected breaking the reflection law, i.e., the reflection angle θ rd is not equal to the incidence angle θ id ), the reflection phase should linearly vary compensating the phase mismatch between the incident and reflected waves.The corresponding local reflection coefficient, defined as the ratio between tangential electric fields of the incident and reflected waves at each point of the metasurface, can be written as Γ(x) = exp[j(sin θ id − sin θ rd )kx].Here, k = ω √ µε is the wavenumber in the background medium, and axis x is directed along the tangential component of the wavevector.This reflectarray method (e.g., [10], [11]) is also known as the generalized reflection law in the optics community [12], and it allows controlling the direction of the reflected wave by engineering the reflection phase produced by each small area of the metasurface.
Fig. 1.Efficiency of anomalous reflectors characterized by a local reflection coefficient with a linear phase gradient.The efficiency is calculated as η = 4 cos θ i cos θr/(cos θ i + cos θr) 2 .Derivation of this formula and a detailed discussion can be found in [13].
Metasurfaces designed using this method perform highefficiency anomalous reflection when the difference between the incident and reflected angles is small or when the energy is sent back into the same direction as the incident wave (the retroreflection scenario, θ rd = −θ id ).This behavior is illustrated in Fig. 1, where one can see that a local reflection coefficient with a linear phase gradient ensures perfect performance around the specular and retroreflection directions.However, as it was explained and demonstrated in [14]- [16], the efficiency of phase-gradient anomalous reflectors with large transformations of the wave propagation direction drops considerably due to the impedance mismatch between the incident and reflected waves [13].Alternative, theoretically perfect solutions have been found using nonlocal structures [15], [17], diffraction gratings [18], auxiliary evanescent fields [19], or power-flow-conformal surfaces [20].These solutions offer possibilities for creation of extreme anomalous reflectors with high efficiency.The reflective properties of these advanced metasurfaces cannot be modelled by a local reflection coefficient varying over the reflecting plane.
Current knowledge about metasurface-based anomalous reflectors allows us to understand such important parameters as the efficiency and the energy distribution of the parasitic reflections for the design conditions.However, in order to characterize the propagation channel in environments with integrated metasurfaces, it is necessary to study the response of anomalous reflectors for angles of incidence that are different from the design incidence angle (θ i = θ id ).In other words, it is necessary to characterize the angular response of metasurfaces for arbitrary illuminations and provide effective design tools for wireless communication systems and other applications.This is necessary because in any realistic dynamic or multipath scenario metasurfaces are illuminated from many directions.
In this paper, we study this problem and make a systematic analysis of the angular response of phase-gradient anomalous reflectors designed using the phased-array principle (the generalized reflection law).
Another important problem that must be solved for further developments of anomalous reflector technology is the estimation of reflected fields from finite-sized metasurface panels.Current literature on metasurface reflectors focuses on the design and properties of infinite planar metasurfaces illuminated by plane waves, and this knowledge is not enough for finding fields that are reflected and scattered by panels of limited size.Known approaches to calculation of reflection from finite-sized anomalous reflectors are based on the local reflection coefficient model, that is, on the assumption that at each point of the metasurface the reflected field equals the incident field with the desired phase shift, e.g., [3], [7], [8], [21], [22].In some works, the amplitude of the reflected field is scaled to ensure that all the incident power is reflected to the desired direction [6], [9], corresponding to perfect operation.The control of the reflection phase is assumed to be achieved by engineering the local surface impedance of the reflector, using the locally periodic approximation.Similar assumptions are made in modeling finite-size metasurfaces for control of transmission [23], [24].
These models do not account for parasitic reflections due to impedance mismatch and the surface periodicity [14]- [16], [25].Another limitation is that the known models can be used only for a single plane-wave illumination at the design incidence angle, and they are not useful if the metasurface is illuminated from many directions (the multi-path scenario).Moreover, the local reflection coefficient model is tantamount to the assumption of locally specular reflection, which is approximately valid only for slowly varying reflecting properties of the metasurface (that is, only when the tilt angles are small).
Here, we propose an approximate analytical method for calculations of reflected fields in the far zone, merging approaches based on the physical optics and on the theory of diffraction gratings.This model takes into account parasitic, multi-beam scattering and is applicable for arbitrary illumination angles.We hope that by using this method it will become possible to extend the applicability of ray-tracing algorithms for estimation of propagation channel to engineered environments with reconfigurable intelligent metasurfaces (RIS).From the more general perspective, the developed theory allows effective calculation of far-zone scattered fields for general periodical structures that operate as diffraction gratings.

II. PROBLEM STATEMENT AND THE DEFINITION OF MODE CONVERSION EFFICIENCY
The first goal of this study is to analyze the response of anomalous reflectors of infinite extent when they are illuminated from an arbitrary direction, so that the illumination angle, θ i , can be different from the angle of incidence for that the metasurface has been designed or tuned, θ id .Our intention is to explore the response of these devices in the general way, without considering specific phenomena associated with a particular implementation.For this reason, we model the metasurfaces as boundaries whose surface impedance varies along the reflection plane [15], [20].The surface impedance is defined locally as the ratio between the tangential components of the surface-averaged electric and magnetic fields.The averaging is made on the scale of the unit cells of the metasurface.
Here, for simplicity, we consider planar metasurfaces whose properties vary only along one direction.Axis x of a Cartesian coordinate system is in the metasurface plane and points along the direction of impedance variation.Axis y is normal to the metasurface plane.Such metasurfaces can tilt the reflected plane wave direction in the xy-plane.The relation between the local reflection coefficient Γ, dictated by the phased-array principle (the generalized reflection law), and the surface impedance Z s can be written as where Z w is the characteristic impedance of the incident plane wave -the ratio between the tangential to the reflector plane components of the electric and magnetic fields of the incident wave.The impedance of the incident plane wave can be written as Z (TE) w0 = Z 0 / cos θ i for TE-polarized waves and Z (TM) w0 = Z 0 cos θ i for TMpolarized waves.Here, Z 0 is the wave impedance of the surrounding medium (usually air).From this expression, it is straightforward to find the surface impedance of phasegradient anomalous reflectors as The surface impedance is a periodic function of the coordinate x, Z s (x) = Z s (x + D), whose period can be found as D = λ/| sin θ id − sin θ rd |, e.g., [15].This period defines the scattering properties of the metasurface and has important influence on the angular response of periodic metasurfaces.As it is known from the Floquet theory, periodical variations of metasurface properties lead to existence of multiple diffracted modes which depend on the illumination angle, the period, and the frequency.Due to the periodicity of the surface impedance, the tangential wave number of different diffracted modes must satisfy where n is the order of the diffracted mode.With this notation, mode n = 0 corresponds to the specular reflection.Consequently, the normal component of the wave vector is It is important to note that only when the value of n ensures that k > k xn , the diffracted mode of index n propagates into the far zone.For propagating modes, the angle of reflection reads The incident power is distributed between the propagating modes depending on the surface impedance of the metasurface.
As we can see from Eq. ( 3), the propagation direction of diffracted modes strongly depends on the illumination angle, and for different illumination angles the propagation directions of reflected waves change.Moreover, this dependence is different for different diffracted modes.The use of the surface impedance model allows us to numerically study the power distribution among all diffracted modes using the mode-matching method.Recently, this method has been used for the analysis of metasurfaces, allowing fast calculations of scattering parameters.For a planar metasurface laying in the xz-plane, we can write the tangential components of the electric and magnetic fields in form of the Floquet series: A n e −jkxnx (4) where A n represents the amplitude coefficient for the n-th reflected harmonic.Notice that the amplitude coefficients are normalized to the amplitude of the tangential component of the incident field.The '∓' sign in Eq. ( 5) is negative for TEpolarized waves and positive for TM-polarized waves.The wave admittance for the diffracted modes can be calculated as Y

(TE) wn
= k zn /ωµ 0 for TE-polarized waves and Y (TM) wn = ωε 0 /k zn for TM-polarized waves.By enforcing the impedance boundary condition where ŷ is the unit normal vector, we determine the amplitudes of all Floquet harmonics.Full description of this method can be found in [26] and a MATLAB code that implements this mode-matching method for metasurfaces modeled by space-varying surface impedance can be found in [27].
Knowing the amplitudes of all the reflected waves, one can easily calculate the conversion efficiency for each propagating plane-wave mode, defined as the ratio between the normal component of the reflected power curried by a certain mode and the normal component of the incident power: η n = P zn /P in z , where P zn = (1/2) Re(E tn × H * tn ), and P in z is the amplitude of the normal component of the Poynting vector of the incident plane wave.For TE-polarized incidence, the power carried away from the metasurface (per unit area) by each Floquet harmonic can be expressed in terms of the corresponding mode-conversion efficiencies for propagating modes 0 for evanescent modes (6) The conversion efficiency for each diffracted mode for TMpolarized waves is defined as for propagating modes 0 for evanescent modes (7) For lossless metasurfaces the sum of the conversion efficiencies of all propagating modes equals unity.

III. ANGULAR PROPERTIES OF ANOMALOUS REFLECTORS
Before continuing with the analysis of the angular response, it is necessary to define the surface impedance for different illumination angles.In this study, for simplicity, we will assume that the surface impedance remains constant for any illumination angle.This means that we neglect spatial dispersion of the metasurface response.This assumption is approximately valid if we choose a proper topology for the actual implementation of the metasurface.For example, one can consider metasurfaces based on metallic grooves acting as local phase shifters or metallic patterns over a highpermittivity dielectric layer backed by a metallic plane.More details about the effect of the topology on the angular response will be presented in Section VI.
We start the study of the angular response with a phasegradient anomalous reflector designed for θ id = 0 • and θ rd = 70 • .This is a typical example of a metasurface designed for a large tilt of the reflected wave direction (with respect to the usual specular-reflection angle).The period of this reflector is D = 1.0642λ.The first step in the analysis is to determine what diffracted modes propagate in the system for different illumination angles.Figure 2(a) shows the directions of propagation of the diffracted modes for different illumination angles.From these results we can see that only diffracted modes with the orders n ∈ [−2, 2] contribute to the angular response in the far zone (even though many other modes can be exited, they do not propagate away from the surface and contribute only to the reactive fields in the vicinity of the reflecting surface).
Next, we study how the energy is distributed between all these propagating modes for different illumination angles.The energy distribution for different illumination angles is represented in Fig. 2(b).For the design conditions, i.e., for the normal illumination, most of the incident power (76%) goes to the first diffracted mode n = 1 [see circular symbols in Fig. 2], realizing the desired anomalous reflection.For this illumination angle, as it is expected from the design conditions, this mode propagates at θ r = 70 • .The rest of the incident power is scattered toward θ r = −70 • direction (mode n = 1) and in the normal direction (mode n = 0).These results agree with the previous analysis of anomalous reflectors done in [13], [15].The energy distribution in Fig. 2(b) shows that when the angle of incidence θ i ∈ [−70, 70] • , most of the incident energy is transferred to higher-order modes (n = 0).It is also important and instructive to consider the case when this metasurface is illuminated from θ i = −28 • (see squared symbols in Fig. 2).For this specific illumination angle, all the energy is efficiently transferred to the diffracted mode n = 1.To understand this behavior, one needs to notice that at this angle of incidence, the propagation direction of this diffracted mode corresponds to the retroreflection scenario.As we can see from Fig. 1, retroreflection is the only nontrivial transformation of the reflected waves that produces perfect performance when we use a linear phase gradient of the reflection phase.Equation (1) can be written for the retroreflection case as Z s (x) = jZ w cot (sin θ retro kx), with θ retro being the angle of retroreflection.As it has been demonstrated in [28], the boundary condition for this surface impedance is perfectly satisfied by the sum of only two plane waves: the incident and the retroreflected ones, ensuring perfect performance.By a simple comparison of this equation and Eq. ( 1), we can see that the retroreflection angle and the reflection angle in the design of the anomalous reflector are related as For the specific case when θ id = 0 • and θ rd = 70 • , the retroreflection angle calculated using this expression equals θ retro = −28 • .As under our assumptions the surface impedance remains the same for different illumination angles, an anomalous reflector designed for θ id = 0 • and θ rd = 70 • performs perfect retroreflection when the surface is illuminated at Another interesting aspect of the angular response appears when the illumination comes from θ i = 28 • (see triangular symbols in Fig. 2).Here, we can see that all the incident power is transfer to the diffracted mode n = −1, meaning that all the power is sent back into the same direction.In other words, an anomalous reflector designed for θ id = 0 • and θ rd = 70 • also performs perfect retroreflection when θ i = 28 • .The explanation of this behaviour is simple.As the metasurface is reciprocal, specular reflection phenomenon is symmetric with respect to the illumination angle sign [see the green line in Fig. 2(b)].This means that, as no power is reflected specularly when the metasurface is illuminated from θ i = −28 • , no power will be sent into the specular direction when θ i = 28 • .The period of this metasurface determines that only one diffracted mode propagates when θ i = ±28 • and, considering that specular reflection is not allowed by reciprocity, the power must go to the diffracted mode n = −1.
The second case of study is an anomalous reflector that is designed for θ id = 0 • and θ rd = 20 • .This is a typical example of a small tilt of reflection direction, when the impedance is a slowly varying periodical function of x, and the period is electrically large.In this specific example, the period of the metasurface is D = 2.9λ.As in the previous case, to understand the scattering properties of this metasurface, it is necessary to consider what diffracted modes propagate in the system for different incident angles.Figure 3(a) shows the reflection angles for different illuminations and for each propagating mode.The most evident difference with the previous example is that, due to the larger period, there are more propagating modes.We can also see that the dependence of the reflection angle on the incident angle becomes more linear.It is particularly interesting to see that for the first diffracted mode, n = 1, the dispersion can be approximated by the linear relationship where θ retro = −9.5 • .This linear approximation is shown in Fig. 3(a) as a dashed line.
The reflected power distribution is represented in Fig. 3(b) for different illumination angles showing important differences in comparison to the previous example.Here, we can see that the efficiency of the anomalous reflector for the design conditions (see circular symbols) increases, and nearly all the power illuminating the metasurface from normal direction goes into the desired direction, θ r = 20 • .However, when the metasurface is illuminated from the retroreflection angle, θ i = −9.5 • , all the power is sent back (see triangular symbols).An important difference appears when we illuminate the metasurface from the opposite direction, θ i = 9.5 • .In this case, all the power is still transferred to the diffracted mode n = 1 that sends the power into θ r = 31.5• (see squared symbols).With this configuration, the existence of multiple diffracted modes allows satisfying the reciprocity condition of the system without enforcing perfect retroreflection for this illumination.From this result, we can conclude that only in the cases where two modes propagate in the system (specular reflection and an additional diffracted mode), retroreflection is enforced by reciprocity when θ i = −θ retro .Figure 4 presents the results of the analysis of the number of propagating modes for anomalous reflectors designed for the normal illumination.
Here we can see that perfect retroreflection due to reciprocity only takes place when θ rd > 42 • .
It has been demonstrated that the first diffracted mode, n = 1, significantly contributes to scattering the impinging power.The angular dispersion of this mode and the amount of power carried by this mode depend on the direction of illumination.Figure 5(a) presents the angular dispersion of the diffracted mode for different illumination angles and different design conditions.Here, we can see that keeping θ id = 0, the angular dispersion curve of the mode becomes wider and more linear.Importantly, the impinging power predominately scatters into the first mode when the angular dispersion is approximately linear, as we can see in Fig. 5(b).This is because the linear phase gradient of the reflection coefficient introduces the proper additional linear momentum to produce this transformation.We see that phase-gradient anomalous reflectors function as diffraction gratings, scattering incident power into all prop-agating modes of the Floquet expansion.The linear phase gradient optimized for a specific pair of the incidence and reflection angles ensures maximization of the power reflected into the desired direction, but ideal conversion efficiency for phase-gradient metasurfaces is possible only for retrodirection and for trivial specular reflection.When the incidence angle changes, reflective properties change dramatically (from highefficiency anomalous reflection to retroreflection).

IV. EFFECTIVE REFLECTION COEFFICIENTS AND FAR-FIELD SCATTERIG
As we have seen in the previous sections, for the analysis and design of phase-gradient metasurfaces, it is convenient to use the local reflection coefficient, defined as the ratio between the tangential component of the incident and reflected electric fields at each point of the metasurface.This reflection coefficient models local properties of the reflector (averaging is made on the scale of one meta-atom).However, this parameter does not allow proper modeling of reflection from the metasurface as a whole and understanding of the effects on the propagation in environments with metasurfaces.The main reason is that gradient metasurfaces are diffraction gratings which scatter into several directions, as was discussed in the preceding section.In this case it is more convenient to define the reflection coefficient of infinite metasurfaces as the ratio between the reflected electric fields of the plane-wave Floquet harmonics and the incident electric field, or, in other words, characterize the macroscopic behavior of the metasurface.
Typically, if we consider plane-wave reflection, the reflection coefficient is defined as the ratio between the tangential components of the electric field of the incident plane wave and the reflected plane wave.In the case of metasurfaces, the existence of multiple reflected plane waves does not allow us to use this definition and introduce a unique reflection coefficient.For this reason, a macroscopic reflection coefficient for modeling metasurfaces can be defined as a combination of the individual reflection coefficients for each individual reflected mode as R(θ i , x) = n r n (θ i )e −jk(sin θrn−sin θi)x e −jk(cos θrn+cos θi)y (10) where r n (θ i ) = |r n (θ i )| exp [jφ n (θ i )] is the individual reflection coefficient for each propagating diffracted mode.The angles θ rn defined by are the reflection angles of the propagating harmonics of index n.Using Eq. ( 4), we can define the value of the reflection coefficient using the amplitude of the reflected field, r n (θ i ) = A n .Here, the reflection coefficient is defined for an arbitrary distance from the metasurface.At the reflecting surface, we let y = 0.
Here, we consider only the propagating Floquet harmonics.Only these modes contribute to the macroscopic, farfield response of metasurfaces of infinite extent (the far-field reflection coefficient for evanescent fields is zero).Evanescent harmonics of the fields contribute to reflections from finitesized metasurfaces due to edge scattering effects, and we will discuss this issue later.We note that macroscopic reflection coefficients in the assumption that the metasurface ideally reflects only one plane wave were discussed earlier in [9].This analysis is a useful tool when we consider infinite metasurfaces that extend over a whole plane.However, for realistic scenarios with finite-size metasurfaces we need to consider the scattering pattern of the metasurface.If we assume that the metasurfaces under study are large compared to the free-space wavelength and that the unit-cell size is small compared to the wavelengths of field variations along the metasurface plane λ n = 2π/k xn for all significant modes, we can use the physical-optics approximation.Namely, we can assume that the currents induced on the finite-size metasurface are the same as on the infinite metasurface.That is, we neglect the perturbations of induced currents near the edges of the metasurface.Importantly, note that we do not assume a local relation between the incident field and the induced currents, which is the other conventional assumption of the physical optics, see, e.g., [29].Instead, we find the induced currents solving the reflection problem for a periodically modulated impedance boundary, which results in the knowledge of the macroscopic reflection coefficients r n .
For simplicity of writing, we assume that the metasurface has a rectangular shape, occupying the area defined by −a < x < a and −b < z < b with 2a and 2b being the total sizes of the metasurface.The geometry of the problem is shown in Fig. 6.If the metasurface is illuminated by a TE-polarized wave defined by the electromagnetic field E i = E 0 e −jk(sin θix−cos θiny) ẑ, the reflected electric field at a certain point at the metasurface plane (position vector r with the coordinates x , z ) reads Here, we use the physical-optics approximation, assuming that the reflected field is not disturbed near the metasurface edges.
The components of the magnetic field associated with the reflected electric field can be written as with Y 0 = k/ωµ being the admittance of the background media.Note that this relation is written for plane waves propagating along the directions defined by the reflection angles θ rn , because we do not assume locally specular reflection.
Using the Huygens principle, we introduce equivalent electric and magnetic surface currents that create these reflected fields.The corresponding surface current densities read j er = ŷ × H r and j mr = −ŷ × E r .On the shadow side we assume, within the physical optics approximation, that the total field is zero, since the metasurface is impenetrable.Thus, the fields scattered by the metasurface in the shadow direction are equal to the negative of the incident fields.The corresponding equivalent Huygens shadow surface current densities read j es = −y × (−H i ) and j ms = ŷ × (−E i ).Since the electrical thickness of the metasurface is negligible, we can sum up these two Huygens currents and assume that both flow at the same surface.The total Huygens current densities correspond to the total fields at the illuminated side: In earlier work [6], it was assumed that phase-gradient reflectors can be viewed as infinitely thin sheets of only electric surface currents with appropriately modulated phase, but that assumption contradicts to the assumption of an impenetrable reflector.Indeed, if the metasurface carries only an electric surface current, the tangential electric field is continuous across the metasurface.Thus, if the field behind the metasurface is zero, the only allowed value of the reflection coefficient is −1, and no phase modulation is possible.For this reason, consistent models require introduction of both electric and magnetic surface currents.Magnetic surface current density can be expressed in terms of the electric current using the characteristic impedance of the corresponding reflected mode.Note that the conventional physical optics expressions (for example, for calculations of far fields of aperture antennas) assume the free-space impedance relation between the electric and magnetic fields at the radiating surface, which in our case corresponds to the assumption of locally specular reflection.This assumption is valid only for slowly varying surface impedances, and its use cannot be justified for majority of anomalously reflecting metasurfaces.
For the equivalent Huygens surface currents in form ( 15), the reflected electric field at any point in space (position vector r) is given by the following integral (e.g., [29,Eq. (3.87)]) where is the scalar Green function, and the integration is over the metasurface plane.This formula can be used for calculations of scattered field at any distance from the metasurface, as long as it is large compared to the unit-cell (meta-atom) size.As is well known from the theory of diffraction, calculations of fields in the far zone can be significantly simplified.In this work, we use the simplified far-field model that is valid under the following assumptions: |r| λ, |r| L, and L 2 /|r| λ.Here, L = max(2a, 2b) is the largest size of the metasurface plate.That is, the distance to the observation point is large compared to the wavelength and to the size of the metasurface plate.Under these assumptions, the reflected field can be found as 18) where Z 0 = 1/Y 0 is the impedance of the surrounding space and r is the unit vector along r, pointing from the center of the metasurface plate to the observation point.Details of the derivation of this formula and justifications of the simplifying assumptions can be found in [29,Sec. 3.4.1].
We see that the far-zone field is a spherical wave, whose fields are orthogonal to vector r (to the direction of scattering).The amplitude and angular dependence of the scattered field are given by the surface integral that can be viewed as the spatial Fourier transform of the surface current distribution.Calculation of this integral is trivial because the integrand is a sum of simple exponential functions.The result is obviously a sum of products of two sinc functions (sinc(x) = sin(x)/x).
It is convenient to express the final result using spherical coordinates, writing r = sin θ cos φẑ + sin θ sin φx + cos θŷ.Let us consider the scattered fields in the incidence plane (φ = π/2).In this case, the θ-component of the reflected field is zero, one of the two sinc functions equals unity, and the expression for the z-component of the scattered electric field can be simplified as where a ef = (sin θ−sin θ i )a, and a efn = (sin θ−sin θ rn )a represents the effective size of the metasurface for each reflected propagating mode.Furthermore, S = 4ab is the area of the metasurface.As a check, for a PEC plate r 0 = −1 for any θ i , and the other terms in the summation equal zero.Furthermore, the reflection angle θ r0 = θ i (specular reflection).In this trivial case, formula (19) reduces to the well-know result of the physical-optics model of scattering from electrically large PEC plates (e.g., Eq. ( 4) in [6]).The first term in Eq. ( 19) represents the contribution of the currents generating shadow radiation.In the ideal case of an infinite impenetrable metasurface, the shadow radiation will be zero in the backscattering direction, and it will ensure complete cancellation of the incident field in the forward direction.However, due to the finite size of actual metasurfaces, the contribution of the shadow radiation to the total radiation should be included because these currents create fields in the whole space.For a clear representation, it is convenient to normalize the pattern so that its maximum value for the PEC plate of the same size is unity.With this normalization, the patterns of the reflected and the shadow fields can be written as and the total scattering pattern is defined as F sc (θ) = F r (θ) + F sh (θ). Figure 7 shows the patterns for the reflected, shadow, and total scattered fields produced by an anomalous reflector designed for θ id = 0 • and θ rd = 10 • under normal illumination.The size of the metasurface is assumed to be 10λ × 10λ.
As expected, total backscattering radiation is predominately in θ = 10 • and additional secondary lobes appear due to the finite size of the metasurface.From the figure, it is clear that the contribution of the shadow radiation to the total backscattering is only noticeable in the side lobes, where it is significant.We note that the reflection pattern of finite-sized anomalous reflectors was recently considered in [22].However, the model in that paper does not take into account the diffraction-grating nature of such metasurfaces (only one reflected plane-wave harmonic is assumed, defined by a local reflection coefficient) as well as the shadow-current contribution to scattering; compare with Eq. (32) in [22].
Formula ( 19) is general, applicable for both propagating and evanescent harmonics of the reflected fields.Indeed, higher-order, evanescent components of the equivalent currents are also exponential functions of x , and integrals of these components also give sinc functions of the same form as for the contributions from the propagating Floquet harmonics.The difference is that for the evanescent modes the arguments of these sinc functions are large compared to unity.Obviously, for the propagating harmonics the value of | sin θ rn | is smaller than unity, and the argument of the sinc function can take zero value, corresponding to the main-lobe maximum.However, for evanescent harmonics of order n, the argument of the sinc function does not cross zero for any observation angle θ because sin θ i + n λ D > 1.Thus, contributions of evanescent modes of the induced surface currents to scattering in the main-lobe directions scale as 1/(ka) as compared to the contributions of the propagating harmonics.We also note that, for high-order modes, n λ D is a large number as compared to unity.Thus, the corresponding sinc functions weakly depend on the angles of incidence and observation.This means that contributions from evanescent current harmonics affect mainly the side-lobe level and shape, and they are predominantly of diffuse character.Finally we note that, within the physical optics approximation, we can take these modes into account only if the variations of the induced currents over the metasurface plane remain slow on the micro-structure scale, see Section VI.This assumption cannot be justified in the frame of the surface-impedance model.
In all numerical examples presented in this paper we include the contributions of propagating field harmonics only, neglecting the small effects of additional diffuse scattering due to edge inhomogeneities of the evanescent fields.

V. ANGULAR BANDWIDTH OF ANOMALOUS REFLECTORS
Anomalous reflectors are designed to change the direction of reflected waves breaking the reflection law, θ id = θ rd .In a sense, for the design conditions, an anomalous reflector acts as a virtual tilted mirror rotated by the angle θ t = (θ id + θ rd )/2.Under an arbitrary illumination, a tilted mirror reflects the incident waves to the angle θ r = θ i + 2θ t .This linear relation can be understood as the angular dispersion of the ideal tilted mirror.
Because of the functional analogy between the two problems, one can define the angular bandwidth of an anomalous reflector as the range of incident angles within which it behaves as an ideal tilted mirror.As we have learnt from the analysis in Section III, anomalous reflectors predominantly couple the energy into the first diffracted mode n = 1.So, in order to analyze the bandwidth of anomalous reflectors we should focus on the angular dispersion of the first diffracted mode.
Figure 8(a) reports the angular dispersion of the first-order diffracted mode compared to the linear dispersion of the equivalent tilted mirror for two different configurations.First, let us consider an anomalous reflector designed to produce a small transformation of the incident waves, θ id = 0 • and θ rd = 10 • .In this case, we can see that the angular dispersion of the diffracted modes can be approximated by the linear Using this definition, the field scattered by the metasurface when θ i = 0 • , 15 • , 30 • is represented in Fig. 8(b).From this analysis one can clearly see that, as far as the angular dispersion of the first-order diffracted mode is linear, the metasurface behaves as an efficient virtually tilted mirror.
Angular dispersion of the diffracted modes strongly depends on the periodicity of the system that is directly related with the desired transformation to be performed by the anomalous reflector.As we can see from the second example shown in Fig. 8(a), where θ id = 0 • and θ rd = 80 • , the dispersion of the first-order diffracted mode cannot be approximated by a linear relation.In this example, the anomalous reflector will only behave as a virtually tilted mirror for the design conditions and the reciprocal scenario.Figure 8(c) shows the scattering properties of this anomalous reflector when θ i = 0 • , 15 • , 30 • for a square-shaped metasurface of the same size (10λ × 10λ).Here, we can see a completely different behaviour than in the previous example.Even a small deviation of the incident angle from the design value leads to a dramatic change in reflector performance, which becomes totally different from the anomalous reflection with the desired tilt by 80 • .This result confirms the strong dependency of the angular bandwidth on the design conditions of anomalous reflectors.

VI. INFLUENCE OF THE ACTUAL TOPOLOGY TO THE ANGULAR RESPONSE
In the preceding analysis, we assumed that the surface impedance remained unchanged for different illuminations angles.As we mentioned before, this assumption is only valid for some actual structures and cannot be considered as a general rule.In this section, we present an example of a topology that satisfies this condition and validates this assumption.
In particular, we will numerically confirm that one can implement a surface impedance dictated by Eq. ( 1) without angular dispersion using closed-end grooves in a metal plate [see in Fig. 9(a)].These meta-atoms have been proposed to implement power flow-conformal metasurfaces [30].We assume that the groves are infinitely long and uniform and the thickness of the walls is much smaller than the width of the groves, so that we can neglect the effects of wall thickness on the surface impedance.For TM-polarized electromagnetic waves, the response of the meta-atom is defined by the surface impedance Z EM s = jZ 0 tan(k ), with being the depth of the groove.As an example, we will design and numerically study an anomalous reflector targeted for θ id = 0 • and θ rd = 40 • .The operation frequency of the metasurface is f = 8 GHz.The surface impedance required to implement this functionality is represented in Fig. 9 To make the metasurface thickness small compared to the free-space wavelength, the grooves can be filled with a highpermittivity dielectric.In this case, each unit cell can be be viewed as a combination of an electric-current element (due to the current at the bottom termination) and a magnetic-current element (due to the out-of-phase currents at the vertical walls).Effective induced magnetic current is absent only in the case of zero depth of the groove.As was discussed in Section IV, this obviously corresponds to the reflection phase 180 • .
The numerical study of the structure is done using the commercial software Ansys HFSS.The response of the structure is represented in Fig. 9(c) by the circular symbols.For comparison, the angular response for a metasurface implementing the continuous surface impedance and calculated using the Floquet harmonic expansion has been included in the figure (see solid lines).We can see excellent agreement between both results verifying the angular independence of the implemented surface profile.We note that for other realizations of phasegradient reflectors (most commonly, as planar patch arrays at a grounded dielectric substrate when the substrate permittivity is not large) the surface impedance depends on the incident angle, which makes the angular dependence of the metasurface response more complicated.

VII. CONCLUSIONS
In this work we have investigated the angular response of phase-gradient metasurfaces designed for manipulating the direction of propagation of reflected electromagnetic waves.Due to their extraordinary properties for wave manipulation and possibilities for simple integration in real-world environments, characterizing the scattering properties of these metasurfaces under different illumination conditions is necessary for the development of smart radio environments and other applications of advanced metasurfaces.
Using the surface impedance model, we have studied the response of anomalous reflectors for arbitrary illumination angles.The results show that for anomalous reflectors designed using the linear phase gradient of the local reflection coefficient, the incident power is predominantly coupled to the diffracted mode n = 1.Only when the design requirements demand strong transformations of the direction of propagation, significant portion of the power is reflected specularly.We have also delved into the relation of the number of diffracted modes and the angular response for reciprocal metasurfaces.We have verified the results obtained with the surface impedance model by designing an actual structure that confirms the angular response predicted by the impedance model.
We have proposed an analytical modeling approach to predict the far-zone reflective and scattering responses of finitesize anomalous reflectors for arbitrary illumination angles.This model can be used as the key element in future analytical and numerical models of radio wave propagation in complex environments including metasurface panels on the surfaces of walls or facades.
To find the reflected fields in the far zone, we merge some ideas of the physical optics (but using macroscopic, essentially nonlocal expressions for the induced currents) and the approaches of the theory of diffraction by gratings (because anomalous reflectors are periodical structures that can create multiple diffracted beams).The developed method is based on the use of introduced macroscopic reflection coefficients that relate the amplitude of the incident wave and the amplitudes and phases of the fields of each diffracted mode.Using the developed formalism, it is simple to calculate the far-field response of the metasurface when it is illuminated from any direction and estimate the link budget for arbitrary positions of the transmitter and receiver.
The developed analytical formulas for far-field scattering from finite-size metasurfaces are valid not only for phasegradient reflectors but also for nonlocal metasurfaces that offer theoretically perfect anomalous reflection.In those advanced designs [15], [17]- [20], the amplitude of reflection from the design direction θ i = θ id into the desired anomalously reflected mode is close to the ideal value defined by Eqs. ( 6) and ( 7) for 100% efficiency.The amplitudes of all other propagating Floquet harmonics are suppressed to negligible levels.In this case, practically all incident power is reflected as a single plane wave propagating in the desired direction, if the metasurface is infinite.For other angles of incidence, there is scatterings into all propagating modes, and the corresponding amplitudes need to be found numerically or experimentally.Importantly, the phase of the main reflected harmonic can be engineered, which opens up a possibility to control the side-lobe structure of scattering from finite-size metasurfaces because scattering into side lobes depends on interference with the shadow radiation, as illustrated in Fig. 7.

Fig. 2 .
Fig. 2. Analysis of the scattering properties of a flat infinite anomalous reflector designed for θ id = 0 • and θ rd = 70 • .(a) For a periodicity given by D = λ/| sin θ id − sin θ rd |, relation between the angle of incidence θ i and the angle of reflection of the diffracted modes θr.(b) Energy distribution between the reflected propagating harmonics for different incident angles.The power is normalized to the amplitude of the normal component of the incident power flow.

Fig. 3 .
Fig. 3. Analysis of the scattering properties of an infinite anomalous reflector designed for θ id = 0 • and θ rd = 20 • .(a) For a periodicity given by D = λ/| sin θ id − sin θ rd |, relation between the angle of incidence θ i and the angle of reflection of the diffracted modes θr.The dashed line corresponds to the linear approximation (9).(b) Power distribution between the propagating harmonics for different incident angles.

Fig. 4 .
Fig. 4. The number of propagating diffracted modes for anomalous reflectors designed for θ id .Red line marks the incident angles for perfect retroreflection θ i = θretro.Dashed red line marks the incident angles for the reciprocal scenario of perfect retroreflection θ i = −θretro.

Fig. 5 .
Fig. 5. Angular response of the first diffracted mode, n = 1, for different θ rd .(a) Angular dispersion for different design conditions.(b) Normalized scattered power by the diffracted mode.

Fig. 6 .
Fig.6.Geometry of the problem.A metasurface of rectangular shape and a finite size (2a × 2b) is placed on the xz-plane.The illumination is considered to be a TE-polarized plane wave, and the incidence plane is the xy-plane.

Fig. 7 .
Fig. 7. Analysis of the scattering properties of an anomalous reflector designed for θ id = 0 • and θ rd = 10 • under normal illumination.The radiation patterns are represented in dB: F i (dB) = 20 log 10 (|F i (θ)|), where i = sc, r, sh refer, respectively, to the total scattered field, the field scattered by the equivalent surface currents defined by the reflected fields Er, Hr, and the shadow currents defined by the incident fields.
(b) with a solid line.To have a fine surface profile and avoid errors due to the discretization of the continuous surface impedance, we use 15 elements per period to implement the surface impedance.The depths of the grooves to produce the required surface impedance are 10.625 mm, 11.875 mm, 13.125 mm, 14.375 mm, 15.625 mm, 16.875 mm, 18.125 mm, 0.625 mm, 1.875 mm, 3.125 mm, 4.375 mm, 5.625 mm, 6.875 mm, 8.125 mm, and 9.375 mm.The implemented surface impedance is represented in Fig. 9(b) with a dashed line.

Fig. 9 .
Fig. 9. Analysis of the scattering properties of an infinite anomalous reflector designed for θ id = 0 • and θ rd = 40 • and implemented using metallic grooves.(a) Normalized surface impedance for the anomalous reflector.Solid line represents the continuous surface impedance dictated by Eq. 1. Line with symbols represent the surface impedance implemented by the metallic grooves.Inset picture represent the scheme of the structure that implements the anomalous reflector.(b) Energy distribution between the propagating harmonics for different incident angles.Solid lines represent the results considering a continuous surface impedance.Symbols represent the results obtained from the numerical simulation of the actual structure.