# A Singularity Free MoM-Type of Formulation Using the Dipole-Moment-Based Approach.

1. INTRODUCTION

Formulating integral equations via the use of Green's functions is a well- established and universally accepted method [1-3] in the context of MoM, and it has been a staple for CEM problems in the past. But MoM requires special treatment at low frequencies where it runs in to difficulties, and it switches to loop-star basis functions to mitigate the problem. Furthermore, MoM needs to deal with the singular and/or hyper-singular behaviors of the Green's functions, and designs special techniques for integrating them when generating the matrix elements. Additionally, both the frequency domain techniques, namely Finite Element Method (FEM) and MoM, experience difficulties when handling multiscale geometries, because the associated matrices for these problems can become ill-conditioned. To mitigate these problems, we introduce a universal MoM-type formulation, which bypasses the use of Green's function to overcome the disadvantages of the conventional frequency domain techniques alluded to above.

2. DIPOLE MOMENT CONCEPT

For the case of a sphere illuminated by a plane wave, the scattered fields can be determined analytically. Consider a PEC (Perfect Electric Conductor) sphere of radius a, which is immersed in free space, and is illuminated by a plane wave [E.sub.x] = [E.sub.o][e.sup.-jkz]. The field scattered by this PEC sphere in the far field region resembles the far field radiated from a combination of an x- directed electric dipole and a y-directed magnetic dipole, whose moments are given by [4]:

[Il.sub.x] = [E.sub.o]4[pi]j/[eta][k.sup.2] [(ka).sup.3] (1a)

[Kl.sub.y] = [E.sub.o] 2[pi]/j[k.sup.2] [(ka).sup.3] (1b)

Along the same lines, we can show that the equivalent dipole moments for a lossless dielectric sphere of radius a, whose relative permittivity and permeability are [[epsilon].sub.r] and [[mu].sub.r] respectively, can be written as:

I[l.sub.x] = [E.sub.o] 4[pi]j/[eta][k.sup.2] [(ka).sup.3] [[epsilon].sub.r] - 1/[[epsilon].sub.r] + 2] (2a)

I[l.sub.y] = [E.sub.o] 2[pi]j/j[k.sup.2] [(ka).sup.3] [[mu].sub.r] - 1/[[mu].sub.r] + 2] (2b)

Equations (2) can be readily modified for a lossy medium by replacing the real valued [[epsilon].sub.r] and [[mu].sub.r], with their complex permittivity e and permeability [bar.[mu]]. It's important to note that the magnetic dipole moment goes to zero for non-magnetic media ([[mu].sub.r] = 1).

Hence, the dipole moment representation of a scatterer generates the same far fields as those scattered by the original object. However, what has not been realized in the past--and what can be proven analytically [5]--is that for a sphere whose radius is electrically small, the dipole moment fields exactly match the original ones scattered by the sphere, all the way up to its surface, and not just in the far field.

3. DM FORMULATION

3.1. Formulation for PEC Objects

When formulating a problem that involves only PEC objects, the first step is to represent the original scatterer by using a collection of PEC spheres. Next these spheres are replaced by their corresponding Dipole Moments (DMs) and a set of them are used to form a suitable set of macro- basis functions. We then evaluate the electric fields generated by these macro basis functions and compute the reactions between them and the testing functions, which are also the same as the basis functions (Galerkin method), to generate the elements of the MoM matrix. The right-hand side of this matrix is obtained by applying the boundary condition on the total tangential E-Field, by testing it with the same functions as those used to generate the matrix elements.

[E.sup.tan.sub.inc] + [E.sup.tan.sub.scat] = 0 (3)

Hence, with the incident E-field polarized along [??], the matrix equation for a thin PEC rod oriented along [??] and modeled by using N macro basis functions, will have the form:

[mathematical expression not reproducible] (4)

In the above Equation (4), I[l.sup.n.sub.z] represents the effective dipole moment of the nth macro basis function--directed along [??], [E.sup.n.sub.z-inc] represents the tangential incident field component at the location of the nth macro basis function and [E.sup.mn.sub.z] represents the scattered field component along [??] on the mth macro basis function by the nth macro basis function. The above matrix Equation (4) is solved for Il's, i.e., the co-efficient of the macro basis functions, and used to compute the induced currents.

3.2. Formulation for Dielectric Objects

The first step in the formulation of the dielectric scattering problem essentially follows along the same line as in the case of PEC objects, in that we again represent the original scatterer as a collection of small-size dielectric spheres. As before, we then go on to replace these spheres by their corresponding DMs, and use them to form a set of macro-basis functions. At this point, we differ from the PEC case and generate the MoM matrix not by imposing a boundary condition but by applying a consistency condition (5) on the tangential E-Field, which reads:

[[epsilon].sub.o] ([[epsilon].sub.r]-1) ([E.sub.inc] + [E.sub.scat]) = F(Il) (5)

where F, called as the consistency factor, is derived by analytically solving the problem of a dielectric sphere with a small radius. This factor, obtained in the manner described above, is given by:

F [approximately equal to] -3j/4[pi][omega][a.sup.3] (6)

4. NUMERICAL RESULTS

4.1. PEC Objects

