© R. Palmeri and T. Isernia, published by EDP Sciences, 2021

1 Introduction

The ever-growing interest in artificial materials, including metamaterials/metasurfaces and photonic crystals, is related to their very intriguing properties and possibilities allowing very numerous and even exotic applications in all ranges of the electromagnetic spectrum including optics [1], microwaves [2] and millimeter-waves frequencies [3].

In this contest, devices based on electromagnetic bandgap (EBG) structures result very attractive and are exploited for a wide class of applications, ranging from radiation [4] and propagation [5] applications, to particle accelerators [6]. Amongst them, in this contribution we focus on EBG waveguides.

EBG structures [7] are periodic arrangement of inclusions usually adopted to improve the performance of devices or to realize ad-hoc functionalities for a given device. In particular, EBG structures have the fascinating property of suppressing field propagation within a frequency band (depending on the structure itself) while allowing it in other bands. If a defect is introduced in the periodic arrangement by removing some inclusions, such a defect will support an electromagnetic mode and hence it is possible to generate localized or guided fields [7].

Despite the very large number of analysis tools for this kind of structures, both at ‘macroscopic’ (like homogenization techniques) and ‘microscopic’ (as finite-element-methods) level, as well as at a ‘mesoscale’ level [8], formalized synthesis and design tools are quite limited in literature. Basic approaches are founded on trial-and-error cyclic processes, making extensive use of numerical simulations that are occasionally combined with genetic-algorithms-based optimizations [9]. Obviously, because of the huge number of involved parameters and the corresponding high number of trials, this is a far-from effective way to proceed. More rigorous design approaches are those based on inverse design, such as the inverse homogenization problems [10,11] which are mainly devoted to synthesizing effective electromagnetic properties of unit cells. Unfortunately, the arising dielectric profiles exhibit a complex distribution which is difficult to fabricate [10], so that fabrication-oriented tools are needed. In this respect, topology optimization [12] based design tools have been widely adopted during the years. As the name suggest, this kind of techniques aims at modifying the topology of the device at hand, that is, the geometrical and e.m. characteristics of the object being designed, in such a way some assigned specifications are satisfied. Developed design tools based on topology optimization [13–17] or the like [18,19] are able to realize binary devices wherein just two materials (one of which is usually air) are involved. Other types of design tools based on group theory [9] and the emerging concept of deep neural network [20] have been also proposed.

Recently, the design of artificial materials through an inverse design approach has been introduced within the framework of inverse scattering problems for antennas [21] and cloaking [22,23] applications, showing very promising performance compared to homogenization-based approaches, but also evidencing limitations as far as the computational burden is concerned. In fact, the adopted tools in [21] require very small mesh elements to correctly represent small scatterers and hence the field they scatter.

As an alternative to above techniques and with the aim of overcoming relates issues, we have recently proposed, for a completely different problem concerned with antennas (see below), a novel inverse design approach based on inverse scattering co-adjuvated by the use of the so called Scattering Matrix Method (SMM) [24]. The SMM is a very accurate and fast analysis tool for structures made of multiple dielectric or metallic scatterers. In fact, it allows to exploit analytical (or off-line numerical) results related to the ‘response’ of each single inclusion to have an accurate modelling of the overall device. Notably, it operates somehow at a ‘mesoscale’ level, as the modelling is not based on the microscopic level of numerical simulations acting on a dense grid, and avoids the macroscopic but possibly approximated way of thinking based on homogenization techniques. Notably, the basic SMM method can be eventually improved (see [8] for more details) by introducing aggregations of elementary scatterers and computing the scattering matrix of the larger scatterer which is obtained.

As we are going to summarize in the next section, the scattering matrix, through proper expansions of the incident field (acting as input) and scattered fields (acting as output), allows to take into account in a very effective fashion the ‘response’ of each single scatter at hand. Then, one is left with a system of equations governing (just) the interactions amongst the different scatterers, while interactions internal to the single scatterers constituting the device are taken into account through the corresponding scattering matrix. Of course, the number of equations and the dimension and structure of the scattering matrices depend on the number of scatterers as well as on their geometrical and electromagnetic characteristics.

