Near Field Models of Spatially-Fed Planar Arrays and Their Application to Multi-Frequency Direct Layout Optimization for mm-Wave 5G New Radio Indoor Network Coverage

Two near field models for the analysis of spatially fed planar array antennas are presented, compared and applied to a multi-frequency wideband direct layout optimization for mm-Wave 5G new radio (NR) indoor network coverage. One model is based on the direct application of the radiation equations directly derived from the A and F vector potentials. The second model is based on the superposition of far field contributions of all array elements, which are modelled as rectangular apertures with constant field. Despite the different assumptions made to develop both models, the degree of agreement between them in the computation of the radiated near field is very high. The relative error between the models is equal or lower than 3.2% at a plane 13λ from the array, and it decreases as the near field is computed further away from the array. Then, the faster model is employed in a general direct layout optimization procedure to shape the electromagnetic near field for application in an indoor femtocell to provide coverage with constant power in a private office. Results show that a magnitude ripple better than 1.5 dB can be achieved in an enlarged coverage area covering the whole n257 band of the 5G NR, corresponding to a 10.7% relative bandwidth.


Introduction
Wireless communications within the 3rd Generation Partnership Project (3GPP) 5G New Radio (NR) paradigm have introduced new concepts and technologies, including ultradense networking, massive machine communications, all-spectrum access, full-duplexing, and massive MIMO [1,2], for applications in the industry 4.0 [3], smart cities [4], Internet of Things (IoT) [5], precision agriculture [6], and self-driving cars [7], among others [8]. These applications require higher data rates, higher channel capacity, and lower latency than those achieved with previous generations of mobile networks [9]. To that end, the 3GPP 5G NR technical specification [10] identifies two frequency ranges (FR) for the implementation of the 5G networks, namely, FR1 (0. 41-7.125 GHz) and FR2 , each of them divided into several bands. From the point of view of frequency allocation, the FR1 spectrum is already very crowded, by now providing services that range from cellular and satellite communications, to wireless local area networks [9]. Thus, the mm-wave ranges provided by the FR2, with their readily available wide spectrum, have become increasingly interesting for the deployment of 5G NR networks. In particular, there is special interest in the three lower bands in FR2, n257 (26.5-29.5 GHz), n258 , and n261 (27.5-28.35 GHz) since they have been released by the most countries, including the United States, the European Union, China, Japan, and South Korea [9]. However, the use of frequencies in the FR2 presents some challenges, such as path loss due to higher losses and shadowing due to blockage by physical objects that act as barriers for the propagating wave [11]. To overcome these issues, the use of femtocells with increased line of sight connectivity [12,13] has been proposed as a means of enhancing the cellular coverage while offloading the macro cell traffic in places with higher density.
Spatially fed arrays have been proposed within the intelligent reflective surface (IRS) paradigm [14][15][16], specifically reflectarray antennas for both indoor [17] and outdoor scenarios [18]. Transmitarrays have also been proposed as a solution to provide coverage in indoor femtocells [19,20]. In order to improve the coverage in indoor scenarios for transfer of information, power, or both simultaneously, near field beam-shaping techniques for the enlargement of coverage areas have been proposed. For instance, in [17] a phaseonly synthesis (POS) based on the generalized intersection approach [21] was carried out to design two reflectarray antennas, one providing coverage in an outdoor scenario by considering a squared-cosecant beam in elevation and sectorial beam in azimuth, while the other reflectarray provides coverage in an indoor scenario. In [19], a novel near field POS algorithm was developed in which the far field resulting from the truncated near field according to specifications is computed as an intermediate step before recovering the field at the aperture by means of the fast Fourier transform (FFT). It was employed to design a 3D-printed dielectric transmitarray antenna that provides a very narrow coverage area. This algorithm was later extended [22] to include simultaneous near and far field constraints. Another near field POS is carried out in [20] to enlarge the coverage area provided by a multi-fed transmitarray antenna.
The techniques employed in the previous works are very efficient for the near field beam-shaping of spatially fed arrays at a single frequency. However, they present some shortcomings. One common characteristic of the previous works is that they perform a POS, where a phase-shift for each array element is obtained at a single frequency such that the array radiates the desired near field. From those phases, the chosen unit cell is tuned in such a way that it produces the required phase-shift [23]. One obvious limitation of such an approach is the narrow bandwidth nature of the obtained arrays. Indeed, since the POS and layout design only consider a single frequency, the in-band behaviour of the array rapidly deteriorates [24]. Moreover, due to the intermediate step of computing the far field in the near field synthesis of [19,22], that algorithm is restricted to a single near field plane that is parallel to the array aperture.
In order to overcome those limitations, this work proposes a general near field beamshaping technique for the wideband optimization of spatially fed arrays. Specifically, the new analysis technique allows computation of the near field in a volume with an arbitrary local coordinate system rotated by three angles with regard to the array, as well as to perform a wideband direct layout optimization by taking into account the full electromagnetic response of the unit cell at several frequencies as provided by the use of a full-wave electromagnetic tool assuming local periodicity directly in the optimization loop. To that end, first two near field models for spatially fed arrays are presented. One model is based on the radiation equations directly derived from the A and F vector potentials [25]. The volumetric integrals are particularized to planar arrays discretized into a number of constituent elements. The resulting surface integrals are evaluated by means of the mid point quadrature, evaluating their convergence for different discretization grids. The second model is based on the superposition of the far fields radiated by each array element, which are modelled as small rectangular apertures with constant field on them. The far field is obtained with Love's principle of equivalence for aperture antenna [26]. In addition, a formulation is given to have the possibility of computing the near field in an arbitrary local coordinate system which is rotated and translated with regard to the array local coordinate system, which might be useful depending on the application. Both near field models are then compared by considering the near field radiated by a transmitarray antenna, showing a high degree of agreement between the two models. Then, the faster model is applied to a novel technique for the wideband near field beam-shaping based on the generalized intersection approach (GIA) [21]. This technique performs a direct layout optimization employing a method of moments based on local periodicity (MoM-LP) directly in the optimization loop to account for the real electromagnetic response of the unit cell at several frequencies of interest, including mutual coupling between the array elements. In this way, a general near field multi-frequency optimization algorithm is obtained. As an example of application, it is applied to an indoor scenario in which it is desired to provide coverage in a private office with a 5G NR femtocell working in the n257 band, which comprises the frequency range 26.5 -29.5 GHz, corresponding to a 10.7% relative bandwidth. The optimization procedure with the GIA is divided into several stages to facilitate convergence towards a suitable solution. The final optimized array presents a maximum magnitude ripple in the enlarged coverage area of 1.5 dB over the whole n257 bandwidth, showing the suitability of the proposed methodology for the wideband near field beam-shaping of spatially fed array antennas.