For the first example we consider a circular PEC loop with a diameter of 600 mm and a thickness of 18.6 mm as shown in Figure 1. The loop is fed with a voltage gap source. Figure 2 compares the frequency variation of the feed current calculated by using DM approach with those obtained using NEC [6].

From Figure 2 we can see that the comparison of the feed current is good all the way down to 50 MHz. It is important to point out that the strength of the DM approach lies in the fact that we can go down further in frequency all the way up to DC, without any special treatment.

Next we consider a PEC sphere, whose diameter is [lambda]/60 with [lambda] at 10 GHz. It is illuminated by a plane wave, incident from [??], and polarized along [??], as shown in the Figure 3. Figure 4 compares the scattered [E.sub.z] field at x = [lambda]/46 calculated by using the DM approach as described in Section 3.1 with those obtained from Mie Series [4] for different frequencies of incident plane wave.

As we can see from Figure 4, the comparison with the Mie series results is good all the way down to very low frequencies, which is one of the major advantages of the DM formulation, as we have pointed out before.

4.2. Dielectric Objects

For the first Dielectric example we consider a square-shaped plate with [[epsilon].sub.r] = 6, which is [lambda]/40 on the side and whose thickness is [lambda]/400. The plate is illuminated by a plane wave traveling along the negative-z direction, as shown in the Figure 5. The backscattered field, calculated by using the DM approach described in Section 3.2, is presented in Figure 6, which also compares these results with the corresponding ones from a commercial MoM package.

We can find that the comparison of the DM results with those from a commercial solver is good. However it should be pointed out, that DM approach has the potential to solve very thin lossy scatterers without any need for modifications in the formulation, while the conventional MoM needs special treatment.

As a second example, we consider plasmonic nano-spheres randomly spread in 3D space, with a diameter of [lambda]/20 at 300THz and [[bar.[epsilon]].sub.r] = -47.5378-- 1.1383.) (corresponding to Gold). It is illuminated by a plane wave, incident from [??], and polarized along [??], as shown in the Figure 7. The scattered fields are calculated by using the DM approach described in Section 3.2 but with the dipole moments Ils calculated using the complex permittivity in the consistency condition (5) and compared with those generated using a commercial FEM code.

From Figure 8 we can see that the scattered fields generated by using the commercial FEM code shows numerical artifacts near the origin, while the fields calculated using the DM approach is more physical.

5. KEY FEATURES OF DM APPROACH

5.1. Quasi-Static Approximation

When the size of problem geometry is very small when compared to the wavelength at the frequency of interest; i.e., r [right arrow] 0[lambda], then the 1/[r.sup.3] term in the scattered field equations will dominate. Hence the expressions for [E.sub.r] and [E.sub.[theta]] in the limit of r [right arrow] 0[lambda] and using the expression for Il from Equation (1a), can be rewritten as:

[mathematical expression not reproducible] (7a)

[mathematical expression not reproducible] (7b)

[E.sub.[phi]] = 0 (7c)

The above equations shows that the electric fields are real and time- independent. These expressions resemble the fields of a static charge dipole. Since this quasi-static approach produces real and timeindependent fields, we can use this to hybridize DM approach with Finite Difference Time Domain (FDTD) to solve a variety of multi-scale problems [7, 8].

The quasi-static approximation can be used for problems for which near-field calculations are of interest, since it is predominantly dictated by the 1/[r.sup.3] term; also, this approach is computationally less expensive.

5.2. Order of Singularity

The boundary condition or the consistency condition is always applied on the surface of the geometry. So for a problem geometry comprising of only PEC objects, the field values at the geometry's surface r = a can be calculated from Equation (7) as follows:

[E.sub.r] = 2[E.sub.o] cos [theta] (8a)

[E.sub.[theta]] = [E.sub.o] sin [theta] (8b)

[E.sub.[phi]] = 0 (8c)

Hence, in the DM approach the fields are always bounded and the order of the singularity is zero. The same can be proven to be true for dielectric and inhomogeneous objects [5].

6. PERFORMANCE ENHANCEMENT

Though the DM approach described in previous sections is accurate and captures all the physics, it is not the most efficient from numerical point of view. This is because the number of spheres used to represent a three-dimensional object can grow very rapidly if the diameter of the sphere is small, as is often the case. For instance, for a thin-wire scatterer, the diameter of the spheres used to represent it is the same as that of the wire. Hence, for the example shown in Figure 5, the number of constituent spheres needed to form the plate can be quite large, even when the length of the plate is relatively small in comparison to the wavelength. In the following sections, we introduce techniques to enhance the performance of DM approach and to adapt this approach to solve a number of representative electromagnetic scattering problems. 7

7. HIGHER ORDER MACRO-BASIS FUNCTION

Our strategy for reducing the number of unknowns significantly and to make it comparable to that needed in the conventional MoM formulation, is to use macro-basis functions (MBF). These basis functions belong to a level higher than that of the dipole moments used to model the geometry in the initial step of the DM approach. The low-level dipole moments associated with such macro-basis functions are represented by a single unknown, with the variation of the dipole moments following the shape of the associated macro-basis function. In order to further improve the performance by reducing the number of unknowns, we introduce higher-order basis functions (HBFs). Towards this end, we use a set of macro-basis functions and form a set of suitable higher-order macro-basis functions and represent them by using a single unknown. The coefficients of the macro-basis function follow the shape of the higher-order macro-basis function. The MBFs can be categorized as sub-domain basis functions, as opposed to entire domain basis functions. Some of the commonly used low-level basis functions used in the concept of MoM formulation are triangles, pulses and Rao-Wilton-Glisson (RWGs) [9], but we will use them as MBFs for the DM approach.