With the aim of pursuing the design of EBG devices, the SMM formulation can be conveniently turned into an inverse problem. More in detail, the role of data and unknowns is inverted so that the scattering matrices of each inclusion can be considered the new unknowns of the problem, while expansion coefficients of the scattered fields must be determined in such a way to obey assigned field specifications, depending on the device at hand.

Very recently, we have proposed our Inverse SMM (I-SMM) tool in the framework of innovative antennas design. In fact, we have successfully applied it in [25] for the same design problem considered in [21]. Note that in the antenna applications above, the specifications were derived starting from an optimal synthesis of the coefficients of a representation of the scattered and total field, and the pursuing of such a field is then the goal of the inverse scattering step, so that such a field substitute the field measurements which are the usual data of inverse scattering for reconstruction problems.

In this contribution, we consider a completely different problem with a completely different kind of design goal. As a matter of fact, EBG guiding structures are considered rather than antennas, and the actual goal is the optimization of the Voltage Standing Wave Ratio (VSWR) along the structure in such a way to avoid undesired reflections in case of bends, rather than the realization of a given field.

Notwithstanding the completely different problem and the different kind of goals, it is shown that scattering matrices, as well as the kind of (state and data) equation which are the core of inverse scattering procedure can still be used.

As well known, the electromagnetic response of a class of EBG devices can be understood by analysing the bands diagram which depends in turn on the kind of lattice, polarization, and the used defects (if any). On the other side, such a ‘basic’ way of thinking does not exploit a number of degrees of freedom which can be gained by changing (for example) the permittivity or dimensions of (some of) the inclusions. For example, which is the basic problem we deal with in the following as a first demonstration of potentialities, one could use a few degrees of freedom regarding a limited number of inclusions to improve the behaviour of an EBG waveguide. In fact, by starting from a basic structure, inverse design procedures are exploited to modify (part of) the lattice for an improving of the VSWR (i.e., saying it in other words, for avoiding waves coming back along the defect). A similar problem has been dealt with in [13–17]. In particular, we apply the I-SMM to improve guiding performance of straight and bent EBG waveguides. Our procedure resembles the one in [26], wherein the multiple scattering theory (MST) is adopted. Both approaches are based on a harmonic expansion for fields and their formulation is general, so they can be formalized for different kind of inclusions; the main difference is that, opposite to the SMM, the MST exploits the T-matrix to link the scattered field to the electromagnetic properties of the scatterer. As far as the design procedure is concerned, an important difference can be highlighted amongst our approach and the one based on the MST. In fact, in [26] a global optimization is performed for the solution of the inverse problem. As a consequence, to limit the computational burden, several a-priori information must be used. In particular, in [26] the authors not only set arrangement and dimension of scatterers, but also their electromagnetic parameters, the remaining degrees of freedom being just the presence or not of the scatterers themselves. In our approach instead the inverse problem is solved through a conjugate gradient based algorithm that is able to handle both geometric and electromagnetic features, assumed that a favourable starting point can be provided for the optimization step.

The remainder of the paper is structured as follows. In Section 2 the SMM is briefly recalled while in Section 3 its inverse formulation is introduced. Section 4 is devoted to present and assess the procedure for designing optimized EBG waveguides. Conclusions follow.

For the sake of simplicity, 2-D scalar fields are considered and the time harmonic factor e ^jωt is assumed and dropped.

2 Basics of the Scattering Matrix Method (SMM)

