Abstract
We present an efficient method for the analysis of different objects that may contain a complex feeding system and a reflector structure. The approach is based on a domain decomposition technique that divides the geometry into several parts to minimize the vast computational resources required when applying a full wave method. This technique is also parallelized by using the Message Passing Interface to minimize the memory and time requirements of the simulation. A reflectarray analysis serves as an example of the proposed approach.
1. Introduction
In the last years we have seen a rising interest in computer simulation tools in several engineering and science fields, due to technological breakthroughs and the development of novel efficient numerical algorithms. Many researchers are involved in the development of new methods focused on the acceleration and improvement of conventional electromagnetic field solvers. The logical consequence of this fact is a great increase in the demand for computer simulation tools.
Addressing Method-of-Moments- (MoM-) based electromagnetic solvers, we have seen a good number of efficient approaches in the last years. Many of these are related to the reduction of the computational requirements of the legacy MoM by only storing the coupling between basis and testing functions that are located in the near field region. There are very well-known methods, like the Multilevel Fast Multipole Algorithm (MLFMA) [1], that compute the matrix-vector products efficiently by aggregating a number of active currents to some specific points, translating them and disaggregating those multipoles in order to take into account the effect of the active basis functions over the testing functions. Other efficient techniques reduce the total number of unknowns by defining extended sets of macrobasis functions [2–6]. Fast matrix compression and low-rank approximation approaches can also be mentioned as efficient methods in terms of memory and CPU time requirements [7, 8]. In addition to the MLFMA there are a number of techniques that also take advantage of the fast evaluation of matrix-vector products in the iterative solution of large problems, among which we can mention the Complex Multipole Beam Approach (CMBA) [9], the Impedance Matrix Localization (IML) technique [10], the Adaptive Integral Method (AIM) [11], or the Multilevel Matrix Decomposition Algorithm (MLMDA) [12]. However, there are some cases for which the use of a full-wave method still poses a bottleneck, even considering these fast approaches, due to the size of the problem, slow convergence, or other issues. In these cases it is sometimes possible to resort to high-frequency methods [13] or hybrid combinations [14] if the nature of the geometrical properties of the problem and the excitation lie within the scope of such techniques. Unfortunately this is not usually the case when analyzing reflectarray antennas. Domain decomposition techniques (DD) allow to limit the maximum size of the system matrix to be solved [15] regardless of the size of the objects involved. It is necessary to include interactions between the domains in order to reach the final solution, which is usually done in an iterative or recursive way [16–18]. There are some works worthwhile to be mentioned that study the convergence problems that may prevent reaching a stop situation [19, 20]. Some typical strategies to avoid this undesirable divergence rely on the definition of a relaxation factor [19, 21], make use of an extension of the original domains [17, 18, 21] and/or use preconditioners [17].
Due to the interest raised by the analysis of reflectarray antennas and similar geometries, a number of works reported in the last years deal specifically with this kind of problems. For instance, [22] presents an efficient technique based on the extension of an adaptive integral method that allows the full-wave analysis of electrically large multilayered printed arrays that have one or more planar metallizations and vertical conductors. The technique proposed in [23] computes the generalized scattering matrix of a dielectric interface with periodic metallizations, which is based on a spectral domain Moment Method, but assuming multiple incident Floquet-harmonics. A full-wave analysis using the finite integration technique is applied in [24], where the mutual coupling between the feeding horn and the elements of the reflectarray is considered. The analysis reported in [25] is based on an application of the Moment Method. In [26] the authors make use of a fast simulation tool based on MoM-MLFMA to analyze a three-layered reflectarray. The work reported in [27] combines MoM-MLFMA with the use of macrobasis functions to address large problems, including reflectarray antennas.
In this work we propose a combination of the MLFMA with a new domain decomposition technique that provides fast analysis of arbitrary structures in complex environments and contains a novel procedure for selecting the interactions between domains that will have a noticeable contribution to the final results. This new domain decomposition method has been validated and compared with the MoM-MLFMA approach. The simulation results of both approaches are compared with measurements of a reflectarray antenna that has been previously designed, fabricated, and measured [28] in the Research Advanced Antenna Technology Communications Research Centre, Canada.
The geometry of the bodies considered in this work is defined using parametric surfaces [29, 30], as well as the mesh required in order to apply the MoM-based approach. We obtain a continuous quadrangular mesh by using the mesher described in [31]. The use of parametric subpatches enables us to conform them to curved surfaces and maintain a high accuracy degree reducing at the same time the sampling rate (which in traditional approaches is considered to be around 10 samples per wavelength). It represents an advantage when compared with other approaches that use a faceted meshing, which loses the real shape of the original bodies. Once the geometry has been discretized, high-order modified rooftop and razor blade functions are used as basis and testing functions, respectively. In addition, the proposed approach also allows the analysis of structures embedded in dielectric slabs. In this case, volumetric rooftops and razor blade functions are defined over curved volumetric elements. For very thin dielectric slabs a surface representation is considered instead of the previously mentioned volumetric approach, and the MoM impedance matrix elements are corrected by a term that depends on the size of the subdomain, the frequency, the electric permittivity, and the thickness of the slab.
Regarding the convergence issues that can be derived from the application of domain decomposition methods, the approach proposed in this work relies on the definition of extensions for each domain. After computing the updated currents over the extended version we only retain those which are included in the original, nonextended domain. We reach a double goal with this strategy: the improvement of the convergence properties in the iterative process, and the correction of the artificial edge behavior of the current. Extension values of about 0.4 to 0.7 wavelengths typically yield good results. The negative side of this strategy is the extra CPU time required to solve the extended problem. However we have found that it only represents a small fraction of the total simulation time for the typical domain sizes considered.
The MPI is used for the parallelization of the approach. This paradigm provides communication functionalities between a set of processes. Therefore, it is possible to address this kind of problems using a multiprocessor computer or cluster, increasing the computational efficiency [32, 33].
The paper is organized as follows. The main features of the proposed approach are addressed in Section 2. Some remarks about the selection of active and passive domains are given in Section 3. Section 4 presents the results obtained when analyzing a specific reflectarray antenna, and Section 5 contains the conclusions derived from this work.
2. Background of the Technique
The domain decomposition approach combines the MoM and MLFMA together with an iterative process used to calculate the interactions between different domains, which will be described in the next section. The first step will be, therefore, a partitioning process of the geometry into several domains. Thus, the total current over the original geometry can be expressed as where denotes the number of domains. Subsequently, a full wave analysis is individually applied over each domain using the MoM and MLFMA, obtaining the primary currents, which are those obtained by considering each block isolated from the rest and the excitation vector given by the external field. Several small problems are solved instead of a very large one, which is the main advantage when compared to the analysis of the whole structure.
The artificial edges derived from the partitioning of the objects in terms of domains may affect noticeably the final results. The solution to this potential burden is to extend the area of each domain and compute the induced currents while considering the extended domain. Figure 1 illustrates this decomposition procedure and one example of an extended domain.
Once the currents have been calculated, those that are located inside the extension are discarded. By doing this we avoid the negative effect of the singular behavior of the currents near the artificially introduced edges. According to our experience an extension with a side length between 0.7 to 1.0 wavelengths is enough to obtain good accuracy.
After computing the primary currents an iterative process is performed to calculate the influence of the currents obtained inside one domain over the others. This stage of the process is required in order to update the primary currents calculated with the MoM. In order to compute the interactions between different parts of the geometry, a selective process is necessary to avoid unnecessary calculations. Let us denote as active domains the ones that contain previously computed currents which will radiate over other domains of the geometry (these ones called passive domains). It can be inferred that a given domain D can serve as active and passive domain at the same time, since the current distribution calculated over D in the kth iteration may produce interactions with other parts of the geometry in the (k + 1)th iteration, while other active domains may contribute to update the current distribution on D as well. A ray-tracing approach determines which active domains contribute to the updating current process of a certain passive domain, and therefore only a few computations are required to obtain the final currents. Multiscale problems can take advantage of this methodology because only the resulting fields radiated by an active domain are required to feed the passive domain.
3. The Iterative Selection of Active and Passive Domains
In order to increase the efficiency of the approach, it is very important to discard interactions between active and passive domains when these interactions are negligible. This strategy reduces significantly the CPU time required. Moreover, realistic structures with complex small details can be efficiently analyzed because no geometrical approximation is considered. The immediate benefit of this method is a significant decrease in the required CPU time and memory storage resources.
The total current over a given domain (let us consider domain-i) can be seen as the superposition of the currents computed by a number of interactions between that domain and the rest of the domains contained in the geometry: where the term represents the primary current computed over domain-i. The other terms of the expression indicate the currents due to the interaction between domain-i and the rest of the domains, where stands for the highest order of the interactions with domain-i considered as the passive domain.
In order to compute the interactions between different domains, it is necessary to calculate the fields radiated by the previously computed currents on the rest of the domains over each passive domain, which can be iteratively performed. Thus, the impress field over domain-i due to the kth order effect can be written as where the operator computes the field over domain-i due to the previously calculated currents over the rest of the geometry, which serve as active domains. The MLFMA approach is very well suited to obtain . Now it is possible to obtain the current over domain-i due to the kth order effect, , by solving the conventional MoM equation for the isolated domain. This current can be efficiently obtained by utilizing MLFMA with an iterative solver such as the Biconjugate Stabilized Gradient method (BiCGStab).
Expression (3) needs to be modified in order to select those pairs of active and passive domains that will entail a noticeable contribution to the final current, and discard the negligible interactions. In practical applications this procedure allows us to reduce the amount of primary memory and CPU time maintaining the same degree of accuracy. In order to take into account this consideration, the modified expression can be seen as follows: where we have introduced the function that takes a zero value if the interaction between domain-i and domain-n is negligible and takes the value of if such interaction entails a noticeable contribution to the current distribution over the passive domain. For the selection of the passive domains which interact with an active one it is important to consider the possibility of total or partial eclipse situations, as well as the shape of the pattern radiated by the active currents. We utilize the following mechanisms to accomplish this task.
First, using the multipole aggregation of the current over the active domain we determine the main directions in which such domain radiates, discarding the passive domain if it is not located in the path of the main set of output directions. This set of directions is determined by setting a threshold over the strongest field value.
In addition to the previous procedure we also consider a ray-tracing algorithm, very successfully used in high-frequency methods, known as the Angular Z-Buffer (AZB) [34, 35]. With the AZB the space is divided into angular regions (anxels), and the surfaces of the geometry are classified depending on which anxels contain them. This information is stored in the AZB matrices in a preprocessing step. As an example, let us consider a domain of the problem called . The AZB matrix associated to contains a number of cells which are assigned to different angular ranges, since one of the dimensions of such matrix makes reference to the theta position of the angular range and the other dimension indicates its phi position. For the generation of the AZB matrix, different ray tubes are shot over (depending on the nature of the external excitation) and the output directions are obtained. The output ray tubes are then checked for intersections with other domains. If some rays belonging to the same tube intersect with different domains, these are processed in order to analyze the possibility of eclipse effects. Only the nonhidden domains are included in the AZB matrix of . Further information can be found in great detail in [36]. This technique is very well suited to be used in combination with the previous approach based on the multipole expansion of the currents over the active domain, because the former gives the initial tube rays to be considered for the search of passive domain candidates to interact with the active domain.
Thus, for a given direction (from the main output radiation directions of the active domain) the AZB yields the passive domains where there can be a strong coupling, discarding the negligible interactions. This can be extended to consider multiple effects between domains, as shown in Figure 2 where , , and are domains contained in the geometry.
4. Analysis of a Reflectarray
In this work we have modified and analyzed a reflectarray antenna for applications such as radar or long-distance communication. Some preliminary results of the numerical approach presented here can be seen in [37]. The design of the broadband reflectarray is presented in [28] and provides high radiation efficiency, approximately 60%. Regarding its construction, flat reflectors are often defined by several layers of periodical metallic patches and dielectric slabs over a ground plane [38]. Although microstrip reflectarray antennas are narrow in bandwidth, they are well suited for dual-frequency operation [39]. In addition to the original reflectarray considered here, there are other works in the literature presenting dual-band reflectarrays [40–43]. Currently, reflectarray antennas are widely used as an alternative to conventional reflector antennas due to several advantages such as reduced weight and volume, mechanical robustness, compatibility with active devices, low cost, and short manufacturing time.
Figure 3 shows the original reflectarray and the shape of the unit cells. The individual elements are designed to scatter the incident field while impressing the appropriate phase shifts in such a way as to form a plane wave front that propagates in a prescribed direction or, more generally, to create the required radiation pattern. The dimensions of the whole square structure are 40 cm × 40 cm, and the reflective surface is located over the XY plane.
The modifications we have introduced consist of a new compact feed system comprising an offset cassegrain structure where the horn is located under the plane of the reflectarray, for which it is necessary to open a hole on the reflectarray to feed the subreflector. This configuration can be seen in Figure 4. The feeding antenna is a pyramidal horn polarized along the X-axis. Figure 4 shows the location of the feed antenna over the reflectarray, and Figure 5 displays its far-field radiation pattern.
(a)
(b)
The aim of the modification provided in this work is the generation of a more compact layout for this antenna, while maintaining its original features up to a certain degree. The Cassegrain structure allows this goal, with good radiation properties as the final result.
The reflectarray is composed of a dielectric layer over a metallic ground plane, and the metallic strips are located over the dielectric. The period of the reflective elements is defined as 12 × 12 mm. The substrate has a thickness of 3.175 mm and a permittivity of 2.17. The simulation has been performed considering a superficial thin layer approximation for the dielectric region.
The reflectarray was originally designed to radiate in the direction defined by θ = 16° and φ = 0° at 14 GHz. The horn is located at x = y = z = 0 and feeds a subreflector whose center is x = y = 0, z = 0.21 and radius is 0.106 m. The focal length of the subreflector is 0.215 m, and the offset is −0.02 m. We have computed the radiation patter of the geometry using two different approaches, namely, the full wave technique given by the MoM-MLFMA and the domain decomposition approach described in this work. With the latter method we have considered 6 domains in the geometry: the first one comprises the horn, the second domain contains the subreflector, and the reflectarray is partitioned into 4 domains. The CPU-time required by the MoM is 4806 seconds, while the Domain decomposition approach took 987 seconds to obtain the results. The total number of unknowns when analyzing the reflectarray antenna at 14 GHz is 271636 unknowns. Around 27% of that amount (73720 unknowns) belongs to the dielectric part. However, the efficiency of the developed technique compared to the conventional full-wave analysis increases when the electrical size of the structure is even higher.
Figure 5 depicts the comparison between the radiation patterns obtained by these two approaches at 14 GHz. The direction of the main beam is now given by θ = 18° and φ = 0° at both frequencies. A good agreement can be appreciated, with some discrepancies for the lower levels, due to the fact that, even considering the offset, there is still some blockage owing to the location of the horn and subreflector. A 3D view of the pattern obtained by the proposed approach can be seen in Figure 6, and the current distribution obtained is shown in Figure 7. In Figures 8 and 9 we compare the radiation pattern obtained with the modified configuration with the original configuration. It is possible to appreciate the differences due to the effect of the subreflector blockage and the hole required for the new horn placement, but the overall behaviour of this more compact configuration can be of interest for a number of practical applications.
(a)
(b)
(c)
5. Conclusion
An efficient approach based on the MoM-MLFMA combined with a domain decomposition scheme has been presented and validated comparing to the full-wave solution. This approach considerably reduces the memory and time requirements by solving several small problems instead of a single large one. The predictions given by both methods show a good agreement. The MPI paradigm was incorporated to increase the performance in the simulations.
Acknowledgments
This work has been supported in part by the Comunidad de Madrid Project S-2009/TIC1485, by the Castilla-La Mancha Project PPII10-0192-0083, by the Spanish Department of Science, Technology Projects TEC 2010-15706 and CONSOLIDER-INGENIO no. CSD-2008-0068, and by University of Alcalá Project UAH2011/EXP-014. The authors thank Dr. Reza Chaharmir for the original model and measurements of the reflectarray antenna reported in this paper.