Statement of the Problem
Let us consider a spatially fed planar array, namely a reflectarray or transmitarray, comprised of N elements with a regular periodicity of p x × p y . The goal is to calculate the near field radiated by the array in a volume discretized in M points in front of the antenna and to shape the near field such that it complies with certain specifications for a given application. To that end, Figure 1 shows a sketch of the near field model for the two considered arrays. From the point of view of the analysis of the array, both spatially fed arrays are analysed in the same way, considering their own local coordinate system (ACS, array coordinate system) defined by (x a ,ŷ a ,ẑ a ). The location (i.e., centre) of the i-th array element is given by r i = (x i , y i , z i = 0), which is a point in the surface of the planar array, and r = (x, y, z) is a point in space in front of the array, i.e., with z > 0. Both r i and r are defined in the ACS, and R = r − r i .ˆˆˆì The goal is to obtain the near field radiated by the arrays in a volume in front of the antenna. The analysis of the arrays is carried out in the array coordinate system (ACS) defined by (x a ,ŷ a ,ẑ a ), while the radiated field will be obtained in the near field coordinate system (NFCS) defined by (x nf ,ŷ nf ,ẑ nf ). The transformation between ACS and NFCS may be consulted in Appendix A.
On the other hand, the radiated field is referred to a rotated-and translated-coordinate system (NFCS, near field coordinate system) noted as (x nf ,ŷ nf ,ẑ nf ). The relation between the ACS and NFCS is given by matrix T, which is derived in Appendix A. Since the analysis of the array is carried out in the ACS, the calculated radiated field will have to be transformed so that it is referenced to the NFCS.
In the following section, two near field models will be developed. One is based on the exact radiated field equations, while the other is based on the superposition of the far field radiated by each array element. In both cases, the radiation element is modelled as a small rectangular aperture with constant field on it.

