Cagniard-DeHoop Technique-Based Computation of Retarded Partial Coefficients: The Coplanar Case

Received June 22, 2020, accepted August 5, 2020, date of publication August 13, 2020, date of current version August 24, 2020. Digital Object Identifier 10.1109/ACCESS.2020.3016316 Cagniard-DeHoop Technique-Based Computation of Retarded Partial Coefficients: The Coplanar Case MARTIN ŠTUMPF 1, (Senior Member, IEEE), GIULIO ANTONINI 2, (Senior Member, IEEE), AND ALBERT RUEHLI3, (Life Fellow, IEEE) 1Department of Radioelectronics, Faculty of Electrical Engineering and Communication, Brno University of Technology, 616 00 Brno, Czech Republic 2UAq EMC Laboratory, University of L’Aquila, 67 100 L’Aquila, Italy 3MST EMC Laboratory, Missouri University of Science and Technology, Rolla, MO 65409, USA Corresponding author: Martin Štumpf ([email protected]) This work was supported in part by the Czech Science Foundation under Grant 20-01090S. ABSTRACT Efficient computation of partial elements plays a key role in the Partial Element Equivalent Circuit (PEEC) method. A novel analytical method for computing retarded partial coefficients based on the Cagniard-DeHoop (CdH) technique is proposed. The methodology is first theoretically developed and then illustrated on the computation of a surface retarded partial coefficient pertaining to two coplanar rectangular surface elements. An efficient way for incorporating loss mechanisms in the time domain (TD) via the Schouten-Van der Pol theorem is proposed. Illustrative numerical examples demonstrating the validity of the introduced solution are given. INDEX TERMS Time-domain analysis, partial element equivalent circuit, Cagniard-DeHoop method. I. INTRODUCTION The PEEC method is an integral-equation technique (e.g. [1], [2]) that is capable of translating an EM scattering problem to an equivalent-circuit representation [3]. Electric and magnetic field couplings are represented separately with partial potential coefficients (P) [4], [5], partial inductances (Lp) [6], [7], as well as resistances. The PEEC method is well suited for both high and low frequency problems down to dc [8], [9]. If a volume-equivalence formulation is assumed [10], the free-space full-wave Green’s function (e−ikR/Rwith k=ω/c0, with c0the speed of light in the free space) is to be considered in the computation of partial elements. If the surface-equivalence principle is adopted, the full-wave Green’s function within homogeneous lossy dielectric/magnetic materials is to be used along with its curl [11] in the computation of interaction integrals describing the electric and magnetic field coupling. Potential integrals over triangular domains in the static limit (k→0) have recently been analyzed in [12]. When the volume-equivalence is used to derive the PEEC formulation, the exponential term in the full-wave Green’s The associate editor coordinating the review of this manuscript and approving it for publication was Su Yan . function corresponds to a retardation. In this case, the exponential term was moved outside of the integrals by using the distance between the center of the source and the observation point [13]. Such an approximation is possible for moderate frequency problems [14]–[16] and allows to avoid to integrate the frequency dependent matrix elements at each frequency. However, this does not lead to highly accurate results for very high frequencies. A better approximation can be obtained by resorting to a Taylor series expansion of the Green’s function [17]. Although more accurate to catch the frequency-dependence of the Green’s function and more effective than the brute-force numerical integrationbased approach, the resulting partial elements are easily computed in the frequency domain but their TD counterpart is not straightforward and would require auxiliary equations to implement higher order terms of the Taylor’s expansion. In this paper, we shall provide a fundamentally novel analytical approach for computing retarded coefficients of a PEEC model based on the CdH technique [18]. The methodology is illustrated on the computation of a retarded partial coefficient pertaining to two coplanar rectangular surface elements. It demonstrated that in such a case the surface retarded potential coefficients can be derived analytically in the TD in terms of elementary functions only. Moreover, a general VOLUME 8, 2020 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ 148989 M. Štumpf et al.: CdH Technique-Based Computation of Retarded Partial Coefficients: The Coplanar Case FIGURE 1. Two interacting rectangular coplanar surface elements. procedure based on the Schouten-Van der Pol theorem [19, p. 1056] enabling the inclusion of general loss mechanisms is proposed. II. PROBLEM FORMULATION A PEEC model is represented through a set of partial elements, the value of which is found upon evaluating spatial integrals over the domains of interacting discretization elements. For the sake of simplicity, we shall limit ourselves in this work to the interaction of two coplanar rectangular surface elements (see Fig. 1). Other configurations such as rectangular surfaces on distinct parallel or orthogonal planes or parallelepipeds will be considered in future reports. To localize the position in the problem configuration, we employ coordinates {x,y,z}with respect to an orthogonal Cartesian reference frame with the origin Oand the standard basis {ix,iy,iz}. Consequently, the position vector is r= xix+yiy+ziz. The time coordinate is denoted by t. Finally, ∗tdenotes the time-convolution operator, H(x) represents the Heaviside unit-step function and δ(x) denotes the 1-D Diracdelta distribution. With reference to [3, eq. (D.1)], we shall study a surface retarded partial coefficient expressed through a 2-D integral ˆ Pmn(s)=1 SmSnZr∈Am dSZr0∈Anˆg(r−r0,s)dS0(1) where Am= {−1x/2<x−xm< 1x/2,−1y/2< y−ym< 1y/2,z=0}and An= {−1x/2<x−xn< 1x/2,−1y/2<y−yn< 1y/2,z=0}, where 1x>0 and 1y>0 denote the spatial discretization steps in the xand y-direction, respectively. Furthermore, sis the (real-valued and positive) Laplace-transform parameter, Sm,n=1x1yis the surface area of domain Am,n, and ˆg(r−r0,s)=exp[−s|r−r0|/c] 4π|r−r0|(2) is the free-space Green’s function pertaining to a homogeneous, isotropic and loss-free medium described by its (scalar and real-valued) electric permittivity and magnetic permeability µwith the corresponding EM wave speed c=(µ)−1/2>0. As demonstrated in Sec. IV, conductive electric and linear hysteresis magnetic losses can subsequently be incorporated via the Schouten-van-derPol theorem [19, p. 1056] by replacing swith [(s+α) (s+β)]1/2, where αand βare non-negative constants (see [19, sec. 26.5]). III. PROBLEM SOLUTION The retarded partial coefficient as expressed through Eq. (1) will be next reformulated via a spatial Fourier representation. To that end, the Green’s function in the s-domain is represented via the wave slowness representation as ˆg(r−r0,s) =s 2iπ2Zi∞ κ=−i∞ exp[−sκ(x−x0)]dκ ×Zi∞ σ=−i∞ exp{−s[σ(y−y0)+γ(κ, σ )|z−z0|]} 2sγ(κ, σ )dσ (3) where γ(κ, σ )=(1/c2−κ2−σ2)1/2with Re(γ)≥0 (4) has the meaning of the wave slowness normal to surfaces Am,n. Making use of the representation (3) in Eq. (1), we get ˆ Pmn(s)=s 2iπ2Zi∞ κ=−i∞ i2 0(sκ1x/2) ×exp[−sκ(xm−xn)]dκZi∞ σ=−i∞ i2 0(sσ1y/2) ×exp[−sσ(ym−yn)] 2sγ(κ, σ )dσ(5) where i0(x) denotes the modified spherical Bessel function of the first kind. Upon expanding the product of the modified spherical Bessel functions into their exponential factors it is found that ˆ Pmn(s) can be expressed in terms of the generic integral analyzed in Appendix VI. Accordingly, we may write ˆ Pmn(s)=S−1 mS−1 nhˆ I(xm−xn+1x,ym−yn+1y,s) −2ˆ I(xm−xn+1x,ym−yn,s) +ˆ I(xm−xn+1x,ym−yn−1y,s) −2ˆ I(xm−xn,ym−yn+1y,s) +4ˆ I(xm−xn,ym−yn,s) −2ˆ I(xm−xn,ym−yn−1y,s) +ˆ I(xm−xn−1x,ym−yn+1y,s) −2ˆ I(xm−xn−1x,ym−yn,s) +ˆ I(xm−xn−1x,ym−yn−1y,s)i(6) Transforming the latter to the TD, we end up with Pmn(t)=S−1 mS−1 nhI(xm−xn+1x,ym−yn+1y,t) −2I(xm−xn+1x,ym−yn,t) +I(xm−xn+1x,ym−yn−1y,t) −2I(xm−xn,ym−yn+1y,t) 148990 VOLUME 8, 2020 M. Štumpf et al.: CdH Technique-Based Computation of Retarded Partial Coefficients: The Coplanar Case +4I(xm−xn,ym−yn,t) −2I(xm−xn,ym−yn−1y,t) +I(xm−xn−1x,ym−yn+1y,t) −2I(xm−xn−1x,ym−yn,t) +I(xm−xn−1x,ym−yn−1y,t)i(7) where the TD function I(x,y,t) is specified in the Appendix via Eq. (21). It is worth mentioning that Eqs. (6) and (7) are valid for arbitrary placed coplanar regions and, thus, they hold not only for mutual partial elements but also for the self-interaction term which is typically the most challenging part to be computed because of the singularity. Other geometrical configurations will be considered in forthcoming works. IV. INCLUSION OF LOSSES Loss effects can be incorporated via the Schouten-Van der Pol theorem [19, p. 1056] for the replacement of swith (s+α)1/2(s+β)1/2as described in the present section. This can be seen upon observing that losses can accounted for by substituting in Eq. (1) the free-space Green’s function of the dissipative 3-D scalar, modified Helmholtz equation (cf. Eq. (2)) ˆg(r−r0,s)=exp[−(s+α)1/2(s+β)1/2|r−r0|/c] 4π|r−r0|(8) where α=σ/ε,β=0/µ are non-negative constants representing the effect of conductive electric and linear hysteresis magnetic losses, respectively. Consequently, assuming that Pmn(t) is available (see Eq. (7)), the retarded coefficient affected by the losses, say P† mn(t), can be represented via ˆ P† mn(s)=Z∞ τ=0exp[−(s+α)1/2(s+β)1/2τ]Pmn(τ)dτ(9) Since the original of exp[−(s+α)1/2(s+β)1/2τ] is known in closed form, the desired TD function immediately follows as P† mn(t) =Pmn(t)exp[−(α+β)t/2] +[|β−α|/2]exp[−(α+β)t/2] ×Zt τ=0I1|β−α|(t2−τ2)1/2/2Pmn(τ)τdτ (t2−τ2)1/2 (10) where I1(x) is the modified Bessel function of the first kind. Apparently, if α=β, then P† mn(t)=Pmn(t)exp(−αt), as expected. Thus, in line with the result represented by Eq. (10), the inclusion of rather general losses amounts merely to evaluating an additional integral with the bounded range of integration. V. NUMERICAL EXAMPLE A popular approximation in PEEC implementations is the midpoint one that amounts to replacing |r−r0|in Eq. (1) by rmn =[(xm−xn)2+(ym−yn)2]1/2. This simplification makes it possible to obtain an approximation of Pmn(t) at once. We next discuss the corresponding results concerning both loss-free and dissipative Green’s functions as given by Eqs. (2) and (8), respectively. A. LOSS-FREE COEFFICIENT Application of the midpoint approximation to Eq. (1) with (2) yields a simple s-domain expression that can be transformed to the TD as Pmn(t)≃δ(t−rmn/c) 4πrmn (11) Equation (11) will next be compared with the exact solution represented by Eq. (7). For this purpose, we define the following parameter I=4πrmn Z2rmn/c τ=0 Pmn(τ)dτ(12) which, under the midpoint approximation, takes the value one. In the following examples, we choose 1x=1y=1.0m, which is supposed to be relatively small (with respect to c×(incident wave pulse time width) or the wavelength at the frequency of analysis). Without any loss of generality, the center of Amis located at the origin, that is, xm=0 and ym=0. In the first example, we take xn=yn=51x, so that surfaces Amand Anare relatively far apart, namely, rmn = 5√21x≃7.071x. The corresponding coefficient (scaled by 4πrmn/c) is shown in Fig. 2a. Apparently, the resulting pulse is relatively narrow, thus resembling the Dirac-delta impulse. Also, the corresponding value of Ias found using the trapezoidal rule is approximately 1.0017, which is very close to the value one pertaining to the midpoint approximation. In the second case, we choose xn=yn=1x, so that the surfaces are relatively near each other, namely, rmn = √21x≃1.411x. In such a case, the computed pulse is not relatively narrow anymore, which indicates inadequacy of the midpoint approximation. Finally note that the corresponding value of parameter Iis 1.0592. B. DISSIPATIVE COEFFICIENT The midpoint rule applied to Eq. (1) with the dissipative Green’s function (8) yields the approximate coefficient P† mn(t)≃δ(t−rmn/c) 4πrmn exp[−(α+β)t/2] +|β−α| 2 rmn/c (t2−r2 mn/c2)1/2exp[−(α+β)t/2] ×I1|β−α| 2(t2−r2 mn/c2)1/2H(t−rmn/c) 4πrmn (13) which consists of the impulsive (instantaneously-reacting) part and the relaxational one forming the tail of the TD coefficient. Apparently, if both α=0 and β=0, which represents the loss-free case, then Eq. (13) boils down to (11). VOLUME 8, 2020 148991 M. Štumpf et al.: CdH Technique-Based Computation of Retarded Partial Coefficients: The Coplanar Case FIGURE 2. Computed loss-free retarded coefficients. Again, for the sake of comparison, we define (cf. Eq. (12)) I†=4πrmn Z2rmn/c τ=0 P† mn(τ)dτ(14) which, under the midpoint approximation, can be expressed as I†≃exp[−(α+β)(rmn/2c)] +|β−α| 2 rmn cZ2rmn/c τ=rmn/c exp[−(α+β)τ/2] ×I1|β−α| 2(τ2−r2 mn/c2)1/2dτ (τ2−r2 mn/c2)1/2(15) The integral in Eq. (15) can be readily evaluated numerically. Again, it is noted that I†takes the value one for the loss-free case. The dissipative coefficient will be evaluated for two distances rmn =5√21xand rmn =√21x, with 1x=1y= 1.0m, again. First, we consider relatively weak electric-type dissipation represented by rmnα/c=0.5 and rmnβ/c=0. In this case, the midpoint approximation of parameter I†as given by Eq. (15) yields I†≃0.80055. The corresponding results are shown in Fig. 3. Their comparison with the lossfree TD coefficients as computed in the previous section V-A clearly show the effect of the non-vanishing losses. It is further observed that the parameter I†for rmn =5√21x and rmn =√21xtakes the values I†=0.80224 and I†=0.85988, respectively, thus justifying the use of the FIGURE 3. Computed dissipative retarded coefficients for rmnα/c=0.5 and rmnβ/c=0. midpoint approximation (13) for relatively distant surfaces, again. Finally, Fig. 4shows the computed TD coefficients for relatively strong dissipation described by rmnα/c=2.5 and rmnβ/c=0.1. Apart from the significant attenuation of the pulse amplitudes, we may observe the change of the resulting pulse shapes, in particular, in their late-time parts. The slow decline of the pulse tails is attributed to the relaxational part of the TD coefficients that is proportional to |β−α|(see Eqs. (10) and (13)). Note that the value of I†as predicted by the midpoint approximation via Eq. (15) is I†≃0.40467 in this case. C. CAPACITIVE COUPLING IN THE FREE SPACE In the last example, the capacitive coupling between regions mand nis considered. They are assumed as squares of size 1 cm on the same plane. The coefficient of potential describing their full-wave electric field coupling in the free space is computed as PFW mn (t)=Pmn(t)/ (16) which corresponds to the CdH-based exact solution. It admits three approximations: •Quasi-static (QS): PFW mn (t)≃PQS mnδ(t) •Quasi-static-delayed (QS −delayed): PFW mn (t)≃ PQS mn δ(t−rmn/c): 148992 VOLUME 8, 2020 M. Štumpf et al.: CdH Technique-Based Computation of Retarded Partial Coefficients: The Coplanar Case FIGURE 4. Computed dissipative retarded coefficients for rmnα/c=2.5 and rmnβ/c=0.1. FIGURE 5. Self induced potential on region ndue to a charge pulse injected in the same region n. •Center-to-center (CC): PFW mn (t)≃1/(4πrmn) δ(t−rmn/c). where PQS mn is the static coefficient of potential [3]. Region n is assumed to be excited by a charge pulse qn(t) of amplitude 1 pC with rise and fall times of 6.67 ps. The self induced potentials are computed using the proposed CdHbased approach, the QS and the QS −delayed approximations through convolution 8mn(t)=PFW mn (t)∗tqn(t),n= 1,m=1,2. For the QS −delayed approximation, it has been assumed rmn =1x=1y=1 cm. Figure 5 shows the self induced potential on region m. It is worth FIGURE 6. Mutual induced potential on region mdue to a charge pulse injected on region nwhen the two regions are 1.4142 1xapart along the diagonal. FIGURE 7. Mutual induced potential on region mdue to a charge pulse injected on region nwhen the two regions are 14.142 1xapart along the diagonal. to observe that the causality is correctly reproduced, while the QS and QS −delayed are inaccurate. Then, region nis moved along the diagonal at a distance rmn =1.41 1xand rmn =14.142 1xin order to compute the near-field and far-field coupling. Figures 6and 7show the mutual induced potentials. As expected, when the two regions are close, the effects of the finite size of the two regions are clearly visible and captured by the proposed CdH-based approach, while the approximated solutions QS,QS −delayed and CC fail. When the two regions are far apart, the QS −delayed and CC approximated solutions are quite close to the proposed CdH-based one because the geometrical dimensions are less important and the propagation delay is well approximated by using the center-to-center distance. On the contrary, the QS approach results to be non-causal since this approximation neglects the propagation delay. VI. CONCLUSION A fundamentally new methodology based on the CagniardDeHoop technique for calculating TD retarded coefficients VOLUME 8, 2020 148993 M. Štumpf et al.: CdH Technique-Based Computation of Retarded Partial Coefficients: The Coplanar Case has been proposed. For the case of two coplanar discretization surfaces, it is demonstrated that this technique yields an exact, closed-form analytical TD expression for the TD surface retarded coefficient of potential. The inclusion of dissipation mechanisms is accomplished with the aid of the SchoutenVan der Pol theorem. Future research directions will address an extension of the methodology to more general problem configurations including noncoplanar and orthogonal elementary surfaces. APPENDIX THE ON-PLANE CASE The TD counterpart of the retarded coefficient under consideration (see Eq. (1)) is given by Eq. (7) in terms of TD functions I(x,y,t). This calls for the inversion of the following integral representation ˆ I(x,y,s)=s 2iπ2Zκ∈K0 exp(sκx) s2κ2dκ ×Zσ∈S0 exp(sσy) s2σ2 dσ 2sγ(κ, σ )(17) for x∈R,y∈Rand {s∈R;s>0}, where K0and S0are the integration paths extending along Re(κ)=0 and Re(σ)=0, respectively, that are indented to the right with semi-circular arcs with centers at the origins and vanishingly small radii (see Fig. 8). FIGURE 8. Complex slowness planes. (a) σ-plane with the CdH-path for y<0; (b) κ-plane with the CdH-paths for x<0. To obtain the TD original of Eq. (17), we shall pursue the CdHapproachcloselydescribedin[20].Accordingly, we first deform S0into the CdH-path, say L∪L∗(here ∗denotes the complex conjugate), where L= {σ(u)= −u(κ)sign(y)+i0}(18) for {1≤u<∞} with (κ)=(1/c2−κ2)1/2, thus representing the loop around the branch cuts extending along {(κ)<|Re(σ)|<∞,Im(σ)=0}. The deformation is permissible by virtue of Jordan’s lemma and Cauchy’s theorem. Consequently, introducing uas the new variable of integration, the inner integral with respect to σcan be expressed as s 2iπZσ∈S0 exp(sσy) s2σ2 dσ 2sγ(κ, σ )=yH(y) 2s(κ) +1 2πs22(κ)Z∞ u=1exp[−su|y|(κ)] du u2(u2−1)1/2(19) where we have combined the contributions from Land L∗ and accounted for the non-zero contribution of the pole singularity at σ=0. Substitution of Eq. (19) back in (17) leads to ˆ I(x,y,s)=ˆ U(x,y,s)+ˆ V(x,y,s) (20) where ˆ U(x,y,s) and ˆ V(x,y,s) are transformed to the TD in the following subsections. Employing these results, we finally end up with I(x,y,t)=|x||y| 4πc|x| 2|y|+|y| 2|x|+c2t2 2|x||y| −c2t2 x2−11/2 −c2t2 y2−11/2 +tan−1"c2t2 x2−11/2# +tan−1"c2t2 y2−11/2#−π 2H(ct −r) +x|y| 2πcc2t2 y2−11/2 −tan−1"c2t2 y2−11/2#H(x)H(ct −|y|) +y|x| 2πcc2t2 x2−11/2 −tan−1"c2t2 x2−11/2#H(y)H(ct −|x|) +xy 2cH(x)H(y)H(ct) (21) which is the main result of the present work, with r=(x2+ y2)1/2>0. This result can be used in Eq. (7) to obtain the desired retarded PEEC coefficient. Its limiting cases follow as I(x,0,t)=[(ct −|x|)2/8π]cH(ct −|x|) +(xct/2π)cH(x)H(ct) (22) and a similar expression applies to I(0,y,t). Finally, upon taking x→0, we get I(0,0,t)=(c2t2/8π)cH(ct) (23) which is useful to computing the coefficient pertaining to overlapping surface elements (cf. Eq. (7)). 148994 VOLUME 8, 2020 M. Štumpf et al.: CdH Technique-Based Computation of Retarded Partial Coefficients: The Coplanar Case GENERIC INTEGRAL ˆ U The first integral to be handled has the following form ˆ U(x,y,s)=1 2πsZ∞ u=1 du u2(u2−1)1/2 ×1 2πiZκ∈K0 exp{−s[−κx+u|y|(κ)]} s2κ2 dκ 2(κ)(24) To obtain the TD counterpart of ˆ U(x,y,s), we deform K0to the CdH path that is defined via −κx+u|y|(κ)=τ(25) where τis the real-valued and positive (time) parameter. The deformation is permissible, again, as the condition for Jordan’s lemma to apply is safely met. Upon solving (25) for κ, we find that the equality is met along the hyperbolic arcs, say G∪G∗, parametrized via G=κ(τ)= − x r2(u)τ+iu|y| r2(u)τ2−r2(u)/c21/2(26) for all τ≥r(u)/c, where r(u)=(x2+u2y2)1/2>0. Further, introducing τas the variable of integration and combining the contributions from Gand G∗, we arrive at ˆ U(x,y,s) =ˆ U0(x,y,s)+1 2π2s3Z∞ u=1 du u2(u2−1)1/2 ×Z∞ τ=r(u)/c Re1 κ2(τ)[κ(τ)]dτ [τ2−r2(u)/c2]1/2 (27) where in the integrand we take the values along Gand ˆ U0(x,y,s) denotes the contribution from the pole at κ=0. The latter is given by ˆ U0(x,y,s)=xH(x)c2 2πs2Z∞ u=1exp(−su|y|/c)du u2(u2−1)1/2 (28) Next, in the second term of Eq. (27), we change the order of integration according to Z∞ u=1duZ∞ τ=r(u)/c dτ→Z∞ τ=r/c dτZU(τ) u=1du(29) in which r=r(1) and U(τ)=(c2τ2−x2)1/2/|y|. Subsequently, the second integral with respect to uis carried out analytically, which after lengthy but straightforward algebra yields ˆ U(x,y,s) =c3y2xH(x) 2πs2Z∞ τ=|y|/c exp(−sτ) c2τ2 dτ (c2τ2−y2)1/2 +3c4x2|y| 4πs3Z∞ τ=r/c exp(−sτ) cτ dτ (c2τ2−x2)3/2 −c4x4|y| 2πs3Z∞ τ=r/c exp(−sτ) c3τ3 dτ (c2τ2−x2)3/2 +c4|x|y2 4πs3Z∞ τ=r/c exp(−sτ) cτ ×3c2τ2−3x2−2y2+2x2y2/c2τ2 c2τ2−x2 dτ (c2τ2−y2)3/2 (30) where the first integral on the right-hand side of Eq. (30) is, in fact, the pole contribution (28), where we substituted τ=u|y|/cfor all τ≥ |y|/c. The integrals in (30) have the form of Laplace integrals, which makes it possible to employ Lerch’s uniqueness theorem [21, appendix] and get the original U(x,y,t) at once. Owing to the factors s−2and s−3in front of the integrals, the thus constructed U(x,y,t) is composed of convolution-type integrals. Since they are amenable to analytical solution, we finally end up with U(x,y,t)=|x||y| 4πc|x| 2|y|+|y| 2|x|+c2t2 2|x||y| −c2t2 x2−11/2 −c2t2 y2−11/2 +tan−1"c2t2 x2−11/2# +tan−1"c2t2 y2−11/2#−π 2H(ct −r) +x|y| 2πcc2t2 y2−11/2 −tan−1"c2t2 y2−11/2#H(x)H(ct −|y|) (31) which is used to derive the main result represented by Eq. (21). GENERIC INTEGRAL ˆ V The second integral to be transformed to the TD has the following form ˆ V(x,y,s)=yH(y) 2 1 2iπZκ∈K0 exp(sκx) s2κ2 dκ (κ)(32) In the first step, again, the indented integration contour K0is deformed to G∪G∗, which now represents the loop encircling the branch cuts {1/c<|Re(κ)|<∞,Im(κ)=0}in the complex κ-plane (see Fig. 8b), that is (cf. Eq. (26)) G= {κ(τ)= −τ/x+i0}(33) for all τ≥ |x|/c. Combining the contributions of integration from Gand G∗, while introducing the time parameter as the variable of integration, we end up with ˆ V(x,y,s)=cxH(x)yH(y)/2s +c3x2yH(y) 2πs2Z∞ τ=|x|/c exp(−sτ) c2τ2 dτ (c2τ2−x2)1/2 (34) VOLUME 8, 2020 148995 M. Štumpf et al.: CdH Technique-Based Computation of Retarded Partial Coefficients: The Coplanar Case where the first term on the left-hand side represents the contribution from the pole singularity at κ=0. The inversion of the first term is straightforward, while the transform of the second term yields a convolution-type integral, again. Carrying out the latter analytically, we arrive at V(x,y,t)=cxyH(x)H(y)H(t)/2 +|x|y 2πcc2t2 x2−11/2 −tan−1"c2t2 x2−11/2#H(y)H(ct −|x|) (35) which is finally used to derive the main result represented by Eq. (21). REFERENCES [1] R. Cicchetti, V. Cicchetti, A. Faraone, and O. Testa, ‘‘Analysis of thin truncated cylinder scatterers using incomplete hankel functions and surface impedance boundary conditions,’’ IEEE Access, vol. 8, pp. 72997–73004, 2020. [2] K. Wang, J.-J. Laurin, Q. Zhang, Q. Zhang, and K. Wu, ‘‘Threedimensional scattering from uniaxial objects with a smooth boundary using a multiple infinitesimal dipole method,’’ IEEE Access, vol. 8, pp. 80842–80854, 2020. [3] A. E. Ruehli, G. Antonini, and L. Jiang, Circuit Oriented Electromagnetic Modeling Using the PEEC Techniques. Hoboken, NJ, USA: Wiley, 2017. [4] A. E. Ruehli and P. A. Brennan, ‘‘Efficient capacitance calculations for three-dimensional multiconductor systems,’’ IEEE Trans. Microw. Theory Techn., vol. 21, no. 2, pp. 76–82, Feb. 1973. [5] A. E. Ruehli and P. A. Brennan, ‘‘Capacitance models for integrated circuit metallization wires,’’ IEEE J. Solid-State Circuits, vol. SSC-10, no. 6, pp. 530–536, Dec. 1975. [6] A. E. Ruehli, ‘‘Inductance calculations in a complex integrated circuit environment,’’ IBM J. Res. Develop., vol. 16, no. 5, pp. 470–481, Sep. 1972. [7] C. R. Paul, Inductance, Loop Partial. New York, NY, USA: Wiley, 2010. [8] D. Romano, G. Antonini, M. D’Emidio, D. Frigioni, A. Mori, and M. Bandinelli, ‘‘Rigorous DC solution of partial element equivalent circuit models,’’ IEEE Trans. Circuits Syst. I, Reg. Papers, vol. 63, no. 9, pp. 1499–1510, Sep. 2016. [9] D. Romano, I. Kovacevic-Badstubner, M. Parise, U. Grossner, J. Ekman, and G. Antonini, ‘‘Rigorous DC solution of partial element equivalent circuit models including conductive, dielectric, and magnetic materials,’’ IEEE Trans. Electromagn. Compat., vol. 62, no. 3, pp. 870–879, Jun. 2020. [10] C. A. Balanis, Advanced Engineering Electromagnetics. New York, NY, USA: Wiley1989. [11] D. Gope, A. E. Ruehli, C. Yang, and V. Jandhyala, ‘‘(S)PEEC: Timeand frequency-domain surface formulation for modeling conductors and dielectrics in combined circuit electromagnetic simulations,’’ IEEE Trans. Microw. Theory Techn., vol. 54, no. 6, pp. 2453–2464, Jun. 2006. [12] D. R. Wilton, J. Rivero, W. A. Johnson, and F. Vipiana, ‘‘Evaluation of static potential integrals on triangular domains,’’ IEEE Access, vol. 8, pp. 99806–99819, 2020. [13] A. Bellen, N. Guglielmi, and A. E. Ruehli, ‘‘Methods for linear systems of circuit delay differential equations of neutral type,’’ IEEE Trans. Circuits Syst. I, Fundam. Theory Appl., vol. 46, no. 1, pp. 212–215, Jan. 1999. [14] W. Pinello, A. C. Cangellaris, and A. Ruehli, ‘‘Hybrid electromagnetic modeling of noise interactions in packaged electronics based on the partialelement equivalent-circuit formulation,’’ IEEE Trans. Microw. Theory Techn., vol. 45, no. 10, pp. 1889–1896, Oct. 1997. [15] P. J. Restle, A. Ruehli, S. G. Walker, and G. Papadopoulos, ‘‘Full-wave PEEC time-domain for the modeling of on-chip interconnects,’’ , vol. 20, no. 7, pp. 877–887, Jul. 2001. [16] J. Nitsch and F. Gronwald, ‘‘Analytical solutions in nonuniform multiconductor transmission lines theory,’’ IEEE Trans. Electromagn. Compat., vol. 41, no. 4, pp. 469–479, Nov. 1999. [17] L. Lombardi, G. Antonini, and A. E. Ruehli, ‘‘Analytical evaluation of partial elements using a retarded Taylor series expansion of the Green’s function,’’ IEEE Trans. Microw. Theory Techn., vol. 66, no. 5, pp. 2116–2127, May 2018. [18] A. T. De Hoop, ‘‘A modification of Cagniard’s method for solving seismic pulse problems,’’ Appl. Sci. Res., Sect. B, vol. 8, no. 1, pp. 349–356, Dec. 1960. [19] A. T. De Hoop, Handbook of Radiation and Scattering of Waves. London, U.K.: Academic Press, 1995. [20] M. Štumpf, Time-Domain Electromagnetic Reciprocity in Antenna Modeling. Hoboken, NJ, USA: Wiley, 2019. [21] M. Štumpf, Electromagnetic Reciprocity in Antenna Theory. Hoboken, NJ, USA: Wiley, 2018. MARTIN ŠTUMPF (Senior Member, IEEE) received the Ph.D. degree in electrical engineering from the Brno University of Technology (BUT), Brno, Czech Republic, in 2011. After his Ph.D. research, he spent a year and a half as a Postdoctoral Fellow of the ESAT-TELEMIC Division, Katholieke Universiteit Leuven, Leuven, Belgium. He is currently an Associate Professor with the Department of Radioelectronics, BUT. During three months in 2018, he was a Visiting Professor at the UAq EMC Laboratory, University of L’Aquila, Italy. He has authored the books Electromagnetic Reciprocity in Antenna Theory (Wiley–IEEE Press, 2017), Pulsed EM Field Computation in Planar Circuits: The Contour Integral Method (CRC Press, 2018), and Time-Domain Electromagnetic Reciprocity in Antenna Modeling (Wiley–IEEE Press, 2019). His main research interests include the mathematical modeling of electromagnetic wave phenomena with an emphasis on EMC and antenna engineering. GIULIO ANTONINI (Senior Member, IEEE) received the Laurea degree (cum laude) in electrical engineering from the University of L’Aquila, L’Aquila, Italy, in 1994, and the Ph.D. degree in electrical engineering from the University of Rome ‘‘La Sapienza,’’ in 1998. Since 1998, he has been with the UAq EMC Laboratory, University of L’Aquila, where he is currently a Professor. He has coauthored the book Circuit Oriented Electromagnetic Modeling Using the PEEC Techniques (Wiley–IEEE Press, 2017). His scientific interest is in the field of computational electromagnetics. ALBERT RUEHLI (Life Fellow, IEEE) received the Ph.D. degree in electrical engineering from the University of Vermont, Burlington, VT, USA, in 1972, and the Honorary degree from Lulea University, Lulea, Sweden, in 2007. He is currently an Emeritus of IBM Research and an Adjunct Professor in the EMC area with the Missouri University of Science and Technology, Rolla, MO, USA. He has authored or coauthored more than 250 technical articles. He is a member of the SIAM. He is the editor of two books ‘‘Circuit Analysis,’’ ‘‘Simulation and Design’’ (North Holland 1986, 1987), and co-author of a book ‘‘Circuit Oriented Electromagnetic Modeling Using the PEEC Techniques’’ (WileyIEEE Press, 2017). 148996 VOLUME 8, 2020