7.1. Numerical Results

As a first example, let us consider a PEC sphere with a diameter of [lambda]/60 at 10 GHz. It is illuminated by a plane wave, incident from [??], and polarized along [??], as shown in Figure 9. Figure 9 compares the scattered [E.sub.z] field along the negative-x axis, computed by using DM approach; DM approach using MBFs; DM approach with HBFs; and, with those obtained from Mie Series [4]. For this problem, triangular basis functions were used, both as MBFs and HBFs, in the context of the DM approach.

Table 1 compares the number of unknowns required to solve the problem using (i) DM approach alone; (ii) with MBFs and (iii) with HBFs. Table 1 shows that the use of HBFs significantly reduces the number of unknowns without sacrificing the accuracy, as may be seen by referring to Figure 9.

For the next example, we consider a square-shaped dielectric plate with [[epsilon].sub.r] = 6, which is [lambda]/40 on the side and whose thickness is [lambda]/400 (Figure 5). The plate is illuminated by a plane wave traveling along the negative-z direction. Figure 10 compares the backscattered field calculated by using: (i) the DM approach only; (ii) the DM approach with MBFs; and (iii) a commercial MoM solver. Table 2 compares the required number of unknowns in each of these cases. For this problem the rooftop basis function was used as MBFs in the context of DM approach.

The results show, once again, that the comparison of the fields is good and the use of MBFs greatly reduces the number of unknowns without compromising the accuracy of the results. As we can see from the previous examples, it is relatively easy to choose these macro-basis functions. In order to better capture the current behavior near the corners and edges, we can increase the number of MBFs near the corners and edges to refine the level of discretization. Also, we can reduce the number of HBFs as we go down in frequency, since the current distribution will not vary rapidly at lower frequencies. 8

8. CHARACTERISTIC BASIS FUNCTIONS

Characteristic basis functions (CBFs) [10] are higher-level entire domain macro- basis functions, and their use helps reduce the size of the matrix rather significantly. The CBFs are tailored for the geometry at hand and often just one or two CBFs are sufficient to solve the problem if the object is relatively small in size. The use of CBFs enables one to solve electrically large problems, because the size of the reduced matrix is often orders of magnitude smaller than one required in the original MoM formulation to achieve the same level of accuracy. The CBFs are physics-based and lead to well- conditioned matrices [11], because their redundancy is removed via the use of Singular Value Decomposition (SVD).

It is a common practice to generate the CBFs by solving for the current distributions by using a number of independent excitations, which the angle of incidence and polarization of the illuminating wave is varied. Next, a matrix is generated by using the resulting current distributions, as its columns and a SVD of this matrix is performed. The threshold for the singular values is typically chosen to be 1% of the highest singular value. Finally, we use the vectors corresponding to these singular values to construct the CBFs [12].

Consider a PEC rod, whose length and diameter are [lambda]/10 and [lambda]/400; respectively, as shown in the Figure 11. The CBFs for this problem were calculated as explained above and are plotted in Figure 12. As we can see from this figure, only two out of the twelve originally generated solutions survived the SVD thresholding.

To test the method just described, we consider the case of the above PEC rod illuminated by a plane wave, incident from the broadside direction. Figure 13 compares the backscattered field calculated by using: (i) the DM approach only; (ii) DM approach using CBFs; and, (iii) commercial MoM solver.

From Figure 13 and Table 3 we see that both the DM results are in agreement with those obtained by using a commercial MoM solver. However, the DM/CBF approach reduces the number of unknowns by a large factor, namely 50, in this example.

9. FAST MATRIX GENERATION

Another way to improve the performance of DM approach is to adapt the Fast Matrix Generation (FMG) technique [13], proposed for the generation of matrix elements in the context of the conventional MoM formulation. To adapt this technique for the DM approach, we consider only the 1/r term in the field expression while calculating the interaction between basis functions separated by a distance greater than or equal to [lambda]/10 and while we consider the 1/[r.sup.2] and 1/[r.sup.3] terms when the distance is less than [lambda]/10 separation.

Consider a square-shaped dielectric plate, shown in Figure 14. The backscattered field, calculated by using the DM approach with and without the use of FMG algorithm, is presented in Figure 14, which also compares these results with the corresponding ones from a commercial MoM package. Table 4 compares the CPU time required by these two different approaches.

Figure 14 shows a good comparison of the electric field calculated by using the DM approach, with and without the FMG, and also with those from a commercial solver. Even though Table 4 shows only a slight time advantage of FMG over the DM, it has been shown elsewhere that this advantage grows rapidly as we increase the electrical length of the geometry. Here we are simply demonstrating that we do not sacrifice the accuracy when we use FMG along with the DM approach.

10. CLOSED-FORM FIELD EXPRESSIONS