Let consider the geometry depicted in Figure 1 wherein the ℓ-th cylinder with cross section C _ℓ is located at (r _ℓ, θ _ℓ) with respect to a global reference system Oxy. Let also E _inc be a given incident field impinging on the set of scatterers. By exploiting the cylindrical harmonics expansion of fields [27], the field scattered by the set of inclusions in the regions of space outside the single scatterers reads as:(1)wherein N is the total number of scatterers, k _b is the wavenumber of the background medium and is the m-th order Hankel function of second kind. Notably, representation (1) is expressed in the reference system of the ℓ-th scatterer (see Fig. 1).

As soon as the observation point is a little bit outside an inclusion C _ℓ, the corresponding internal summation can be safely truncated¹ to a finite number, say M _ℓ, which also depends on the single inclusion maximum dimensions [28].

Accordingly, the total field can be represented as:(2)

From equation (2) it is possible to note three contributions: the first one is the primary incident field; the second addendum is a summation of scattered fields from all inclusions but the ℓ-th one, that accordingly acts as a secondary incident fields on it; finally, the last term simply represents the field scattered by the ℓ-th cylinder as defined in equation (1).

In equation (2), the fields representations are still expressed in the reference system linked to the pertaining cylinder. In order to have a field representation in a unique global reference system (e.g., the ℓ-th one), the Graf's formula [29] for the Hankel function can be used. We finally obtain what follows (see [24] for more details):(3)

a _ℓ,m coefficients being the following compact expression:(4)wherein Q _ℓ,m encodes terms arising from the harmonic expansion of the primary source, while T _ℓ,i,m,q takes into account mutual coupling effects amongst scatterers [24], that is, it is related to the secondary incident field. For the sake of completeness, analytic expressions of T _ℓ,i,m,q and Q _ℓ,m taken from [24] are given in Appendix.

If b _ℓ is the vector encoding the field representation coefficients while a _ℓ is the vector encoding the source representation coefficients, and according to the fact that a finite number of parameters can be used to represent both fields and sources [28], for any fixed scatterer one can conveniently design a finite dimensional ( M _ℓ ×  M _ℓ) scattering matrix S _ℓ encoding all relevant information of the scatterer, such that:(5)

Hence, after multiplying equation (4) by the scattering matrix, we finally achieve the linear system in the unknown b _ℓ:(6)

that formally describe the scattering model at hand.

Once equation (6) is solved, the scattered field in the region of space outside scatterers can be finally evaluated through equation (1) [24].

Fig. 1 Geometry sketch for the problem at hand.

3 The Inverse-SMM (I-SMM)

The aim of an inverse design problem is to infer optimal properties on certain parameters, within a set of constraints, to accomplish assigned specifications at best. In the present case, the parameters to be optimized are the geometric and electromagnetic features of inclusions in such a way the scattering from the optimized set is able to satisfy some given specification on the scattered or total field.

If a nominal assigned field to be pursued is given, in system (6) both the scattering matrices S _ℓ and the coefficients b _ℓ are unknowns, and these latter must be determined such to fit the given assigned field as well, that is, such that:(7)

Accordingly, in the I-SMM design problem, S _ℓ and b _ℓ must be optimally retrieved to contemporarily satisfy the model (6) and the data constraints (7). Notably, such an inverse problem is not trivial to solve because of its non-linearity that could cause the optimization problem (based on a gradient scheme) to be trapped into local minima. However, as long as a good matching is gathered with both the model and the data, sub-optimum solutions can be accepted.

A further issue to be taken into account concerns the existence of a solution for the problem, especially in terms of a final actual device. More in detail, even if scattering matrices S _ℓ can be found as solution of the inverse design problem, the existence of corresponding physical reliable scatterers is not ensured. For this reason, a-priori information on the scatterers are required to gather a feasible solution of the problem. This allow, for instance, to look for diagonal matrices if cylindrical scatterers with circular symmetric cross-section are dealt with. Also, one can look for a single S matrix if the scatterers are identical. Moreover, we can directly search for the parameters actually determining the entries of the matrices if the electromagnetic features of the inclusions are supposed to be of a given kind, namely if we are dealing, for instance, with homogeneous dielectric or metallic inclusions, single layer or stratified inclusions, and so on. In this case, in fact, analytic expressions for the scattering coefficients are available [25].