Field at the Aperture
Before presenting the two near field models, the field at the aperture needs to be calculated, since it is used as the source from which the radiated field is obtained. First, a feed whose phase center is placed at r f = (x f , y f , z f ) in the ACS generates a tangential incident field on the surface of the array, which can be written for the i-th element as: with i = 1, . . . , N. In order to alleviate the notation, the dependence on r i of the field components will be dropped from here on. From the incident field, the field at the aperture may be calculated as: where M i is a 2 × 2 matrix of complex numbers which depends on the considered array. For reflectarrays, M i is the matrix of reflection coefficients: and the field at the aperture E ap would correspond to the reflected field. In the case of transmitarrays, M i is the matrix of transmission coefficients: and the field at the aperture would correspond to the transmitted field. For both cases, M i is computed using a full-wave analysis tool assuming local periodicity. In addition, M i depends on the working frequency, unit cell topology, periodicity, substrate characteristics, and angles of incidence of the plane wave employed in the analysis [27]. Once the tangential electric field at the aperture has been obtained with (2), the tangential magnetic field may be calculated with: where ω = 2π f , f is the frequency, µ 0 is the free-space permeability, k 0 is the free-space wavenumber, and (θ i , ϕ i ) is the angle of incidence of the plane wave for the i-th array element (see Figure 1). Theẑ component of E ap , necessary to solve (5), is obtained by solving the plane wave equation k · E ap = 0. Then, the components of the reflected tangential magnetic field are: with: In this way, the tangential field at the aperture is fully characterized by: where E ap ( r i ) is obtained from (2) and H ap ( r i ) from (5) and may be readily employed as a source for the calculation of the radiated field. Sometimes, the electric ( J) and magnetic ( M) current densities are required. Since planar arrays are considered, the notation J s and M s is employed. They can be obtained from the tangential field at the aperture as follows: wheren =ẑ a is the unit vector normal to the surface of the planar array. Please notice that the surface current densities are also tangential vectors, and thus the longitudinal components are zero, i.e., J z = M z = 0.

Near Field from the Radiation Equations
The first near field model is based on the radiation equations obtained as superposition of the contributions of the vector potentials A and F [25]: where In (12) and subsequent equations, the dependence of the field components on the observation point r has been dropped to alleviate the notation. There are a total of 12 field components in (11), for which an integral of field sources must be evaluated. An example of one such component is: where R is the distance from a source point in the volume V to the observation point, and β is the phase constant that in this case is equal to the free-space wavenumber k 0 . The rest of the integrals may be consulted in [25] (pp. 285-286). Since the arrays are planar, (13) can be reduced to: where S is the surface of the array. The array is divided into N unit cells, where the elements of matrix M i , and thus the field at the aperture and surface currents, are computed. Thus, we express the integral as a sum of the integral over each array element: where S i is the surface of the i-th unit cell, and: To finally obtain E F z in (15), the double integral may be evaluated by any quadrature method, for instance, the mid point quadrature. With this technique, the integral is approximated by a double sum: where surface S i is divided into N x N y rectangles of length ∆ x = p x /N x and width ∆ y = p y /N y , and: with: Thus, the function is evaluated at the center of each rectangle. By substituting (17) into (15), we finally obtain: The same process can be applied to the other 11 integrals in [25] (pp. 285-286), taking into account that the planar arrays are placed at the plane z = 0 and that the longitudinal components of the surface currents are zero (J z = M z = 0). For the sake of completeness, the final expressions for all the integrals particularized for spatially fed planar arrays are gathered next: where η is the free-space impedance, and:

