Optimization Approaches in Sideband Calculations and a Non- Iterative Harmonic Suppression Strategy in 4D Arrays

3.1. A general overview on metaheuristic optimization

Metaheuristics may be defined as the nature-inspired algorithms which mimics a natural phenomenon in order to achieve a certain goal. To be specific, the main purpose of this kind of approach is to find a close-optimal parameter set over a specifically constructed mathematical space so-called a search space with some limiting constraints for the problem at hand. Most of the metaheuristic algorithms operate on a point-based search mechanism in order to keep the dimension space unchanged. In other words, these kinds of algorithms operate on a search space sampling strategy which forms a parallel optimal point search mechanism over the search space. These search space samples form a potential solution set called, usually, the population. If there is no pre information about search space behavior, usually, these optimizers start their search from a uniformly distributed search space sample set so-called the initial population.

Here, the search space mentioned above is a user-defined mathematical function space, usually called as cost function, fitness function, and so on, which varies depending on the problem construction. Definition of the search space is one of the key issues for these kinds of algorithms since it directly defines the problem. A poorly defined search space may lead inappropriate results even the algorithm works properly. A search space may be constructed in several ways but constructing it on an error function such as mean square or absolute error is a popular approach since the optimal point in an error function is predictable or known. However, using an error-based function is not a necessity, and some other well-defined functions may be used as the function that is to be optimized. If an error-based cost function is formed for a problem at hand, typically the problem becomes a minimization process, and the algorithm tries to find the possible solutions close to a desired one with a minimal error. A typical error-based cost function may be defined as

where x¯stands for the parameter vector to be optimized, y¯is the constraint set and fx¯is the cost function may be defined by fx¯=∑igidx¯−giax¯21/2. Here, in the cost function, g, d and a represent the specification that is to be optimized, the desired level and the actual level, respectively.

After constructing and sampling the search space, that is, defining the cost function and constructing the initial population, all point-based metaheuristics move these initial samples to different points by their predefined operators. Actually, at this point, it should be noted that mainly the algorithms differ from each other by their relocation processes. Each algorithm uses a different method to relocate the initial samples, and usually this relocation process gives the algorithms name. By the relocation of the initial samples, a newly generated sample points are created which stands for new solutions. Comparing these new solutions with old ones and keeping the suitable ones and discarding the others cause an overall improvement in the solution set. In other words, after creating new points by relocation process, comparing them with already known possible solutions and keeping the better ones provide an overall convergence to the optimal point. This comparing and selecting process may be conducted in several ways such as using a roulette wheel or a binomial selection. In any way, the basic idea behind this process is comparing the solutions and keeping the better ones in order to improve the solution quality at hand. The relocation and selection process continues until a satisfied convergence rate is achieved or an iteration limit is reached. By this way a close optimal point may be found after the overall process in an iterative manner.

The basic idea behind the metaheuristics may be defined in main four steps as

1. Step 1: Constructing a search space

2. Step 2: Forming a randomly generated initial sample set

3. Step 3: Moving the samples in search space and create new solutions

4. Step 4: Comparing the new set with old one and keeping the suitable ones

As stated above, Step 3 and Step 4 are repeated in an iterative manner to achieve a satisfactory result. Here it should be noted that, because of their randomized nature, none of the metaheuristics guarantees to find an optimal solution. But it can be said that elitist methods (the methods that always keeping the current best solution in the solution set) with enough randomization of the solution set (this process is usually necessary to avoid local optimal points) guarantees to find a close optimal solution with enough convergence rate as the time tends to infinity.

3.2. Application to the harmonic suppression problem

