2009 International Workshop on Nonlinear Circuits and Signal Processing
NCSP'09, Waikiki, Hawaii, March 1-3, 2009

A Spice-Oriented Frequency Domain
Analysis of Electromagnetic Fields of PCBs
Akiko Kusaka†, Takaaki Kinouchi†, Yoshihiro Yamagami†, Yoshifumi Nishio†, Akio Ushida‡†Department of Electrical and Electronic Engineering,
Tokushima University, Japan‡Department of Mechanical and Electronic Engineering,
Tokushima Bunri University, Japan

Abstract
For designing PCBs (printed circuit boards), it is very important to find out the locations gradating strong EMF (electromagnetic fields), where the voltages will have large peakvalue. In the peak search of the frequency characteristics, it
is diﬃcult to find the peak for large scale circuits. If we carryout frequency analysis as changing the frequency, we have toincrease the frequency gradually. However, if the changing
step is large, we may miss the peaks. Also there is no guarantee that the peak values exist at the frequencies. We needto find out both the exact peak voltages and locations on thePCB.
In this article, we propose an algorithm to find the peak of
the frequency characteristics by using the Sine-Cosine circuit[1]-[3] and nonlinear limiter to control the step size of SPICE.
1. Introduction
Nowadays, sizes of LSIs and PCBs become smaller andsmaller, and the operating frequency higher and higher. In
these cases, PCBs are modeled by linear LRCG large scaleplane circuits as shown in Fig. 1. They can be solved by timedomain and/or frequency domain [4][5] and we concentrateon the latter technique in this study. In this case, the frequency
response curve may have many sharp resonant points for highQ circuits. It is diﬃcult to find the exact peak points.
The resonant and anti-resonant points appear at poles andzeros of driving and transfer functions, respectively. We can
not find out them with analytical method when the circuit size
becomes larger. They are also found by tracing the frequency
response curve |V(ω)|, where the resonant and anti-resonant
points correspond to the highest and lowest values. In thiscase, when we use AC analysis of Spice simulator, we may
pass over to find out all the resonant points especially for thehigh Q circuits. Furthermore, it is impossible to find out theexact resonant voltages. Therefore, we need to propose a newalgorithm based on the HB method where the determining
equation is formulated by the coupled Sine-Cosine circuits[1]-[3]. Both the resonant and anti-resonant points correspond to the solutions satisfying f 0 (ω) = 0 on the frequency

Lrj(t)Crg

L

Figure 1: LRCG plane model of PCB.

response curve, so that we use a diﬀerentiator and a nonlinear limiter to find the f 0 (ω) = 0 points. The circuit is solvedby the variable step size transient analysis of Spice. Thus,
we can develop an eﬃcient simulator to trace the frequency
response curves, and to find the exact peak voltages.
2. Reactance circuit
Now, we review the circuit properties of a general linear
LRC circuit as shown Fig. 1. The circuit equation is given by Y11
 Y
 21 ..
 .
Yn1

Y12
Y22..
.
Yn2

 

... Y1n   V1   J   ... Y2n   V2   0  =  ..   ..   .. .   .   . 0Vn... Ynn(1)where Vs are nodal voltages for the input J, and the mutual
admittance is given by
Yi j = jωCi j + Gi j +Then, we have

- 526 -

Vk =

∆1k J
∆

1jωLi j(2)

(3)

where∆=Y11Y21..
.

Y12Y22..
.

Yn1

Yn2... Y1n... Y2n..
.

Then, the current iC is given by... YnndvC= − ωCVLC sin ωt + ωCVLS cos ωtdt

(9.2)

Thus, the coeﬃcients of sin ωt, cos ωt are described byNow, we consider the case r = g = 0 in Fig. 1 which corresponds to the reactance circuit. It is known that the voltage(3) can be written byn a0k ∑
ja2i,k ω
Vk (ω) = ++ ja∞,k ω J(5)22jω i=1 ω2i − ωObserve that the resonant points correspond to the poles ofboth the driving and the transfer functions. They arise at thesame frequencies. They also satisfy
0 < ω1 < ω2 < ω3 < . . . < ∞iC = C(4)

ILS = −ωCVLC , ILC = ωCVLS(9, 3)Namely, the capacitor is replaced by the coupled voltagecontrolled current sources in the Sine-Cosine transformation.
The circuit corresponding to the determining equation of theHB method is driven by a constant voltage E = J and ω isgiven by (7).
Node voltage(6)Sine-CosineCircuit

Note that r and g in the high frequency domain are small compared to the reactance values, so that the frequency response
behaves like as that of reactance circuits, and it has sharp resonant points.

ω3. Peak points tracing algorithm