Near Field as Superposition of Far Field Contributions
To develop this model, we have to consider each unit cell in the array as a unit of radiation and calculate its far field. In such a case, the total near field will be the superposition of the far fields radiated by each element: where E FF,i and H FF,i are calculated as the far field radiated by an aperture [26]. The electric far field is: and its components according to Love's equivalent principle: where the angles θ i and ϕ i are defined in Figure 2. P and Q are known as the spectrum functions, and they are the Fourier transform of the field at the aperture: where u = sin θ i cos ϕ i and v = sin θ i sin ϕ i . Since the aperture is element i of the array, we consider a constant field, which comes out of the integral, yielding: The double integral in (27) can be solved analytically yielding [28]: where sinc x = sin x /x is the unnormalized sinc function. Once the electric far field has been obtained, the magnetic far field may be obtained using the plane wave relation: The far field of an array element obtained with (24) and (29) is expressed in spherical coordinates, while we are interested in obtaining the near field in a Cartesian basis. Thus, to perform the superposition of (23), the far field of an array element must be expressed as: where the Cartesian components are found with: The same procedure is applied to the H FF field.

Change of Coordinates
The near field obtained with either (11) or (23) is referenced to the ACS. Moreover, to compute the near field following the equations derived in the previous subsections, the point in space r = (x, y, z) must also be expressed in the ACS. However, the goal is to obtain the near field in the NFCS, since the array may be rotated with regard to the desired local coordinate system (for instance, if the reflectarray is placed on a wall, and it is desired to provide indoor coverage on a table).
In general, the coordinates of the points where the near field is to be computed will be initially given in the NFCS. Thus, in order to compute the near field with (11) and (23), these coordinates need to be transformed from the NFCS to the ACS as detailed in Appendix A. Then, after the near field is computed in the ACS, it needs to be transformed to the NFCS as follows: where T is defined in Appendix A, and E NF ( r) ACS and H NF ( r) ACS are the same as in (11) and (23).

Comparison of the Near Field Models
To compare the near field predicted by both models, let us consider a center-fed square planar transmitarray comprised of 60 × 60 elements in a regular grid of periodicity 3.84 mm 2 . The feed is modelled as a cos q θ function [29] with q = 22 and placed at r f = (0, 0, −180) mm in the ACS (see Figure 1). The working frequency is 39 GHz. For the analysis, the phases of the direct transmission coefficients in (4) are set to: where d i is the distance from the feed phase center to the i-th array element, and (θ 0 , ϕ 0 ) is the pointing direction of the far field pencil beam defined by the phases of τ xx and τ yy . For the present example, (θ 0 , ϕ 0 ) is set to (0°, 0°). The radiated near field will be calculated with the two models in planes parallel to the array with z constant and values z = 0.1 m, 1 m. Thus, the matrix of change of coordinates from the ACS to the NFCS is the identity. In addition, in the case of the model based on the radiation equations, it will be evaluated for several values of N x and N y , specifically for N x = N y = 1, 3, 5 to assess the convergence of the integrals. Finally, each plane is divided into a grid of 150 × 150 points (22,500 in total) where the near field is computed. Figure 3 shows the main cuts in x for y = 0 for the near field component E x in magnitude and phase for both z planes. The agreement between both models is in general very good. In the case of the magnitude, the model based on the radiation equations with N x = N y = 1 presents some discrepancies in the plane z = 0.1 m with the rest of the models, albeit for very low values, 60 dB lower than the peak magnitude. In the case of the phase, the model based on the superposition of the far field radiated by each array element presents a small discrepancy in the plane z = 0.1 m. In this case, the shape of the phase is the same as in the other models, but it has a shift of approximately 1.3 deg. Both discrepancies disappear when the field is calculated further away from the array.   Figure 4 shows the magnitude and phase of the field components E x and E y as calculated by the models based on far field superposition (top row) and the radiation equations with N x = N y = 5 (bottom row) for the plane z = 0.1 m. As can be seen, visually the results are practically the same.
A better way of assessing the similarities between both models is to compute the relative error in the calculated near fields. To that end, the following expression is employed: where · is the Euclidean norm; E is one component of the near field computed with the model based on the radiation equations with N x = N y = 5, which is considered as reference to compute the relative error in percentage; andẼ is the same component of the near field computed with any other model. Table 1 shows the relative error for both models. It shows the relative errors for the three components of the near field for both planes. In the case of the model based on the radiation equations, even a value of N x = N y = 1 shows a sufficiently low error. As the near field is computed further away from the array, the error decreases. The use of N x = N y = 1 means that the double integrals in (15) are just evaluated at the center of the array element, i.e., considering only one sample, effectively reducing the triple summation in (21) to a single sum, considerably accelerating computations without barely affecting accuracy in the prediction of the near field. In the case of the near field model based on the superposition of the far fields radiated by each array element, the relative error is slightly larger, and it also decreases as the near field is computed further away from the array. However, even for a plane at z = 0.1 m, which is 13λ at the working frequency, the error is very low.   Table 1. Relative error of the field components calculated with the different near field models when compared with the model based on the radiation equations with N x = N y = 5, which is taken as reference.
In light of these results, it is clear that both models offer practically the same results in the computation of the near field for spatially fed planar arrays. However, from a computational point of view, the model based on the superposition of far fields is more interesting. Indeed, as shown in Table 2, this model is at least one order of magnitude faster than the model based on the radiation equations. This is particularly important when employing the model in optimization procedures to shape the near field, since the analysis routine where the near field is computed is invoked many times. Then, the model based on the radiation equations could be used, if desired, to assess the validity of the results that were obtained employing the approximate model based on superposition of the far fields radiated by the array elements. Table 2. Computational performance of each algorithm in the computation of the near field in a total of 22,500 points employing an Intel i9-9900 CPU working at 3.1 GHz. Computations are parallelized with OpenMP [30].