As all metaheuristic approaches, definition of the problem is the key issue in sideband suppression problems. In TMA design which the sidebands are not to be used on specific purpose, the sideband radiation levels are preferred to be as low as possible in order to concentrate the total power in main radiation, lowering the total signal interference and unwanted signal radiations. However, suppressing sideband radiations is not the only problem while designing a switched antenna array. The most important problem needs to be solved is shaping main radiation beam with required side lobe level which may be a ‘must’ of design procedure since it affects communication quality directly. Then, other parameters which affect the performance such as harmonic-level suppression, harmonic beamforming and spatial filtering (i.e., adaptive beamforming to use the power effectively and null-steering to reduce interference) need to be handled to have the best fit design for application. Thus, if all the constraints are taken into consideration at the same time, different optimization algorithms become practical once today’s computers’ calculation speed, performance and technology are getting better day by day.

As mentioned earlier, designing a switching array via a metaheuristic algorithm starts with the definition of the problem at hand which will directly appears as the cost function that the algorithm uses. After the problem definition, the switching strategy should be decided which directly affects the parameter vector construction. These two steps should be taken into account carefully. First, a possible mistake in problem definition causes a defected mathematical search space, actually the search space will not be mathematically defected, but wrong statements will lead completely different search space that does not related to the actual problem. The algorithm will produce some results but completely unrelated to the problem actually at hand. Second, an improper selection of the excitation strategy may lead the problem unsolvable. As a simple example to this situation, for a problem involving only the suppression of specific harmonic levels, the selection of a VAS or PS switching scheme makes the problem unsolvable, since the harmonic levels are bounded to each other in these schemes. Hence, the independent suppression of different harmonics is mathematically impossible for this problem using the mentioned time schemes. At this point, it is beneficial to decide which conventional array design parameters such as element orientation or excitation phase to use in addition to the ‘time’ parameter. Since many problems may be solved via only correct switching strategy, this selection is not crucial and the traditional parameters may be excluded from the problem. However, adding some extra degrees of freedom will relax the solution of the problem and in some specific cases such as main beam steering, inclusion of the excitation phase to the parameter vector is inevitable. From this point on, the rest of the procedure is the application of the optimizer to the problem. The key points in the TMA design may be summarized as

• Defining the problem correctly,

• Selecting the switching strategy correctly,

• Deciding whether the conventional parameters are necessary.

The first study dealing with harmonic suppression via an optimizer has been appeared by Yang et al. in 2002 [5]. In Yang’s study, the switch-on durations of a VAS schemed linear array has been tried to be optimized. Unlike conventional array design, in designing a switched array, the harmonic radiations should be taken into account and the question of ‘how?’ is an issue that needs to be considered. At early researches, the calculation of the sidebands was conducted as a sidelobe-like algorithmic procedure, and the result of the procedure is used in metaheuristic algorithms as in Yang’s study [6, 7]. In these researches, the pattern of the sidebands has been sampled in a certain precision, and the maximum of these sample set has been taken as the maximum sideband level (SBL) which has been used as the parameter to be minimized. The basic idea behind this approach is to lower the maximum level of the infinite harmonics below a certain level usually a certain communication threshold. However, her the main problem is: there exist infinite harmonics in number and the question is how many of them should be considered? For this issue, usually the first M harmonics are taken into account in practice and M depends on the problem at hand. However, there also exists a common opinion that conducting operation only on the first harmonic is enough to extract the maximum harmonic level since the Fourier coefficients should follow an envelope of a ‘sinc’ function. This idea is not entirely true since there exist situations that the first harmonic is not the highest especially in forced cases (e.g., the situation in [23]). SBL suppression is a generally applicable generic way to gather information about sidebands; however, calculating sidebands from array factor definition with a certain accuracy is a time- and system resource-consuming process even conducting operation in the first harmonic. If first M harmonics are considered, the time consumption is M fold for SBL extraction. Over and above, in azimuthal asymmetric cases such as planar and volumetric situations, the calculation time is enormously increasing exponentially in proportion to the azimuth resolution. Hence, there exists a solution time-accuracy trade-off in this technique, which always needed to be considered. An illustrative block diagram of this technique may be found in [31].

In the cost function of an optimizer, this SBL information may be directly added. As an example, for a linear array and a problem involving only sideband suppression, the cost function of the metaheuristic algorithm may be written as