At the end of the day, thanks to the fact one is reducing as much as possible the number of unknowns of the problem and by taking advantage from convenient starting points derived from a general physical understanding or homogenization based suggestions, the minimization procedure will ensure actual solutions of the design problem, without resorting to global optimization algorithms, with the inherent advantages in terms of computational burden.

4 An ‘ad hoc’ procedure for optimal EBG waveguides design

The SMM and I-SMM presented above are general formulation for the scattering model dealt with, namely they still apply whatever the arrangement of the scatterers, their electromagnetic (e.g., dielectric, metallic, stratified) and geometric (cylinders with circular or elliptic cross-section) features, whatever the background medium and the primary source². Thereinafter, in order to deal with EBG devices, we consider circular lossless dielectric scatterers with radius ρ _ℓ and constant relative permittivity ϵ _ℓ, embedded in air and arranged in a periodic (triangular or square) lattice with spacing d.

On the basis of ρ _ℓ, ϵ _ℓ,d values and lattice type, the EBG structure presents a gap within a definite frequency band, that represents a forbidden region for fields to propagate. When a line of inclusions is removed from the periodical lattice, a new guiding effect can be obtained. As a matter of fact, the line defect turns into a waveguide for an electromagnetic wave with a frequency within the band gap of the (non-defective) structure. An example of this mechanism is shown in Figure 2, from which it is possible to see that an impinging field at f = 288.46 THz cannot propagate within the full structure because it belongs to the band-gap frequencies, while it propagates when the line defect is realized.

However, the defective EBG structure is able to confine and somehow guide the field, but an unsatisfactory VSWR is obtained. In fact, wave reflections due to the discontinuity in the space produce oscillations of the field inside the guiding region. Notably, this behaviour occurs whatever the kind of waveguide (straight, bent).

Hence, we can exploit the I-SMM framework to optimize (some of) the inclusions in order to improve VSWR and get a better design. In particular, we fix the arrangement of the scatterers and their dimension, which is supposed to be identical; then, we look for the optimal permittivity value of inclusions belonging to the so-called ‘design region’ (DR) allowing a minimum ripple for the field amplitude within the waveguide and, in particular, in the so-called ‘shaping region’ (Λ_shaping ).

Hence, differently from [25], the inverse design problem can be finally formulated as follows:(7.1)

Subject to:(7.2) (7.3)wherein:(8)constraint (7.2) allows to fit the scattering model, while constraint (7.3) enforces a symmetry in the permittivity of the structure.

In the following subparagraphs an assessment of the inverse design procedure is performed for straight and bent waveguides. The considered primary source is a plane wave from the left side of the structure. As we are interested in perspective to reconfigurable devices, the permittivity of the inclusions are the unknowns of our (fixed dimensions) elementary bricks. Of course, a completely similar approach is possible by looking for dimensions (with a fixed permittivity) or both.

As far as the starting point of the optimization problem is concerned, in the specific problem dealt with herein (which can be seen as the problem of optimizing the field which is already propagating within waveguides), the standard defective structure is an obvious (and effective) choice. Also note that constraints on the admissible range values can be considered in the optimization problem to take into account actual feasibility issues.

Fig. 2 Electromagnetic field propagation mechanism for EBG devices. (a) Band-map in logarithmic scale. Amplitude of electric field with λ = 1.04 μm for (b) non-defective and (c) defective EBG structure.

4.1 Straight waveguide