Scenario Definition and Antenna Specifications
As an example of application of the developed models for near field beam-shaping, let us consider the scenario depicted in Figure 5, where the goal is to provide coverage on top of an office desktop. The coverage area is a rectangle of dimensions 0.9 × 0.5 m 2 on the right side of the table. The reflectarray that provides the coverage is placed on the ceiling, at a distance of 2 m from the center of the coverage area, while the feed is placed on a shelf to illuminate the reflectarray. The goal is to achieve a maximum magnitude ripple in the coverage area of 2 dB in the whole 5G NR n257 band (26.5 -29.5 GHz) in dual-linear polarization.ˆ−ˆ=

Unit Cell Characterization
In order to be able to successfully perform a wideband direct layout optimization, the chosen unit cell should be able to provide a large enough phase-shift range as well as enough degrees of freedom for the optimization. In this regard, a unit cell comprised of only one rectangular patch would not provide enough phase-shift and degrees of freedom for its use in a wideband optimization procedure [32]. However, by stacking two rectangular patches, the range of the phase-shift can be increased. Thus, the chosen unit cell consists of two stacked rectangular patches backed by a ground plane as shown in Figure 6a. Both layers employ the same substrate, the Rogers DiClad 870 with relative permittivity ε r = 2.33, loss tangent tan δ = 0.005, and thickness h = 30 mil. A CuClad 6250 bonding film between the two layers with the same electrical characteristics as the DiClad 870 is considered in the simulations. The substrate has been chosen such that a large enough phase-shift range is achieved. For the initial layout design after the POS (see Section 4.3), a fixed relative size of the stacked patches is kept, with T x 1 = 0.77T x 2 and T y 1 = 0.77T y 2 . The scaling factors were found after performing a parametric study to obtain a smooth variation of the phase-shift for T x 1 = T y 1 at several frequencies and angles of incidence. Figure 6b shows the phase response of the unit cell at 28 GHz for normal and oblique incidence. The total phase-shift provided by the cell with the selected substrate is around 500°, and it can be seen that it is stable with the angle of incidence. Conversely, Figure 6c shows the unit cell phase response at three different frequencies for oblique incidence with (θ = 40°, ϕ = 30°). From these curves, it can be seen that the unit cell also provides around 500°of phase-shift at several frequencies, making it suitable for a multi-frequency optimization. Figure 7 shows the proposed wideband near field optimization strategy for spatially fed array antennas. It is divided into several stages in order to facilitate convergence towards a suitable solution. After defining the antenna and near field specifications, a phase-only synthesis is carried out at central frequency. Afterwards, the geometrical features of the unit cell are tuned so each element produces the desired phase-shift obtained from the POS. After this point, the near field radiated by the antenna closely complies with the specifications at central frequency. However, due to the narrow bandwidth nature of planar spatially fed arrays, it does not comply at other frequencies. Thus, this initial layout is used in a multi-stage direct layout optimization (DLO) at several frequencies, finally obtaining an optimized wideband layout. To obtain an array working in dual-linear polarization, the POS is carried out independently for both polarizations, while the DLO is performed considering specifications in both polarizations simultaneously.  Figure 7. Proposed wideband optimization strategy for spatially fed arrays using the generalized intersection approach (GIA). First, a phase-only synthesis (POS) at central frequency is carried out in order to obtain a phase-shift distribution from which the initial layout is obtained. Then, the layout is used in a direct layout optimization (DLO) at multiple frequencies until specifications are met.