To calculate the elements of the interaction matrix, we need to sample the macro-basis function for the purpose of approximating the integration with a numerical summation. Let us consider a triangular current distribution which extends over a length of [lambda]/10, as shown in Figure 15. The fields radiated by this triangular current element is calculated along a parallel line, at an offset of [lambda]/50, using different numbers of samples of the current element. Figure 16 compares the [E.sub.z]- field variation, for different numbers of samples, with the closed-form field expression [14] for the same current distribution.

We note from Figure 16 that we need at least 200 samples to achieve a good match between the closed-form result for the integration and its approximation via numerical summation, when the offset distance is of [lambda]/50. However, as shown in Figure 17, when we move further, say to a distance of 500A, we can achieve a good match between the direct integration and numerical approximation with just 2 samples. This implies that to compute the matrix entries for each macro basis function we may need up to 200 samples, depending on the distance where we apply the boundary condition. It would be even more computationally expensive if we use a fine mesh, or if the object is electrically large. Hence, it would be useful to find a closed-form expression for the field, generated by the current distribution, by representing it in a suitable form as shown below.

A typical current distribution on a wire is piecewise sinusoidal in nature, as represented in Equation (9) for a wire length of [H.sub.1] + [H.sub.2].

[mathematical expression not reproducible] (9)

where [I.sub.m] is the maximum amplitude of the current.

The closed-form expression for the field from this type of current distribution along a bent wire (see Figure 18) can be derived, with separate expressions for the top half and bottom half of the current distribution, as follows:

Top Half (0 [right arrow] [H.,sub.1])

[mathematical expression not reproducible] (10a)

[mathematical expression not reproducible] (10b)

Bottom Half (-[H.sub.2] [right arrow] 0)

[mathematical expression not reproducible] (11a)

[mathematical expression not reproducible] (11b)

where the directions [mathematical expression not reproducible] are unit vectors along the axis of the wire, while [mathematical expression not reproducible] are perpendicular to its axis, as shown in Figure 18. From the above equations we can see that when [[upsilon].sub.1] or [[upsilon].sub.2] is 0; i.e., when the observation point is either along [mathematical expression not reproducible] the electric field [mathematical expression not reproducible] becomes singular. In order to calculate the correct field values for these cases, we model the wire geometry with a sinusoidal current distribution as shown in Figure 19, and use the DM approach with 200 samples. The [E.sub.y] fields are calculated along a observation line parallel to [??], by using the expressions given in Equations (10) and (11), as well as by using the DM approach. Figure 20 plots the fields calculated by using these two approaches and they all are seen to agree well with each other. We can also see from this figure that the value of the field monotonically increases from zero as we move away from the wire axis. Hence, we set the field value to be zero along the wire axis whenever we are interested in the field value at a point located on the axis.

10.1. Numerical Results Using Triangular Basis Functions

It is important to note the fact that the sinusoidal basis function closely resembles a triangular basis function (TBF), provided the support of the basis is less than or equal to [lambda]/10. Hence we can use the closed form expressions given in Equations (10) and (11) to calculate the field values for the most commonly used TBFs, whose supports are less than or equal to [lambda]/10.

For the first example we consider a PEC rod, whose length and diameter are [lambda] and [lambda]/500, respectively. It is illuminated by a plane wave, incident from broadside, as shown in Figure 21. Figure 22 compares the peak amplitude of the induced current in the PEC rod calculated by using the DM approach and closed form expressions, with those obtained from a commercial MoM package for different frequencies of the incident plane wave, ranging from 1 MHz to 10 GHz.

From Figure 22 we can see that the results from the DM approach using closed- form expressions compare well with those from the commercial MoM solver. However, we have found that the commercial MoM solver failed when we go below 10 Hz, while the DM approach is able to handle the problem without any special treatments or modifications. The number of TBFs used to model the wire varies with frequency, starting with 9 TBFs at 10 GHz and progressively decreasing to 1 at 1MHz (or below), as listed in Table 5. We expect this to be the case since the variation in the current distribution on the wire varies less rapidly as we decrease the frequency. Incidentally, although it makes little difference in the accuracy level of the solution whether we use the same number of TBFs over the entire frequency range or decrease their number progressively, we find that the condition number of the interaction matrix improves when we use a variable number of TBFs with frequency, as may be seen from Figure 23.

In order to study the improvement in performance with the use of closed-form expressions, we consider a circular PEC loop with a diameter of [lambda]/23.3 and a thickness of [lambda]/100 at 1GHz. It is illuminated by a plane wave, incident from the [??] direction, and polarized along [??], as shown in Figure 24. Figure 25 compares the variation of the induced current calculated by using the DM approach with that derived by using the DM approach and closed-form expressions. Table 6 compares the times required for these two approaches and we clearly see that the use of closed-form expression speeds up the process by a factor of 3.3 for this problem without compromising the accuracy, as may be seen in Figure 25.

10.2. Numerical Results Obtained by Using Rooftop Basis Functions

One of the most commonly-used basis functions for representing current distributions on surfaces is the Rooftop. It is a two-dimensional basis function comprising a pulse basis function along one direction, and a triangular basis function along the other, as shown in Figure 26. The current density distribution associated with a rooftop is given by:

J(x, y) = [I.sub.m] sin [k(H - [absolute value of u])] rect (x/2H) rect, (y/w) (12)