As a first example, let us consider the geometry shown in Figure 3a. The full EBG structure is composed of N=168 rods with radius ρ =0.075 μm and permittivity ϵ = 12.25 arranged on a triangular lattice with spacing d = 0.45 μm. The central row composed of 21 rods is removed to realize a straight waveguide. The amplitude of the total electric field on domain for the defective structure is shown in Figure 3b. As previously outlined, the field amplitude achieved within the channel is oscillating, corresponding to a Ripple value equal to 0.3821. To improve such a parameter, the design problem (7) is solved by considering rods on the edge of the channel as DR, see Figure 3a. The optimized structure and the corresponding amplitude field are shown in Figures 3c and 3d, respectively. As it can be seen, a proper adjustment of some inclusions' permittivity allows to reduce the Ripple to 0.0173 (with a comparable value of the average field) and hence to improve the guiding effect.

Fig. 3 Optimization of a straight waveguide. (a) Permittivity map of the defective EBG device where the DR is superimposed with light grey rectangle and (b) the corresponding electric field amplitude distribution. (c) and (d) the same as (a) and (b) but for the optimized EBG device through the I-SMM [ϵ opt =8.4,10.8,13.2,12.4,10.1,11.4,13.4,10.7,10.2,13.7,12.7,9.6,12.3,16.4,10.9,9.4,16, 16.7,9.3,11.5,18.1].

4.2 60° bend waveguide

The second example concerns the optimization of a double 60° bend waveguide. In this case, the total number of scatterers is N = 315, while the structure's geometry as well as the cylinders' characteristics are unchanged. By removing 17 rods as in Figure 4a the amplitude field distribution of Figure 4c is achieved, for which the corresponding Ripple, evaluated for horizontal and diagonal channels separately, is summarized in Table 1. By optimizing the device by using the DR highlighted in Figure 4a, we obtain the device in Figure 4b and a reduced Ripple within the channel, see Figure 4d and Table 1.

Fig. 4 Optimization of a double 60° bend waveguide. (a) Permittivity map of the defective EBG device where the DR is superimposed with light grey rectangle and (c) the corresponding electric field amplitude distribution. (b) and (d) the same as (a) and (c) but for the optimized EBG device through the I-SMM [ϵ opt =83,77.7,79.2,78.3,79,79.5,82.4,79.9,80,79.7,79.5,82.7,80.1,79.7,78.8,79.3,79,87.1].

4.3 90° bend waveguide

As last example we deal with a double 90° bent waveguide. The rods of the EBG structure (N = 225) have the same parameters as before but for the lattice that has a square arrangement in this case. By removing 20 inclusions as shown in Figure 5a we realize a waveguide with sharp bends, for which the field distribution of Figure 5b is obtained. A similar device is considered in [16] but with rods arranged in a triangular lattice, which is unusual for this kind of sharp bends. Conversely, in [14] the same authors investigate a single 90° bent waveguide but just the corner area is optimized.

In the present case the DR still corresponds to the rods on the channel's edges. As can be apprised from Figure 5d, that is the field pertaining to the optimized structure of Figure 5c, the design through the I-SMM allows to reduce the Ripple of the guided field. See Table 2 for a quantitative comparison.

Fig. 5 Optimization of a double 90° bend waveguide. (a) Permittivity map of the defective EBG device where the DR is superimposed with light grey rectangle and (b) the corresponding electric field amplitude distribution. (c) and (d) the same as (a) and (b) but for the optimized EBG device through the I-SMM [ϵ opt =12.2,14.3,12.2,12.5,11.9,12.5,17.7,11.3,12.4,7.4,8.3,13.7,14.4,14.2,12.9,13.2,12.8,12.9, 11.8,21.6].

5 Conclusions