Phase-Only Synthesis at 28 GHz
For the POS, the generalized intersection approach (GIA) [33] will be employed. The POS assumes an ideal phase-shifter simplification in matrix M i , in which there are no losses (|ρ xx,i | = |ρ yy,i | = 1 for reflectarrays, or |τ xx,i | = |τ yy,i | = 1 for transmitarrays) and no cross-polarization (ρ xy,i = ρ yx,i = 0 for reflectarrays, or τ xy,i = τ yx,i = 0 for transmitarrays). In this way, only the phases of the direct coefficients are considered as optimization variables. Once the POS has converged to a near field that complies (or it is close to compliance) with the imposed specifications at central frequency, the sizes of the rectangular patches of the unit cell are adjusted such that the unit cell provides the required phase-shift that was obtained in the POS. This procedure can be summarized as follows. First, and for each reflectarray element, two phase shift tables are generated independently for each polarization. Then, the dimensions of the rectangular patches are adjusted using a linear equation. Finally, a zero-finding routine is employed to adjust both dimensions at the same time. During this process, a fixed relative size between the stacked patches is kept, T x 1 = 0.77T x 2 and T y 1 = 0.77T y 2 . Further details on the GIA applied to POS may be found in [21,33].
Since the employed algorithm is a local optimizer, the starting point of the optimization is very important. Thus, the starting phase-shift distribution will be such that it provides a focused near field at the center of the coverage zone. To that end, the conjugate-phase approach [34] is employed, which establishes that the phase of the field at the aperture should be: for i = 1, . . . , N. On the other hand, the phase of the field at the aperture is also equal to the phase introduced by the unit cell plus the phase of the incident field from the feed [32]: where d i is the distance from the feed phase center to the i-th unit cell, and ρ is a direct reflection coefficient (either ρ xx or ρ yy , depending on the polarization). Thus, by combining (35) and (36), we can obtain the phase-shift introduced by each reflectarray element to focus the near field at a point r = (x, y, z): Figure 8a shows the focused near field for linear polarization X when the phases in (37) are implemented. The magnitude of the field is concentrated in the middle of the coverage zone, but it presents an important drop towards the edges of the coverage zone. Indeed, the ripple in this case is 48 dB for polarization X and 38 dB for polarization Y. After the POS is completed, a phase-shift distribution is obtained such that the ripple is considerably reduced at 28 GHz. A layout is obtained from this phase-shift distribution by following the steps described above and in [23] and using the MoM-LP described in [35]. When this layout is simulated, the near field shown in Figure 8b is obtained. As can be seen, the near field presents a more uniform distribution in the coverage area when compared to that of the starting point of the POS shown in Figure 8a. Now, the ripple for polarizations X and Y is 2.66 dB and 2.48 dB, respectively. Although it does not comply with the specifications of a maximum ripple of 2 dB, it is a net improvement over the starting point. However, we are also interested in the performance in the 5G NR n257 band (26.5 -29.5 GHz). To that end, Table 3 gathers the results at seven equispaced frequencies, showing the maximum achieved ripple as well as the percentage of the coverage area that complies with a 1 dB, 2 dB, or 3 dB ripple. As can be seen, the best results are obtained at the frequency of design, 28 GHz. For frequencies 500 MHz above and below 28 GHz, the results deteriorate slightly. However, as the frequency moves further away from the design one, the results quickly worsen, showing the narrow bandwidth behaviour of the reflectarray antenna and the need for a wideband optimization.