In order to model this rooftop using the field expressions previously derived in this section, we represent the rooftop basis function with a number of TBFs that have the same maximum amplitude. We carried out a wide range of numerical experiments with the number of TBFs, and have found that we need to represent a rooftop with 7 TBFs to get accurate results. As an example, let us consider a square-shaped PEC plate, which is [lambda]/2 on the side, and whose thickness is [lambda]/25. The plate is illuminated by a plane wave traveling along the negative-z direction, as shown in Figure 27. The backscattered field, calculated by using the DM approach and closed-form expressions, is presented in Figure 28, which also compares these results with the corresponding ones from a commercial MoM package. To calculate the field directly at the center of the plate by using the DM approach, we take the average of the field over the footprint of the rooftop basis function at the center of the plate. We do this in the DM approach in order to improve its accuracy, instead of using point matching for this case, as we do for other locations on the plate.

Table 7 compares the time required by the DM approach with and without the use of closed-form expressions. It is evident, from Figure 28, that the use of closed-form expression speeds up the process without compromising the accuracy, even when the observation point is close to the surface of the plate. To further improve the performance in terms of the CPU time without compromising the accuracy, we can use 7 rooftops to compute the self term, 5 rooftops to calculate the interaction between the rooftops located within a distance of [lambda]/10; 3 rooftops for the calculation of interaction between the rooftops when this separation distance is greater than [lambda]/10 but less than [lambda]/5; and, a single TBF for separation distances greater than [lambda]/5.

11. SOME EMBELLISHMENTS TO THE BASIC DM APPROACH

Lumped loads are often used either to match the impedance of an antenna or to shift its resonance. The resonance behavior achieved by using lumped loads is often sharp, and requires a fine frequency sampling to capture this resonance behavior. From our experience, many of the commercial solvers fail to capture these resonances and it would be useful to incorporate lumped loads in the DM approach to see if it performs better than the existing codes. Incorporating lumped loads in the DM approach is relatively simple, and is achieved by the addition of the lumped load impedance to certain matrix elements depending upon the spatial locations of the load.

As an example, we consider a power coil comprising of 4 loops, shown in Figure 29, which is modeled by using 192 TBFs. Loop-1 is fed by using a voltage gap source. The frequency range of interest is 6 MHz to 11MHz, and we expect a series resonance to occur around 8 MHz, between the lumped capacitance and the inductance of the loop. Figure 30 compares the peak current in the output coil, i.e., the loop-4, obtained by using the DM approach, and compares it with the one generated by using the NEC code.

Figure 30 shows a good comparison between the amplitude variations; however, there is a mismatch between the amplitude at the resonant frequency, which is attributable to the sharpness of the resonance, and the difference in the methodologies used in the DM approach and the NEC.

11.2. Input Impedance of Nano Antennas

As we move towards the miniaturization of electronic devices, the sizes of the electronic components also follow suit. One of the important components of interest is the small antenna, the evaluation of whose input impedances becomes more challenging as the cross section of the small antenna becomes comparable to its length. A whole host of techniques have been proposed in the literature for computing the input impedance of antennas. Harrington has presented an expression given in (13), for evaluating the input impedance of an antenna by using the current induced on the same.

[mathematical expression not reproducible] (13)

Since expression (13) is variational, it supposedly generates a result for the impedance that is second-order accurate even when the induced current [J.sub.s] inserted in (13) has only first-order accuracy. However, we will show that by using the DM approach we can accurately calculate the input impedance of small antennas. Consider the vertical monopole antenna, shown in Figure 31, which has a length of [lambda]/40 and a thickness of [lambda]/500 at 10 GHz. Figure 32 shows the current variation along the antenna calculated by using the DM approach. Table 8 below summarizes the input impedance results computed by using the DM approach; with Harrington's formula in (13); and with a simplified transmission line model.

Table 8 shows that using the conventional truncated sinusoidal type of representation for the induced current in (13) yields results that deviate considerably from those computed accurately by employing the DM approach.

11.3. Irregular Geometries

Handling irregular geometries can be challenging, since we need to use different mesh sizes and basis functions to accurately model different parts of the geometry. In the conventional MoM as well as FEM, this variation in the mesh sizes can lead to a poorly-conditioned matrix. To handle such geometries using the DM approach, we first calculate the elements of the interaction matrix using the closest possible regular geometry, and then modify the corresponding elements with the ratio of areas of the footprints of the basis functions in the regular and the actual geometries.

Consider a PEC plate with a thickness of [lambda]/25, shown in Figure 33. The closest regular geometry is a [lambda]/2 square plate with the same thickness, i.e., [lambda]/25. The interaction matrix is generated by using the regular square PEC plate using the DM approach as described in the Section 10.2. This interaction matrix is then modified, to handle the actual geometry, by multiplying the difference in the foot prints of the rooftop basis functions before calculating the currents. Figure 34 compares the scattered electric field pattern, at a radial distance of 2[lambda], with [phi] = 45[degrees], calculated by using the DM approach with those obtained from the commercial MoM solver at 5 GHz.

From Figure 34 we find that the results from DM approach show good agreement with those from the commercial MoM solver. However, it was found that when the geometry becomes thinner and the number of irregularities increase, the commercial MoM solver is not able to handle the problem, while the DM approach was able to do so with relative ease.