Figure 2: Peak detector of frequency response curve.

Although, for relative low Q circuits, the AC analysis ofSpice and a cubic spline combining Newton method can be
usefully applied to find out the peak points, it may pass overthe points for high Q circuits. Thus, we propose a new algorithm based on the HB method such that the determining
equation is solved by transient analysis of Spice, where ω isa function of time t as follows:ω = Kt(7)

Now, let us discuss the Sine-Cosine circuit corresponding tothe determining equation of the HB method.
Let the current through an inductance L be
iL = ILS sin ωt + ILC cos ωtE
Differentiator Nonlinear limitter

Since ωs at the peak voltages satisfyd|Vk (ω)|= 0, k = 1, 2, . . . n,dω(10)on the response curve, we need to find the zero points satisfying (10) Hence, |Vk (ω)| need to be firstly diﬀerentiated by adiﬀerentiator. In order to detect the exact peak points, the output is limited and expanded with a nonlinear limiter as shownin Fig. 2, which consists of a limiter and nonlinear diodes.
The output of the limiter is given by
−VmaxkvinvL =  V(8.1)

max:::for vin < −VLfor − VL ≥ vin ≥ VLfor vin > VL

(11.1)

Then, the voltage vL is given bydiL= − ωLILC sin ωt + ωLILS cos ωtvL = Ldt

whereThus, the coeﬃcients of sin ωt, cos ωt are described by
VLS = −ωLILC , VLC = ωLILS(8, 3)

Namely, the inductance is replaced by coupled currentcontrolled voltage sources in the Sine-Cosine transformationof the HB method. In the same way, let the voltage across acapacitor C bevC = VCS sin ωt + VCC cos ωtkVL = Vmax

(8.2)

(9.1)