In this paper the inverse design of artificial-materials-based devices has been addressed. In particular, electromagnetic bandgap (EBG) devices have been considered. Differently from approaches relying on full-wave methods or trial-and-error procedures, inverse design based techniques represent powerful tools to obtain optimal devices realizing desired properties in a more effective fashion. In this contribution, an inverse design tool for EBG devices based on the scattering matrix method (I-SMM) has been proposed. More in detail, the SMM has been adopted to model EBG devices (which are a set of scatterers actually) and hence to develop an ad hoc procedure to determine optimal parameters describing the inclusions composing the device in order to accomplish desired specifications. The proposed I-SMM represents a general framework and makes use of simple geometry structures with multiple materials that could be possibly tuned to accomplish with different specification (i.e., with reconfigurable devices). The proposed approach has been successfully assessed by optimizing the guided electromagnetic field within straight and bent waveguides. In particular, by optimizing the permittivity value of inclusions located on the waveguide's edges, our tool is able to control oscillations of the field while keeping an adequate transmitted power.

As we are interested in perspective to reconfigurable devices, we have focused herein on fixed spacings and fixed dimensions of the inclusions. Interestingly, the model we have used can be both immediately used for the case of fixed permittivities-variable dimensions as well as extended to the case of variable spacings (by acting on the T matrices rather than on the S matrices).

Acknowledgments

This work was supported in part by the Italian Ministry of University and Research under the PRIN research project “CYBER-PHYSICAL ELECTROMAGNETIC VISION: Context-Aware Electromagnetic Sensing and Smart Reaction,” prot. 2017HZJXSZ.

Appendix

The present Appendix is devoted to report the analytic expressions of entries of T and Q matrices involved in the SMM model of equation (6).

Accordingly, in the chosen coordinates system, they read:(A1)

wherein represents the mutual position between the chosen ℓ-th coordinate system (O ^ℓ x ^ℓ y ^ℓ) and the coordinate system of the generic i-th inclusion, while, depending on the primary source, a different Q _ℓ,m formula must be derived.

In particular, if a plane wave E _inc  (x, y) = e ^{−jk _b  (x sin θ _inc −y cos θ _inc )} is considered (θ _inc being the angle with respect to the y-axis) the pertaining expression for Q _ℓ,m is the first one in (A2), while if a line current I _s located at (r _s , θ _s ) is supposed as primary source, that is, , the second relation in (A2) holds, being the mutual position between the chosen ℓ-th coordinate system and the source position. Other kind of sources can be obviously accommodated as well.(A2)

References

All Tables

All Figures

	Fig. 1 Geometry sketch for the problem at hand.
In the text

	Fig. 2 Electromagnetic field propagation mechanism for EBG devices. (a) Band-map in logarithmic scale. Amplitude of electric field with λ = 1.04 μm for (b) non-defective and (c) defective EBG structure.
In the text

Fig. 3 Optimization of a straight waveguide. (a) Permittivity map of the defective EBG device where the DR is superimposed with light grey rectangle and (b) the corresponding electric field amplitude distribution. (c) and (d) the same as (a) and (b) but for the optimized EBG device through the I-SMM [ϵ opt =8.4,10.8,13.2,12.4,10.1,11.4,13.4,10.7,10.2,13.7,12.7,9.6,12.3,16.4,10.9,9.4,16, 16.7,9.3,11.5,18.1].

In the text

	Fig. 4 Optimization of a double 60° bend waveguide. (a) Permittivity map of the defective EBG device where the DR is superimposed with light grey rectangle and (c) the corresponding electric field amplitude distribution. (b) and (d) the same as (a) and (c) but for the optimized EBG device through the I-SMM [ϵ opt =83,77.7,79.2,78.3,79,79.5,82.4,79.9,80,79.7,79.5,82.7,80.1,79.7,78.8,79.3,79,87.1].
In the text

Fig. 5 Optimization of a double 90° bend waveguide. (a) Permittivity map of the defective EBG device where the DR is superimposed with light grey rectangle and (b) the corresponding electric field amplitude distribution. (c) and (d) the same as (a) and (b) but for the optimized EBG device through the I-SMM [ϵ opt =12.2,14.3,12.2,12.5,11.9,12.5,17.7,11.3,12.4,7.4,8.3,13.7,14.4,14.2,12.9,13.2,12.8,12.9, 11.8,21.6].

In the text