where H∙represents the standard Heaviside step function. Here, Λ=maxx¯μm¯x¯θ−Ξd, μ∙is the array factor, Ξdis the desired sideband level and m¯represents the harmonic number vector where∀m∈ℤ+, m≤Mand Mrepresents the harmonic number taken into account. Additionally, it should be kept in mind that in switched array design, the parameter vector x¯should contain the switching parameters. For this kind cost function, the expected value of the cost of the solution is zero. In other words, the optimizer should stop when the cost value reaches zero and the parameter vector x¯that satisfies f=0regarded as a solution to the problem.

Extracting harmonic information to use in an optimizer from harmonic pattern samples is a generic way, but it is not the only way to gather information about sidebands. Calculating the power in sidebands and forming a bound function may also be used in sideband suppression problems. These two methods will be discussed in separate sections.

4.1. Harmonic power calculations in literature

In 2008, Bregains et al. published a paper that contains a closed form equation which defines an asymptotic approximation of the total power in infinite harmonic radiations for a VAS schemed linear array [25]. The main idea behind this approach is the total power of the radiated field, which may be expressed as the integral of the absolute squared array pattern over elevation under far field and κ≫1approximations. Since the infinite series of absolute squared complex Fourier coefficients is convergent which appears in the raw equation, the infinite summation reduces to a specific function involving switch-on durations. After taking the integration over elevation, the result becomes an elegant closed form equation which gives an asymptotic approximation of the total power in infinite harmonics. After this first paper, Poli et al. published a paper that extends this equation to a rectangular grid planar array [26], and Aksoy and Afacan rewrite the equation for arbitrarily distributed planar case [27]. The main idea behind these papers is the integration of the array factor in planar case, which may be expressed as ordinary Bessel functions. After this conversion, the rest of the way to writing the equation is similar to the original Bregains’ paper. One handicap of this equation is that it only holds for the VAS scheme. For a shifted case (i.e., PS scheme), the results of the equation have been considered as inappropriate, since the equation does not contain information about starting instants. In 2012, this situation is added and combined with original equation by Aksoy and Afacan [28]. The idea behind this attempt is writing the infinite summation of the absolute squared complex Fourier coefficients as the Riemann’s zeta function. From this moment on, the equation has been valid for VAS and PS switched linear and planar arrays. But the equation still did not include volumetric cases, and in 2015, Aksoy published a paper to close this gap [29]. The key point in volumetric calculations is solving a definite integral over elevation involving an ordinary Bessel function times a complex exponential function. Writing this multiplication using spherical Bessel and Legendre functions makes this integral easily solvable, and the result reduces to a zeroth order spherical Bessel function which is equal to a ‘sinc’ function with argument of Euclidean distances as in linear and planar cases. Hence, by the publication of [29], the equation has been taken its final and more general form for VAS and PS switched arrays. Here, it must be noted that the studies mentioned so far have been conducted using ideal pulses (i.e., rectangular pulses with no transition regions). In addition to the studies given above, a pulse model having transition region has been studied by Bekele et al. and may be a reference to practical pulse models [32]. The summary of the derivation of the general form of power equation (from this point on, it will be referred as SR equation) is given in the next subsection, and for more information, references given above may be followed.

4.2. Derivation of generalized total harmonic power calculation

In this section, the general form of the harmonic power equation will be briefly derived. Before beginning derivation, some remarks should be noted. First of all, the equation derived here is an asymptotic approximation which the far field and κ≫1approximations are assumed. Since the wave front of a spherically propagating wave may be assumed as planar in Fraunhofer region, this assumption helps to write the array factor in a simple summation consisting of planar wave front oriented phase shifts. Hence, in the region where the maximum phase error is below a certain tolerable level, the wave may be assumed as planar wave. Second, as stated earlier, the Poynting vector of a switched array may be written as in ordinary cases under κ≫1approximation which makes the power density calculations as in conventional arrays.

Under these approximations, assume that Nisotropic radiators are randomly oriented in a three-dimensional space. In this general case, the time average power at harmonics may be represented as

Here, μqθϕis being the array factor of isotropic sources and can be written for volumetric arrays as