11.4. Curved Surfaces

Another object of interest is a faceted surface, which is difficult to handle with the conventional solvers when the geometry has sharp edges. Consider the faceted PEC surface shown in Figure 35. It has a height of [lambda]/4, a width of [lambda]/20 and a thickness of [lambda]/25 at 5 GHz. Figure 36 compares the backscattered field calculated by using the DM approach (modeled with 40 rooftops), as described in the Section 10.2, with that calculated from a commercial MoM solver.

From Figure 36 we see good comparison between the fields calculated by using the DM approach with those obtained from a commercial solver. Whenever we deal with either dielectric objects with curved surfaces, or those with thin curved geometries the commercial solver generates an ill-conditioned matrix, whose solution is questionable in terms of accuracy; however, the DM approach exhibits no such behavior.

11.5. Geometries with Apertures

DM approach can handle geometries with apertures or slits with relative ease without requiring any modification to the approach. Consider a rectangular PEC plate, split at the center, with a gap of [lambda]/40 and a thickness of [lambda]/25 at 10 GHz, as shown in Figure 37. Figure 38 compares the backscattered field computed by using the DM approach as described in the Section 10.2, with that calculated by using a commercial MoM solver.

Figure 38 shows good agreement between the two backscattered fields even as we approach the surface of the plate. If we decrease the split gap size to [lambda]/80, the associated matrix in the commercial solver becomes ill-conditioned, while the associated matrix in the DM approach remains well conditioned regardless of the gap size.

For the next example, we consider a PEC square plate of side length of [lambda]/2 with a square aperture of size [lambda]/20, at a frequency of 10 GHz, as shown in Figure 39. In order to model the aperture, the first step is to calculate the E-fields over the aperture. Hence we model a smaller plate of size [lambda]/4 with the same aperture size using the DM approach without MBFs and calculate the dominant field component in the aperture, namely the [E.sub.y] fields in this case. Using the [E.sub.y] field so obtained, we compute the magnetic field current [mathematical expression not reproducible] over the aperture. Once we know the equivalent magnetic current in the aperture, we can back it with a PEC surface by invoking the Huygens' principle. The total field will be a summation of two sets of fields: (i) scattered by the PEC square plate with the aperture closed, but in the presence of [mathematical expression not reproducible], solved by using rooftop basis function in the context of DM, as described in Section 10.2; (ii) fields radiated by [mathematical expression not reproducible], either when it placed over a ground plane, or by [mathematical expression not reproducible] in free space, if we make the assumption that the plane is infinitely large. However this approximation is justified since the fields radiated by [mathematical expression not reproducible] over the surface z = 0 is concentrated only near the aperture and they rapidly decrease as we move away from the aperture, as shown in Figures 40 and 41.

Figure 42 compares the amplitude of the backscattered field [E.sub.y] calculated by using DM approach with those computed by using a commercial MoM solver.

11.6. Microstrip-Based Structures

As the size of the semiconductor devices goes down, the thickness of the microstrip substrate also becomes smaller. With increasing integration of electronic packages sharing the same substrate, the mutual coupling becomes critical, and often requires an accurate modeling. The thinness of the substrate, along with the finer widths of the signal traces challenges the existing CEM algorithms when used to calculate the electromagnetic response characteristics of such packages. To demonstrate that the DM approach can handle geometries with fine structures, we consider the example of a microstrip transmission line shown in Figure 43, whose length is 2[lambda] and, which has free space as its substrate. The transmission line is modeled by using the DM approach with a voltage source exciting the line at one end (x = -[lambda]), under the infinite ground plane approximation, which is typical. We use rooftops, as described in Section 10.2, to model the current densities [J.sub.x] and [J.sub.y] with [lambda]/10 x [lambda]/100 and [lambda]/20 x [lambda]/20, respectively, at a frequency of 10 GHz. Figure 44 plots the variation of [I.sub.x] along the trace calculated by using the DM approach and shows the standing wave pattern as excepted.

For the next example, we consider a similar microstrip geometry, as shown in Figure 45. The line is illuminated by a plane wave traveling along the negative-z direction and polarized along [??]. Figure 46 compares the scattered field along [??], calculated using the DM approach, as described in the previous example, with those calculated by using a commercial MoM solver.

Figure 46 shows a good agreement between the results generated by the DM approach and those derived from a commercial MoM solver. We point out that the advantage of using the DM approach lies in the fact that it continues to work well when we make the substrate thinner, add losses to it, and use even finer signal traces. The commercial solvers, on the other hand, break down under these circumstances.

12. OBSERVATIONS AND CONCLUSIONS

In this paper we have presented a new physics-based approach for formulating MoM problems that is based on the use of dipole moments (DMs)--as opposed to the conventional Green's functions. The absence of the Green's function, as well as those of the vector and scalar potentials, helps to eliminate two of the key sources of difficulties in the conventional MoM formulation, namely the singularity and low-frequency problems. Specifically, we have shown that there are no singularities that we need to be concerned with in the DM formulation; hence, this obviates the need for special techniques for integrating these singularities.

