Coefficients for calculation of the scaling factor Ψ.
For over four decades the evolution of electronic technology has followed Moore's low where the number of transistors in an integrated circuit (IC) approximately doubles every two years. This naturally increases the number of internal interconnections needed to complete the system. The increase in chip complexity is achieved by a combination of dimensional scaling and technology advances. A variety of chip types exist, including memory, microprocessors and application specific circuits such as System-on-Chip (SoC). Since it is expensive to fabricate the large and simple passive components such as on-chip capacitors and inductors on the same die as the active circuits, it is desirable to fabricate these on separate dies then combine them in System-in-Package (SiP). The main advantage of SiP technology is the ability to combine ICs with other components, including passive lumped elements already mentioned but also antennas, high speed chips for radio frequency communication etc., into one fully functional package. The high complexity of SiP brings many challenges to the design process and physical verification of the system. In many cases the design process relies on detailed 3-D numerical electromagnetic simulations that tend to be slow and computationally demanding [1,2,3] in many cases limited by the available computer memory capacity and computational speed. Therefore, directly including the detail of the dense interconnect networks into the numerical model is demanding due to the amount of memory required to hold the detailed mesh, and numerical penalties associated with small mesh cell sizes relative to the wavelengths of the signals being modelled.
From a package-level point of view, the on-chip interconnects can be seen as a mixture of metal inclusions located in a host dielectric [4-8]. Most of the current studies based on numerical analysis of 2-D or 3-D structures with two constituents show that the effective properties of the mixture strongly depend on the volume fraction, its geometrical profile and spatial orientation in periodic or random arrangements [9,10]. It has been shown that the macroscopic properties of dielectric-only mixtures can be represented by a homogenous dielectric with an effective permittivity that is determined using an empirical mixing model . Metal-dielectric mixtures have been less thoroughly explored, with work limited to treating spherical or ellipsoidal metal inclusions [12,13]. The approach was extended [4-8] to cope with rectangular cuboid metal inclusions representative of on-chip interconnect structures. The use of a single fitting parameter was retained, and is calculated for a wide range of aspect ratios (0.6 – 3), dielectric host materials (1 – 11.7), metal fill factors (0.2 – 0.6) and signal frequencies (1 – 10 GHz) that are likely to be of interest to System-in-Package designers. Here a simplified empirical mixing model defined for a narrowed down range of aspect ratios (1.4 – 3) which accounts for the interconnect geometries is presented. The model is at the same time straightforward and more accessible as well as more accurate. The accuracy improvement is related to the neglected range of low aspect ratios where the scaling factor Ψ has more rapid increment [4, 5].
Interconnects often form regular gratings (e.g. bus structures), hence an infinite metallic grating in a homogenous dielectric host, as shown in Fig. 1., represents a straightforward but broadly applicable model.
Wide parameter space for the dimensions is considered, as follows. The interconnect pattern density is initially limited by the design rules within 20% - 80% metal fill . A statistical analysis of the pattern density of a real chip showed that, typically, the maximum pattern density in actual metal layers does not exceed 60% . Thus metal fill factors f in the range 0.3 - 0.6 is considered. While metal layer height is fixed for any given layer in any given process, track width is less restricted. Aspect ratios (xAR) are continuing to increase as technology develops , hence structures with narrowed down [4,5] values within 1.4 ≤ xAR ≤ 3 are studied. Due to the growing use of low-k dielectrics a host materials with permittivity εe in the range 1 ≤ εe ≤ 11.7 are considered. The interconnect pitch Λ is often measured in micrometres or nanometres, whereas the wavelength λ of the clock signal is typically measured in centimetres. Thus, we can expect to successfully apply an appropriate effective medium approximation because the condition Λ << λ/4 is met . Here, the interconnect pitch is fixed Λ = 100 μm for the sake of clarity in illustrating the method. While we include the effect of the interlevel dielectric, we deliberately neglect the bulk substrate and all but the top layer of interconnects.
The modified Maxwell Garnett mixing rule is used
where εi and εe are the dielectric functions of the inclusion and host material respectively (here, a metal and a dielectric), Ψ is a constant relating the fields inside and outside the inclusions (typically Ψ = 3 for spherical inclusions), f is the filling factor or ratio of the volume of the inclusion to the total size of the unit cell . In earlier work a limited example of such an approach for a fixed value of Ψ applicable to a single structure was presented [6-8]. Further, the approach was expanded by developing a compact equation to calculate the appropriate value of Ψ for a broad range of parameters [4,5]. Here a new empirical model is presented where the considered values of aspect ratio are within 1.4 – 3. The frequency dependent dielectric function of a metal inclusion εi(ω) can be expressed by a Drude model [13,19]
where ω is the frequency of interest, ωp is the plasma frequency and γ is a damping term representing energy dissipation. Despite wide spread use of copper interconnects for the intermediate levels of the interconnect stack, aluminium is often used for the global wiring with which we are concerned, and has ωp = 15 eV and γ = 0.1 eV . Note that the energy is related to the free space wavelength λ0 by ω = 1.24x10-6(λ0)-1. This model was used in both the analytical and numerical calculations.
3. Empirical model
Rigorous coupled wave analysis (RCWA)  was performed for the TM polarisation with the electric field vector ETM coplanar with the grating vector K as shown in Fig. 1. The structure was illuminated by a plane wave with incident angle -89 ≤ θ ≤ 89º and free-space wavelengths λ = 30 cm, λ = 10 cm, λ = 6 cm, and λ = 3 cm. Hence, the adjusted height h of the homogenised layer does not simplify the model the height of the homogenised equivalent layer was kept the same as the grating. The reflection and transmission coefficients for the homogenised structure were calculated using an analytical formula defined for a stratified medium comprising a stack of thin homogenous films .
It is not necessary to make Ψ dependent on the host dielectric as this is already accounted for explicitly in Eq. (1). It was verified by numerical experiment. The scaling factor Ψ in the narrowed down range of aspect ratios was observed to have a linear dependence on this parameter, therefore the general form of the empirical model was chosen:
where the coefficients α, β are determined by linear regression from data obtained from nearly 5000 simulations spanning a four dimensional parameter space. The coefficients are represented as a linear function of metal fill factor by
where k1, k2 are well approximated by
where ν is the frequency (in units of GHz), and factors k11, k12, k21, k22 are presented in Table 1.
The fit of the model was assessed using a linear least square method. Figure 2 illustrates the good agreement between Ψ obtained from a ‘brute force’ fitting algorithm and that from our linear approximation for an example grating structure. The grating has f = 0.5, and the illumination frequency is 5 GHz.
In Fig. 3 the reflection coefficient obtained using RCWA and homogenised model for a sample of 27 different structures are depicted. For frequencies in the range 1 GHz ≤ ν ≤ 10 GHz the error between RCWA results for the detailed structure and those obtained for the homogenised structure is less than 2.5% for reflection coefficient and 0.2% for the transmission coefficient (not shown here) when 0.3 ≤ f ≤ 0.6 and θ ≤ ±30º. The accuracy of the model for reflection coefficient calculations tends to improve with an increase of the metal fill factor. For structures with metal fills 0.4 ≤ f ≤ 0.6 the error varies between 0 – 1.5%. When the model is applied, without modification, to interconnects with a trapezoidal cross section, sometimes found in fabricated structures, the error remains similar as for the empirical model presented in [4,5] and is below 5% for sidewalls with angles of up to 5º and incident angle up to 30º.
The calculated effective permittivity varies according to the particular mixing rules used to analyse a given mixture. However, there are theoretical bounds to the range of calculated effective permittivieties. For the compound of the two dielectrics the effective permittivity calculated from the Maxwell-Garnett mixing rule has to fall in between the following bounds 
These bounds are also called Wiener bounds. The upper limit for the effective permittivity εeff,max is defined for a layered material with boundaries between inclusions and host dielectric parallel to the field vector. The lower bound εeff,min is obtained for the case where the field vector is perpendicular to the boundaries between inclusions and host. Since the Wiener bounds are defined for anisotropic mixtures, stricter bounds, Hashin-Shtrikman bounds, have been defined for the statistically homogenous, isotropic and three dimensional mixtures. The upper and the lower bounds are as follows
where it is assumed that εi > εe. The lower limit corresponds to the Maxwell-Garnett mixing rule whereas the upper limit is the Maxwell-Garnett rule for the complementary mixture obtained by transferring the constituents: εi → εe, εe → εi, f → 1 - f.
It was verified that in a set of about 5000 simulations of the grating structure run to define the empirical model all effective refractive indices (neff = √εeff) are well within the Wiener bounds. Nevertheless, the predicted neff has values close to the lower limit. This is related to the specific alignment of the grating structure (single layer of interconnects) and the angle of incidence wave. Such regular and linearly distributed arrangement of the inclusions with the field vector perpendicular to the grating surface results in an effective permittvity from the bottom range of the possible values defined by Wiener bounds. The upper limit is several orders higher in magnitude, hence even if satisfied, for the purpose of the analyses of this particular grating structure it can be lowered by replacing it with the Hashin-Shtrikman lower limit. It is illustrated in Fig. 4, for a random structure, that the real parts of neff obtained from the empirical model are within the lower limits of the Wiener and Hashin-Shtrikman bounds.
The more strict Hashin-Shtrikman bounds overestimate the obtained values of neff. Hence these limits are based on the Maxwell-Garnett mixing rule for the complementary mixtures and the lower limit is just the classical Maxwell-Garnett rule with εi > εe. Therefore, for the analysed interconnect grating structure it can be assumed that the upper bound for the effective refractive index is the classical Maxwell-Garnett rule whereas the lower bound is the Wiener lower limit.
4. Experimental validation
Experimental validation was carried out by the free space measurement of S-parameters of an air-copper grating structure (Λ = 500μm, f = 0.3, AR = 1) attached to a Rogers 4350 dielectric plate (thickness 762μm, εr = 3.66) illuminated by a plain wave. A pair of horn antennas operated at the X-band (8.2 – 12.4 GHz) frequencies with Teflon’s hemispherical lenses connected to the network analyzer was used. The plane wave illumination focused on a relatively small area was achieved by the special equipment arrangement. A free space calibration method along with smoothing procedure was implemented in order to eliminate systematic errors occurring in the measurement data [22,23]. A 2-D finite difference time domain (FDTD) was defined as shown in Fig. 5, with the case of detailed grating structure in Fig. 5(a) and its homogenised equivalent in Fig. 5(b).
The domain size was 60 cells in x by 1 cell in y and 8203 cells in z direction. The space increment in both directions was set to 5μm and it was ensured that the domain size in the z-direction was at least a half wavelength from each of the absorbing Perfect Matched Layer (PML) boundaries as the behaviour of these boundaries is not reliable in the presence of evanescent fields. The grating structure is periodic in x-direction, with one period of the grating defined in the domain. The structure was illuminated by a wave with frequencies within 8 – 12 GHz in steps of 1 GHz applied at the top of the domain and propagated in the negative z-direction.
The homogenised equivalent structure was obtained by replacing the grating layer with a solid dielectric. The dielectric properties were calculated from the modified Maxwell-Garnett mixing rule. The value of scaling factor Ψ was empirically found as due to the structural difference between experimental settings and the structure studied in order to define the empirical model, the straightforward application of the empirical model underestimated factor Ψ. It was verified that Ψ, when equal to 2.5, gives good approximation of the calculated effective permittivity εeff.
The FDTD and analytical calculations, based on characteristic matrix method, of reflection and transmission coefficients for the gratings with structural period 500μm and its homogenised equivalent are plotted in Fig. 6 and Fig. 7 along with measured return and transmission losses respectively. In order to validate the proposed approach for a wider range of frequencies, 1 - 18 GHz, numerical calculations were performed using finite element method (FEM) . This analysis shows that the results of the two numerical techniques and measured results follow the same trend over a wide range of frequencies and allowed extrapolation of the measured reflection and transmission coefficients outwith the measured domain. Simulated and measured results are in good qualitative agreement and the observed deviation tends to increase simultaneously with frequency increase.
The empirical model which allows a single layer of on-chip interconnects to be accurately replaced by homogeneous material slab in electromagnetic simulations containing integrated circuits was presented. The model is applicable to a wide range of interconnect dimensions (metal fill 20 - 60%, aspect ratio 1.4 – 3 and host permittivity 1 – 11.7) and is accurate to better than 2.5% (0.2%) error for reflection (transmission) when illuminated by plane waves with frequency 1 - 10 GHz, incident at up to 30º off the normal. The error does not increase with respect to the change in the grating profile and it can be applied to trapezoidal gratings with sidewall angles up to 5º. Our approach allows the behaviour of on-chip interconnects to be accurately captured in full vector electromagnetic simulations without incurring the significant computational penalties associated with a finely detailed mesh. The experimental data supports the conception that the metal-dielectric grating structure specified for interconnects can be homogenised. The validation was carried out by comparing the numerical results with experimental data.