where Wnand Gnqrepresent the complex excitations and complex Fourier coefficients, respectively. xn,ynand znstand for the coordinates of the array element in three orthogonal axis in standard Cartesian coordinates, k=2π/λois being the wavenumber and i=−1. In Eq. (18), the absolute squared array factor may be written as μqθϕ2=μqθϕμq∗θϕ,and after simple manipulations, PHturns out to be.

Here, the Fourier coefficients depend on the switching strategy. If an ideal PS time scheme consisting of rectangular pulses is considered, the switching function may be modeled as

where tn1and tn2represent the starting and finishing instants of the corresponding pulse. For this time scheme, the Fourier coefficients may be written as

where qstands for the harmonic number and the infinite summation of GmqGnq∗product may be examined for m=nand m≠ncases.

For m=n, this product may be written as

Hence, the infinite summation of this product turns out to be

where τn=tn2−tn1/Tprepresents the normalized on-time duration of nth element. Since, this series is convergent, using the Riemann’s zeta function where ∑x=1∞1/xn=ζn, the infinite series of coefficient product reduces to

For m≠nby using the fact that cosine is an even function and sine is odd, this infinite series turns out to be

By using relation for summable series, see [25, Appendix], Eq. (25) reduces to

For above equation, there exist 53 different situations in terms of equality and inequality conditions of switching instants. However, considering the fact that tn2>tn1and tm2>tm1should be satisfied always, these possibilities reduce to six. If these six equality and inequality possibility of starting and finishing instants of the pulses are considered, these six conditions may be written in one form given by

In this expression, the bar represents the duration of the intersection between n^th and m^th element’s pulses. For more detailed calculations, [28] may be followed. By substituting this result into Eq. (19), it becomes

where Ωθϕis being

and it should be solved. Solving this integral is not an easy task, and some manipulations should be conducted. It can be started from to write some exponential terms in terms of Bessel functions. To do that, the following fact may be used

By using this conversion and integral definition of ordinary Bessel functions with a=kxm−xnsinθand b=kym−ynsinθ, the Ωθϕcan be written as

Here, J0∙represents the ordinary Bessel function of order zero. By letting c1=kxm−xn2+ym−yn2and c2=kzm−zn2with Ωθϕ=2πΩ′θϕ, Eq. (31) becomes

Here, c1and c2are orthogonal to each other and can be represented as c1=Rcosψand c2=Rsinψwhere R=kR′and ψ=cos−1âz·rm−rn. Here, R′represents the Euclidean distance between corresponding elements (i.e., n^th and m^th elements under process) and rmxmymzmand rnxnynznstand for the position vectors m^th and n^th elements. By using this conversion, Eq. (32) can be written as

Once the Ω′θϕis written as in Eq. (33), the following relation may be used instead

Here, jn·and Pn·represent the n^th order spherical Bessel and Legendre functions, respectively. Multiplying both sides of Eq. (34) by Plcosθsinθand integrate over 0πresults in the following relation

where

In this expression, δn,lrepresents the Kronecker delta function defined by

Hence, for l=0, the equation Eq. (35) becomes

and Ωθϕcan be written as

If this result is substituted to Eq. (19) with Eq. (24) and Eq. (27), total power approximation may be finally written as

Once more, in Eq. (40), R=xm−xn2+ym−yn2+zm−zn2represents the Euclidean distance between corresponding elements, Irepresents the complex excitations and τrepresents the normalized switch-on durations for an Nelement volumetric array. For a more detailed proof, [29] may be followed, and for linear and planar case, see [25] and [27], respectively.

4.3. Using total power in suppression problems

The first usage of total power in harmonic suppression problems was conducted by Poli et al. in 2010 [10]. The main idea behind this technique is based on the idea of the total power reduction in harmonics will concentrate the power to main radiation, hence the harmonics become suppressed. Since the result of the equation does not meet the actual power in harmonics, the ratio of the harmonic power to the actual power is more meaningful in practice. Hence, usually the direct result of the Eq. (40) is not preferred to be used in an optimizer. Instead of using the direct result, the ratio of the harmonic power to the actual power is more meaningful which may be written by