Multi-Frequency Optimization of the Near Field
The multi-frequency optimization of the near field consists in a direct layout optimization using the GIA [21]. To accommodate specifications at a number of frequencies, the cost function (F) in the backward projector of the GIA is modified in the following fashion [36]: where Ω is the region where the near field is computed; N f is the total number of frequencies at which the direct layout optimization is performed; L is the number of points in which the region Ω is discretized; w f is a weighting function that depends on the frequency and observation point r l ; E NF, f ( r l ) is the computed near field at each iteration of the GIA and frequency;˜ E NF, f ( r l ) is the reference near field in region Ω resulting from applying the specification templates at a given frequency; and ∆Ω is the step of the discretization. Further details on the specifics of the optimization algorithm may be consulted elsewhere [21]. Unlike in POS, where the optimizing variables (sometimes also referred to as degrees of freedom or DoF) were the phases of the direct coefficients, the direct layout optimization employs the geometrical features of the unit cell as optimizing variables. In the present case, the dimensions of the rectangular patches are employed, T x 1 , T x 2 , T y 1 , and T y 2 (see Figure 6). Thus, a total of 7744 DoF are available for the optimization. Moreover, the direct layout optimization will be carried out in several stages in order to facilitate convergence towards a suitable solution. First, only variables T x 1 , and T x 2 will be considered at three frequencies, 27 GHz, 28 GHz and 29 GHz. In the second stage, the near field is optimized at the same three frequencies, but this time only considering as variables T y 1 and T y 2 . Finally, in the third stage, all four variables per element are used in the optimization at the same seven frequencies of Table 3. In all stages, the specifications for both linear polarizations are imposed. In addition, the MoM-LP tool is directly employed in the optimization loop to obtain the real electromagnetic response of the unit cell at each considered frequency. In order to accelerate computations, the differential contributions technique [37] is used.
The improvement in the coverage zone can be visually seen in Figure 9, where the normalized magnitude of the electric field for linear polarization X is represented at three different frequencies (26.5 GHz,28 GHz,and 29.5 GHz) before and after the direct layout optimization. The field is only plotted at those points that comply with the requirement of 2 dB of ripple. Although the starting layout does not radiate a near field that fully complies with specifications at 28 GHz, Figure 9b shows that most of the coverage zone does comply, and only at one edge is the ripple higher than 2 dB, demonstrating the successful optimization carried out by the POS. However, due to the intrinsic narrow bandwidth behaviour of the reflectarray, at other frequencies, it does not comply. Nevertheless, after the direct layout optimization, in which the sizes of the two stacked rectangular patches are directly optimized, the antenna is able to provide full coverage in the designated area with a maximum ripple lower than 2 dB.  Table 4 gathers the performance of the optimized layout at seven equispaced frequencies covering the whole 5G NR n257 band. Now, a compliance of 100% is achieved in the whole band for a ripple of 2 dB. Indeed, the worst achieved ripple is 1.5 dB at 29.5 GHz for polarization X. When compared with the results shown in Table 3, the worst ripple has been reduced more than 7.5 dB, and even the compliance for a ripple of 1 dB, although not considered a requirement, has improved from the lowest value of 26.4% for polarization X at 29.5 GHz to 76.7% for polarization X at 28.5 GHz. In the case of the 2 dB compliance columns, the minimum value before the direct layout optimization was 44.3%, which was improved to a 100% compliance after a successful wideband optimization.