Yet another salutary feature of the DM approach is its ability to handle thin and lossy structures, whether they are metallic, dielectric-type, or even combinations thereof. We have found that the DM formulation can handle these types of objects with ease, without running into ill-conditioning problems, even for very thin wire-like or surface-type structures, which lead to ill- conditioned MoM matrices when these problems are formulated in the conventional manner.

The technique is valid over the entire frequency range, from low to high, and it does not require the use of loop-star or other special types of basis functions in order to mitigate the low frequency problem. The DM formulation is universal and can be used for both PEC and dielectric objects, and it requires only a relatively minor change in the formulation when we go from PEC to dielectric scatterers. The approach is also well suited for hybridization with Finite methods, such as the FEM and the FDTD, and such an embellishment renders it suitable for handling multi-scale problems conveniently and efficiently.

We have shown that the use of higher-order basis functions significantly reduces the number of unknowns, without compromising the accuracy and combines the DM with the CBFM technique helps reduce this number even further. The use of closed-form expressions for the interaction matrix elements speeds up the process of matrix generation, regardless of the problem size. For electrically large problems, employing FMG helps to speed up the interaction matrix generation considerably.

We have shown how we can incorporate lumped loads in the DM approach and that it is able to capture sharp resonances even at low frequencies, where the commercial solvers become inaccurate or break down. The DM approach is able to accurately calculate the input impedance of small antennas; fields from irregular geometries; from faceted surfaces; from geometries with slot and slit; and, is able to model microstrip line type of geometries with fine features. In all of the above examples we have solved for, the matrices associated with the DM approach remained well-conditioned throughout the entire frequency range of interest, without the use of special treatments. However, this was not the case with Commercial MoM and FEM solvers, even after special treatments were incorporated in these solvers.

REFERENCES

[1.] Peterson, A. F., S. L. Ray, and R. Mittra, Computational Methods for Electromagnetics, IEEE Press, New Jersey, 1997.

[2.] Harrington, R. F., Field Computation by Moment Methods, The Macmillan Company, New York, 1968.

[3.] Mittra, R., Computer Techniques for Electromagnetics, Hemisphere Publishing Corporation, New York, 1987.

[4.] Harrington, R., Time-harmonic Electromagnetic Fields, IEEE Press, New Jersey, 2001.

[5.] Bringuier, J., "Multi-scale techniques in computational electromagnetics," Ph.D. Dissertation, The Pennsylvania State University, 2010.

[6.] "Numerical electromagnics code," [Online], available: http://www.nec2.org/.

[7.] Panayappan, K., "Novel frequency domain techniques and advances in finite difference time domain method (FDTD) for efficient solution of multiscale electromagnetic problems," Ph.D. Dissertation, The Pennsylvania State University, 2013.

[8.] Mittra, R., Computational Electromagnetics: Recent Advances and Engineering Applications, Springer-Verlag, New York, USA, 2013.

[9.] Rao, S., D. Wilton, and A. Glisson, "Electromagnetic scattering by surfaces of arbitrary shape," IEEE Transactions on Antennas and Propagation, Vol. 30, No. 3, 409-418, 1982.

[10.] Lucente, E., A. Monorchio, and R. Mittra, "An iteration-free MoM approach based on excitation independent characteristic basis functions for solving large multiscale electromagnetic scattering problems," IEEE Transactions on Antennas and Propagation, Vol. 56, 999-1007, April 2008.

[11.] Mittra, R. and K. Du, "Characteristic basis function method for iteration- free solution of large method of moments problems," Progress In Electromagnetics Research B, Vol. 6, 307-336, 2008.

[12.] Mehta, N., "Numerical analysis of frequency selective surfaces," Master's Thesis, The Pennsylvania State University, 2010.

[13.] Kwon, S. J. and R. Mittra, "Impedance matrix generation by using fast matrix generation (FMG) technique," Microwave and Optical Technology Letters, Vol. 51, 204-213, January 2009.

[14.] Jordan, E. and K. Balmain, Electromagnetic Waves and Radiating Systems, Prentice-Hall Series in Electrical Engineering, Prentice-Hall, 1968.

Kadappan Panayappan (1) and Raj Mittra (1, 2) *

(Invited Paper)

Received 14 January 2015, Accepted 30 March 2015, Scheduled 7 April 2015

* Corresponding author: Raj Mittra (rajmittra@ieee.org).

(1) Electromagnetic Communication Lab, Pennsylvania State University, USA. 2 University of Central Florida, USA.

Caption: Figure 1. A PEC circular loop.

Caption: Figure 2. Frequency variation of the feed current for a PEC loop shown in Figure 1.

Caption: Figure 3. A PEC sphere.

Caption: Figure 4. Amplitude comparison of backscattered electric field [E.sub.z] from the PEC sphere in Figure 9.

Caption: Figure 5. A dielectric plate.

Caption: Figure 6. Amplitude comparison of backscattered electric field [E.sub.y] from the dielectric plate in Figure 5.

Caption: Figure 7. Plasmonic nano-spheres.

Caption: Figure 8. Amplitude comparison of scattered electric field [E.sub.y] from the plasmonic spheres in Figure 7.

Caption: Figure 9. Amplitude comparison of scattered electric field [E.sub.z] from the PEC sphere in figure.

Caption: Figure 10. Amplitude comparison of backscattered electric field [E.sub.y] from the dielectric plate in Figure 5.