where P0represents the power in main radiation given by

Without loss of generality for a problem concerning only harmonic reduction, the usage of Eq. (41) in a cost function may be written as

where ϖ=PH´I¯τ¯R¯−Θdwhile Θdis being the desired harmonic to main power ratio.

5.1. HLB concept

As shown in [30] that neither harmonic level reduction nor power reduction ensures the total suppression in terms of both power and communication level. Hence, a combined way may be a more suitable approach if the both level and power reduction is necessary. At this point, the question of ‘is there any simple method exist for the level calculation in order to avoid time consuming processes?’ may be emerged in mind, and it can be said that such studies exist in the literature as an answer to this question.

As mentioned earlier, the SBL method extracts the harmonic level information in an algorithmic way similar to sidelobe calculations. Since there exists no analytical solution to find a maximum point in an unknown sidelobe region, the sidelobe calculations are being conducted such an operation involving sampling the pattern in sidelobe region and finding its maximum. In contrast, the sideband calculations are being conducted in whole visible region, which makes a difference in both SLL and SBL calculations. In other words, since the SBL calculations are being operated over a complete elevation and azimuth space, in some cases, the maximum of a pattern may be extracted analytically. Since the maximum points of all individual harmonic patterns can be calculated, they form a maximum sideband level set, and the maximum point of this set bounds the whole harmonic maxima. Hence, this maximum appears as a bound function covering all individual maximum harmonic levels which are infinite in number. By this way, lowering the total bound ensures that all harmonic levels are below this level.

First attempt to calculate a bound is conducted by Aksoy for linear arrays with a VAS time scheme [23]. In this attempt, the Poynting vector of a TMA has been considered, and the pattern equations are written as the normalized time average power densities. After this first attempt, Aksoy and Afacan published a proof including the planar case using the same idea [24]. By this way, an overall harmonic bound for a VAS switched time modulated linear and planar arrays has been written as

where Ψstands for the normalized upper sideband bound, Inand τnrepresent the excitation amplitude and the normalized on-time durations of the n^th element, respectively. Also, ξnis the dynamic excitation defined by ξn=Inτnand Nrepresents the total element number in the array. Here it must be noted that, to clear possible misunderstandings, the excitation vector represents only the amplitude of the complex excitations and defined in region ∀I∈R+,0<I≤1. Similarly, the normalized switch-on durations are also defined in region ∀τ∈R+, where 0<τ≤1. Since this equation is a closed form expression, the total computation time for harmonic information is tremendously short as compared to the SBL method. Since it is in a shorter form, the computation time is slightly better compared to SR equality. A detailed comparison may be found in [31].

Like all switched array calculations, this technique also uses an asymptotic approach and valid only under far-field and κ≫1approximations and has its own advantages and disadvantages. The major advantage of this technique is that it ensures the harmonic level suppression among infinite sidebands in a very short time. However, currently, it is only applicable to linear and planar arrays with VAS time scheme, which may be considered as a handicap of this technique.

5.2. HLB applications

If the harmonic levels are wanted to be suppressed via the bound function in a metaheuristic approach, the result produced by the bound function may be used directly in the cost function of the optimizer. The aim of the optimization is to find a switch-on duration sequence to satisfy a desired level. Hence, the parameter vector should contain at least switch-on durations τ¯. For a uniform array, the excitation I¯may be taken as a vector whose elements are set to unity. Without loss of generality for a linear case with only sideband level, reduction is necessary for such a cost function that may be written as

where H∙stands for the standard Heaviside step function and λ=ΨI¯τ¯−Ξdwhile Ξdis being the predefined desired harmonic level. Here, the τ¯I¯set stands for a solution for f=0.