Discussion
The proposed wideband near field beam-shaping technique was able to achieve a maximum ripple of 1.5 dB in the whole 5G NR n257 band, corresponding to a 10.7% relative bandwidth. It is worth mentioning that specifications have been met employing a unit cell consisting of two stacked rectangular microstrip patches backed by a ground plane, offering up to four degrees of freedom to perform a direct layout optimization. However, bandwidth could be further improved by employing a more suitable unit cell. For instance, the three stacked patches configuration is known to offer higher bandwidth of operation than only two layers [38]. Alternatively, a unit cell comprised of two layers of parallel and coplanar dipoles is known to be able to offer, after a direct layout optimization, relative bandwidths of up to 20% [36]. When applied to near field beam-shaping, this would allow offering coverage in both the n257 and n258 bands, covering from 24.25 GHz up to 29.5 GHz, which represents 19.5% of relative bandwidth. Additionally, even though the use of a MoM-LP tool directly in the optimization loop penalizes computational performance to achieve wideband results, some techniques could be employed to accelerate the process. In this work, the differential contributions technique [37] was used to accelerate some building blocks of the GIA. However, the unit cell analysis could also be sped up by substituting the MoM-LP with a database of reflection coefficients [39] or the use of machine learning techniques [40].

Conclusions
This work has presented two near field models for the analysis of spatially fed array antennas, namely, reflectarrays and transmitarrays. The first model is based on the direct application of the radiation equations directly derived from the A and F vector potentials. The volumetric integrals are particularized for a planar array discretized into a number of rectangular unit cells, and thus they are reduced to a summation of a double integral over each unit cell. In turn, this integral is numerically evaluated by a quadrature method, for instance, the mid point quadrature. The second model is based on the superposition of the far field radiated by each array element. The far field is computed using Love's principle of equivalence assuming a constant field at each unit cell. Both near field models were compared by computing the field radiated by a transmitarray antenna. For the model based on the radiation equations, the near field was also computed for several discretizations of the unit cell for the numerical computation of the double integrals. Both models provide very similar results. When the model based on the radiation equations is used as reference, the relative error in the radiated field obtained with the other model quickly drops as the field is computed further away from the array. With relative errors lower than 3% in the computed tangential near field at 13 λ or further from the array, the model is effectively validated for subsequent use in optimization algorithms for near field beam-shaping. Then, the model based on far field superposition was applied to perform near field beam-shaping. The considered scenario is aimed at providing indoor coverage in an office for mm-Wave 5G NR networks in the n257 band. The proposed wideband near field opti-mization strategy is composed of several stages to facilitate convergence towards a suitable solution. First, a phase-only synthesis (POS) was carried out to obtain a reflectarray layout that provides coverage at central frequency (28 GHz). Due to the narrow bandwidth nature of spatially fed planar array antennas, poor compliance is achieved at other frequencies within the n257 band. Thus, a multi-frequency optimization was then carried out employing a method of moments based on local periodicity (MoM-LP) directly in the optimization loop along with the near field model. The cost function of the optimization algorithm was modified to include requirements at several frequencies in the near field region of interest. At the same time, the optimization was carried out in several stages, gradually increasing the number of DoF to improve convergence. Eventually, the optimized antenna complied with the imposed specifications, achieving a magnitude ripple equal or lower than 1.5 dB in the whole 5G NR n257 band, corresponding to a 10.7% relative bandwidth, showing the suitability of the proposed methodology for the wideband near field beam-shaping of spatially fed arrays.

Conflicts of Interest:
The author declares no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: igure A1. Rotation of the array coordinate system (ACS) (x a ,ŷ a ,ẑ a ) to obtain the near field coordinate system (NFCS) (x nf ,ŷ nf ,ẑ nf ). This rotation is specified by angles θ, ϕ, and ψ.
Then, if r a = (x, y, z) ACS is a point in the ACS, and r nf = (x, y, z) NFCS is the same point in the NFCS, the matrix of change of coordinates from the ACS to the NFCS is: with: nf ·x axnf ·ŷ axnf ·ẑ â y nf ·x aŷnf ·ŷ aŷnf ·ẑ â z nf ·x aẑnf ·ŷ aẑnf ·ẑ a   . (A5) Conversely, to perform a change of coordinates from the NFCS to the ACS we can calculate: where T T is the transpose of T, which is an orthogonal matrix, i.e., T T = T −1 .