Caption: Figure 11. A PEC rod.

Caption: Figure 12. Calculated CBFs for the PEC rod in Figure 11.

Caption: Figure 13. Amplitude comparison of backscattered electric field [E.sub.y] from the PEC rod in Figure 11.

Caption: Figure 14. Amplitude comparison of backscattered electric field [E.sub.y] from the dielectric plate in Figure 5.

Caption: Figure 15. Triangular current distribution placed along Z axis from--[lambda]/20 to [lambda]/20.

Caption: Figure 16. Amplitude variation of [E.sub.z] at a 2 = [lambda]/50 radiated by the current distribution shown in Figure 15.

Caption: Figure 17. Amplitude variation of [E.sub.z] at a z = 500[lambda] radiated by the current distribution shown in Figure 15.

Caption: Figure 18. A bent wire.

Caption: Figure 19. A wire geometry.

Caption: Figure 20. Amplitude variation of [E.sub.y] along the observation line in A for the wire geometry in Figure 19.

Caption: Figure 21. A PEC rod.

Caption: Figure 22. Frequency variation of peak current for the PEC rod in Figure 21.

Caption: Figure 23. Frequency variation of condition number for the PEC rod in Figure 21.

Caption: Figure 24. A PEC circular loop.

Caption: Figure 25. Amplitude variation of the induced current for the PEC loop in Figure 24.

Caption: Figure 26. A Rooftop basis function.

Caption: Figure 27. A PEC plate.

Caption: Figure 28. Amplitude comparison of backscattered electric field [E.sub.y] from the PEC plate in Figure 27.

Caption: Figure 29. A geometry of a power coil.

Caption: Figure 30. Frequency variation of the peak output current for a power coil shown in Figure 29.

Caption: Figure 31. A vertical monopole antenna.

Caption: Figure 32. Amplitude variation of the current for a monopole shown in Figure 31.

Caption: Figure 33. A PEC plate with a staircase corner.

Caption: Figure 34. Amplitude comparison of scattered electric field [E.sub.[theta]] at [phi] = 45[degrees] from the PEC geometry in Figure 33.

Caption: Figure 35. A faceted PEC surface.

Caption: Figure 36. Amplitude comparison of backscattered electric field [E.sub.z] from the faceted PEC surface in Figure 35.

Caption: Figure 37. A PEC plate with a split.

Caption: Figure 38. Amplitude comparison of backscattered electric field [E.sub.y] from the PEC surface in Figure 37.

Caption: Figure 39. A PEC plate with a square slot.

Caption: Figure 40. Equivalent magnetic current.

Caption: Figure 41. Amplitude variation of the radiated electric field [E.sub.y] from the equivalent magnetic current in Figure 40.

Caption: Figure 42. Amplitude comparison of backscattered electric field [E.sub.y] from the PEC surface in Figure 39.

Caption: Figure 43. A microstrip geometry.

Caption: Figure 44. Variation of Ix along the trace of a microstrip line in Figure 43.

Caption: Figure 45. A microstrip geometry.

Caption: Figure 46. Amplitude comparison of backscattered field [E.sub.x] from a microstrip line in Figure 45.
```Table 1. Comparison of unknowns required for DM approach
using MBF and HBFs for the PEC sphere shown in Figure 9.

Method                 No. of Unknowns   Time in Seconds

DM Approach Only            2322               2.6
DM Approach with MBF         86               0.06

DM Approach with HBF         43               0.04

Table 2. Comparison of unknowns required for DM approach
using MBFs for the dielectric plate shown in Figure 5.

Method                  No. of Unknowns   Time in Seconds

DM Approach Only             8112              81.31
DM Approach with MBFs         50               0.04

Table 3. Comparison of unknowns required
for DM approach using CBFs for the PEC
rod shown in Figure 11.

Method                  No. of Unknowns

DM Approach Only              100
DM Approach with CBFs          2

Table 4. Comparison of simulation time for
DM approach using FMG for the dielectric
plate shown in Figure 5.

Method                 Simulation Time

DM Approach Only        19.22 minutes
DM Approach with FMG    18.98 minutes

Table 5. Number of TBFs used in different
frequency ranges for the PEC rod shown
in Figure 21.

Frequency Range   No. of TBFs Used

1 MHz-1 GHz              1
1 GHz-3 GHz              3
3 GHz-5 GHz              5
5 GHz-7 GHz              7
7 GHz-10 GHz             9

Table 6. Comparison of simulation times using DM
approach with and without closed-form expressions
for the PEC loop shown in Figure 24.

Method                            Simulation Time

DM Approach without Closed-form    0.53 seconds
DM Approach with Closed-form       0.16 seconds

Table 7. Comparison of simulation times using the DM
approach, with and without closed-form
expressions, for the PEC plate shown in Figure 27.

Method                            Simulation Time

DM Approach without Closed-form     69 seconds
DM Approach with Closed-form        57 seconds

Table 8. Comparison of input impedance calculated by using DM
approach, Harrington's approach (13) and a simplified transmission
model for the monopole shown in Figure 31.

Method                          Input Reactance

DM Approach                     - 1.1916e + 003j
Variational Formula (13)           - 2.5216,)
Simplified Transmission Model   - 1.4080e + 003j
```