Moreover, the equality has been derived for the VAS scheme, but it can be used for any switching strategy involving one pulse in one period. On the other hand, the efficiency of the bound function in shifted cases does not appear as good as in the VAS scheme, since the equality does not contain information about starting instants. In practice, the actual harmonic levels are usually less than, or at least equal to, the results produced by the bound function. In VAS scheme, the difference between actual and the bound is usually getting smaller as the iteration counts, but the same statement does not hold true for the PS scheme. Hence, using the bound function in shifted cases causes a random walk on starting instants, and this situation may lead to some impractical results.

6.1. Variable pulse amplitude (VPA)

Instead of using “on–off” switching scheme, a novel approach is proposed by Aksoy in 2014 named as “Variable Pulse Amplitude (VPA),” which is based on switching each element between amplifiers [12]. The main purpose of the newly proposed method is preventing the array silencing, which means all the elements are switched to “off” position and no communication at that instant (e.g., see Figure 1). The binary “On – off” switching with modifications such as shifting and splitting allows the array to be silenced even for a very short time instant and the received or transmitted data may be missed at defined instants that all the elements are switched off in every period of switching. It may be seemed as minor possibility, but it results in additional optimization constraint to be checked whether there is an instant that all the elements are turned-off or not.

6.2. Harmonic suppression method

It is stated in previous sections that metaheuristic optimization algorithms do not guarantee to find the proper solution; however, they generally produce acceptable results. In addition, optimization time may be longer than usual according to host computer performance, optimization constraints and iteration count. On the contrary, a noniterative approach offers a fast and reliable solution with respect to metaheuristic methods, and it is preferable if exists.

In this part of section, a novel approach in order to suppress the sideband levels using variable pulse amplitude excitation scheme will be presented with an explanatory example.

6.2.2. A simple example

Let us assume a 10-element -30 dB Chebyshev array whose elements are located along z-axis with interelement spacing of d=0.5λ. For simplicity, excitation amplitudes of each element are selected equal to unity (i.e., In=1), and no progressive phase shift is applied to the elements (i.e., βn=0). Normalized pulse amplitude reversal time is selected as τn=0.45for each element, but it does not have any effect on maximum harmonic level since normalized pulse difference values (i.e., ∆n) are adjusted according to τnvalues.

Harmonic suppression ratio, γ, should be selected using 20log101/γ=δwhere δrepresents desired the maximum sideband radiation level in dB. It should be noted again that with the adjustment of γ, the maximum level of fundamental harmonic radiation is set to related dB value, thus, further sideband radiation levels are self-suppressed due to the nature of harmonic radiations.

In this simple example, maximum harmonic levels are selected as −30dBand−40dBso that harmonic suppression ratios should be γ=31.6and γ=100, respectively. Once −30dBmaximum level of sideband radiation is considered, main and first harmonic radiation patterns are given together in Figure 2a. Moreover, maximum levels of first 30 harmonic radiations are shown in Figure 2b. According to the calculation steps given in previous subsection, normalized pulse difference and adjacent pulse amplitude values that are used to obtain the results are given in Table 1.

n		1	2	3	4	5	6	7	8	9	10
γ=31.6	Kn1	0.258	0.430	0.669	0.878	1.000	1.000	0.878	0.669	0.430	0.258
	Kn2	0.233	0.389	0.605	0.794	0.905	0.905	0.794	0.605	0.389	0.233
	∆n	0.025	0.041	0.064	0.084	0.095	0.095	0.084	0.064	0.041	0.025
γ=100	Kn1	0.258	0.430	0.669	0.878	1.000	1.000	0.878	0.669	0.430	0.258
	Kn2	0.249	0.417	0.648	0.851	0.969	0.969	0.851	0.648	0.417	0.249
	∆n	0.008	0.013	0.021	0.027	0.031	0.031	0.027	0.021	0.013	0.008

Table 1.

Normalized pulse difference and adjacent pulse amplitude values for γ=31.6 and γ=100.

If −40dBof maximum harmonic level is taken into consideration (i.e., γ=100), Figure 3a and Figure 3b display the radiation patterns for q=0together with q=1and maximum levels of first 30 harmonic radiations, respectively. ∆n,Kn1and Kn2values are also shown in Table 1.