The outputs of diodes are given by{I s exp(λvo ):for vo > 0io =−I s exp(−λvo ) :for vo < 0(11.2)(11.3)

This means that the regions around d|Vk (ω)|/dω = 0 arelargely expanded. Furthermore, the characteristic has large
nonlinearity around the zeros. Thus, the transient analysis
around the zero points is executed with a very small step size,and we can find out precise peak points.

- 527 -4. Illustrative example
4.1. Theoretical analysis of transmission line

4.2. Discrete model of transmission line

We consider a single transmission line of the length d[mm],whose parameters arer[Ω/mm], L[H/mm], g[S /mm], C[F/mm].
The circuit equation with a complex frequency s at x[mm]from the near end is given bydV(x, s)= −(r + sL)I(x, s)dx(12.1)

dI(x, s)= −(g + sC)V(x, s)dx(12.2)

d2 V(x, s)= (r + sL)(g + sC)V(x, s)
dx2(13.1)

Thus, we haved2 I(x, s)= (r + sL)(g + sC)I(x, s)
dx2

Transmission line is usually modeled by discrete RLCG
ladder circuit, where we neglected G. This time, we analyze the two-dimensional circuit in the case of 5×5 as shownin Fig.1, where C = 1[nF], L = 1[mH], e(t) = E cos ωtand E = 100[µV]. In order to obtain the frequency response
curve, the circuit is transformed by HB method into the SineCosine circuit. Note that the detector can be attached toan arbitrary node. We use a simple RC diﬀerentiator with
Cd = 1[pF] and Rd = 1[Ω] in the detector. Now, we show
the simulation results using the transient analysis of Spice.
Figures 3(a), (b) and (c) show the results for the cases of
R = 1[Ω], 10[Ω] and 100[Ω], respectively, where we set thepeaks 1-12 from the left hand side. The sharpnesses dependon the quality factors so that we have changed them by Rs. We
found that although the resonant frequencies are almost samefor all Rs as shown in Table 1, the peak values are largely
diﬀerent.

R=1[Ω]
f [kHz]V[mV]f [kHz]V[mV]f [kHz]V[mV]
R=10[Ω]
f [kHz]V[µV]
f [kHz]V[µV]
f [kHz]V[µV]
R=100[Ω]f [kHz]V[µV]
f [kHz]V[µV]
f [kHz]V[µV](13.2)

Set the boundary conditions V(0, s), I(0, s) at the near endand V(d, s), I(d, s) at the far end. Then, we have[I(0, s)
I(d, s)]=1  coth γ(s)dZ0 (s) − sinh 1γ(s)dwhere[]− sinh 1γ(s)d  V(0, s)V(d, s)coth γ(s)d(14)Z0 (s) =√(r + sL)/(g + sC)

γ(s) =√(r + sL)(g + sC)

Thus, the poles are located at the frequencies satisfying
sinh γ(s)d = 0, namely(r + sL)(g + sC) = ( jnπ/d)2 , n = 1, 2, . . .

(15)

rp0 = − , pn = un ± jvn , n = 1, 2, . . .
L(16)Hence,

where
un = −

Lg + rC,2LC√4LC(gr + (nπ/d)2 ) − (Lg + rC)2vn =.
2LC
Observe that the transmission line has an infinite number ofpoles, where the frequency response curve has peak points.

Table 1 Calculated peak values.
peak 1peak 2peak 3
30.7
108.7
149.2
8.0
1.6
5.5
peak 5peak 6peak 7
220.7
259.7
280.2
3.4
1.9
3.4
peak 9 peak 10 peak 11
356.4
397.9
428.5
0.80
0.37
0.08
peak 1peak 2peak 3
30.4
108.2
149.3
1169.7847.5
943.3
peak 5peak 6peak 7
220.5
259.6
282.0
970.9
309.3
204.4
peak 9 peak 10 peak 11359.5
400.9
427.9
19.6
18.4
14.5
peak 1peak 2peak 3
25.6
104.8
147.2
144.5
101.5
107.7
peak 5peak 6peak 7
220.3
260.0
282.2
106.0
44.9
62.3
peak 9 peak 10 peak 11
353.1
394.0
424.1
29.6
12.2
4.0

peak 4192.0
4.7
peak 8303.6
1.1
peak 12477.7
0.004
peak 4191.9
572.8
peak 8320.4
227.8
peak 12477.3
3.6
peak 4
190.4
69.3
peak 8308.2
52.5
peak 12475.5
3.7

It seems that there are small diﬀerences in the calculated
peak points for R = 1 and for R = 100. It is due to pass overthe point because the peak is too narrow and small to detecteven when the nonlinear limiter is introduced. In this case,- 528 -we need to trace the resonant curve around peak points for
R = 1 with a selection of much smaller initial step size.
5. Conclusions and remarks
In this study, we have proposed an algorithm to find thepeaks of the frequency response curve by combining the SineCosine circuit based on the HB method. The circuit is tracedby the transient analysis of Spice. Since the peak points correspond to the gradient being to zeros, we diﬀerentiate thecurve and find its zero points. In order to find the exact zeropoints, we have developed the peak detector using nonlinear
limiter such that the regions are largely expanded. Thus, thestep size around zero values becomes small when we use variable step size transient analysis.

0.025
0.02
0.015
0.01

V[V]

0.005
0-0.005-0.01-0.015-0.02-0.02501e-009

2e-0093e-0094e-009Acknowledgments

5e-009Time[s](a)This research was partly supported by the fund of Open Research Center Project from MEXT of Japanese Government(2007-2011).

0.00140.00120.001
0.0008References

V[V]0.00060.00040.0002[1] A. Ushida, Y. Yamagami and Y. Nishio, “Frequency responses of nonlinear networks using curve tracing algorithm,”Proc. of ISCAS’02, vol.I, pp.641-644, 2002.

0-0.0002-0.0004-0.0006[2] A. Ushida and M. Tanaka, Simulations of Electronic
Circuits, Corona Publishing Co., Ltd, 2002.

-0.0008
01e-009

2e-0093e-009

4e-0095e-009Time[s](b)[3] T. Kinouchi, Y. Yamagami, Y. Nishio and A. Ushida,“Spice-oriented harmonic balance Volterra series methods,” Proc. of NOLTA’07, pp.513-516, 2007.

0.000160.000140.000120.0001[4] T. Tang and M.S. Nakhla, “Analysis of high-speed VLSIinterconnect using the asymptotic waveform evaluation technique,” IEEE Trans. Computer-Aided Design,
Vol.14, pp.341-352, 1992.

[V]8e-005

6e-0054e-0052e-005[5] Y. Tanji, Y. Nishio and A. Ushida, “Analysis of nonuniform and nonlinear transmission lines via frequency
domain techniques,” IEICE Trans. on Fundamentals,
Vol.E79-A, pp.1486-1494, 1996.

0-2e-005
01e-009

2e-0093e-0094e-0095e-009Time[s](c)Figure 3:
Frequency characteristics of the circuit.
(a) R=1[Ω]. (b) R=10[Ω]. (c) R=100[Ω].K = 2π × 1014

[6] Y. Nakano, Electronic Circuit Theory, Hirokawa publishing Co., 1976 (Japanese).
[7] A.S. Sedra and K.C. Smth, Microelectronic Circuits,
Oxford Univ. Press, 2004.
[8] H. Hasegawa, et al., Numerical Computation Methods,
Ohmsha publishing Co., 1998 (Japanese).

- 529 -