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

Spice-Oriented Nonlinear Circuit Analysis Using Harmonic Balance Method
Takaaki Kinouchi† Yoshihiro Yamagami† Yoshifumi Nishio† and Akio Ushida‡†Tokushima University
2-1 Minami-Josanjima, Tokushima, Japan
Phone:+81-88-656-7470
FAX:+81-88-656-7471
Email: {t-kino, yamagami, nishio}@is.tokushima-u.ac.jp‡Tokushima Bunri University
1314-1, Sido, Sanuki, Kagawa, Japan
Email: ushida@fe.bunri-u.ac.jp

1. Abstract

form of schematic and/or net-list of Spice. The circuit can besolved by Spice simulator and the frequency characteristics
To analyze RF circuits and modulators is very important for are easily obtained.
designing integrated and communication circuits. They are
driven by multiple frequencies, one of which is usually very
high carrier frequency compared to the other. To know the 3. Nonlinear devices of HB combining with MATLAB
transient behaviors, we need to calculate many carrier waveIn this section, we propose our HB devices combining
forms, so that the transient analysis is very time-consuming.
Spiceand MATLAB. Analog integrated circuits usually conHence, we propose an efficient envelope analysis for calcusistofmany kinds of nonlinear devices such as diodes, bipolating the asymptotic behaviors of the amplitudes, which islartransistors
and MOSFETs, whose Spice models are debased on the harmonic balance (HB) method with the slowly
scribedby
theseveral special functions such as exponential,
varying coefficients. In order to develop the Spice-oriented
square-root,
piecewise-continuousfunctions and so on. For
simulators, the Fourier expansions of nonlinear devices suchthese
specialfunctions,theFourier
expansions cannot be exas bipolar transistors are executed with MATLAB, and theecutedintheanalyticalforms.
Therefore,
we need to apFourier modules should be stored in our computer library.
proximate
them
with
the
Taylorexpansions.
For an example
Thus, we can easily formulate the determining equations ofofNPN
bipolartransistor
shownbyFig.1(a),the Ebers-Moll
HB method called Sine-Cosine circuit.
model is given by (b).
2. IntroductionCC

For designing integrated circuits and communication systems, it is very important to analyze the characteristics of nonlinear circuits. In these systems driven by multi-frequency inputs. There are various method to solve the characteristics. In
this paper, we propose an efficient envelope analysis. The frequency components of communication systems usually consist of very high carrier frequency and low signal. Then, thetransient analysis based on the numerical integration technique is very time-consuming, because we need to calculate
many carrier waveforms to know the asymptotic transient behaviors. Our envelope analysis is a modification of HB (harmonic balance) method to the high carrier frequency components, where we assume that the Fourier coefficients are
slowly varied by the low frequency signal. The harmonic balance (HB) is a well-known method for the frequency domain
analysis, which gives good results even for relatively strong
nonlinear circuits. In our HB, which combining MATLAB
and Spice, it is largely improved such that it can be easily
applied to large scale nonlinear circuits, where the determining equation called Sine-Cosine circuit can be given in theBαF id 1

id 2

id 1

αRid 2BE(a) NPN Transistor

E(b) Ebers-Moll modelFig.1 NPN transistor Ebers-Moll model.
The characteristics of two diodes with Ebers-Moll models are
given byvd
(1)
id = IS exp( ),VTwhose parameters are given, for example, IS = 10−12 , VT =0.025. Applying the Taylor expansion, we have
id = IS exp(

vd0
λ2
λ3
) 1 + λvd + v2d + v3d + . . .
VT
2!
3!

!
(2)

for a small variation vd at the operating point vd0 , and λ =- 530 -

exp( V1T ). We assume the input waveform as follows;
vd (t) = V0d +KXVd,2k−1 cos kωt + Vd,2k sin kωt ,Then, the inductor voltage and capacitor current are given byvL (t) = L

(3.1)

k=1

diL
= ωc LIS (t) cos ωc t − ωc LIC sin ωc t
dt

+L

where K is the maximum higher harmonic component to take
into account. Then, the output waveform of (1) is calculated.
However, these calculation is very large and complex. Therefore, we use MATLAB and calculated as follows;dIC
dIS
cos ωc t + L
sin ωc t,
dt
dt

(8.1)

dvC
= ωcCVS (t) cos ωc t − ωcCV IC sin ωc t
dt
dVC
dVS
+C
cos ωc t + C
sin ωc t.
(8.2)
dt
dt
K
X

id (t) = I0d +Id,2k−1 cos kωt + Id,2k sin kωt ,(3.2) Using (8.1) and (8.2), we set the Cosine-component and Sinecomponent of Fig.3 in the manner that they respectively satk=1
isfy the circuit equations. Thus, the Sine-Cosine circuit arein the symbolic form. Namely, (I0d , Id2 , . . . , Id,2K ) are func- obtained as shown in Fig.4.
tions of (V0d , Vd2 , . . . , Vd,2K ), and the input-output relations
ω cC Vare realized by ABMs of Spice;- ω cC VSCIk,C = Iˆk,C (V0C , . . . , V2K,C , V0B , . . . , V2K,B , V0E , . . . , V2K,E )Ik,B = Iˆk,B (V0C , . . . , V2K,C , V0B , . . . , V2K,B , V0E , . . . , V2K,E )Ik,E = Iˆk,E (V0C , . . . , V2K,C , V0B , . . . , V2K,B , V0E , . . . , V2K,E )k = 0, 1, 2, . . . , 2K

iC (t) = C

L(4)

ω cL IS

L

IC- ω c L IC

ISC+ VC -

C
+ V

R

S

Vsinω s t
Sine

CosineMth Sine

Mth Cosine

Cosine circuit-

R

Sine circuit

Fig.4 Sine-Cosine circuit of RLC.
DC component

Mth higher harmoniccomponent

Fundamental component

Fig.2 Symbol of NPN transistor Fourier module.
It is convenient to store the Fourier modules in our computer
library. Replacing each element by the corresponding HBmodel, we can formulate the “Sine-Cosine circuit that is corresponding tothe HB determining equation. Resulting circuit
is efficiently solved by Spice.

They are coupled with the voltage-controlled current sources
and current-controlled voltage sources in each other. From
the relation (6), the voltage V sin ω s t in Fig.4 corresponds tothe envelope of the forced input which is slowly varying because of ωc  ω s . Now, we show that our envelope approachis validity when the difference between two input frequencies
is large. The parameter values are set as follow;
L = 0.1, C = 0.01, R = 0.01, ωc = 50, V = 1,
and the other frequency is changed to ω s = 2, 5, 7, 10.
3

3Envelope analysis 2

RLC envelope analysis 1

4. Envelope analysis
4.1 RLC circuit Now, consider a simple example RLC
circuit Fig.3 driven by

2

2

1

1

0

0

-1

-1

-2

-2

-3

v(t) = V sin(ω s t) sin(ωc t), ωc  ω s .

0

(6)

1

2

3

4

5

,

-3
0

(a) e(t) = sin 2t sin 50t

1

2

3

4

5

(b) e(t) = sin 5t sin 50t
3L

C

i L (t)+ v (t) -Envelope analysis 3

Envelope analysis 4

4

2

2

v(t)

C

1
0

0

R

-1

-2

-2

-4
0

Fig.3 A simple RLC circuit driven by a modulated signal.

1

2

3

(c) e(t) = sin 7t sin 10t

4

5

,

-3
0

1

2

3

4

5

(d) e(t) = sin 10t sin 50t

Fig.5 Envelopes of currents in RLC circuit of Fig.4.
We assume the inductor current and capacitor voltage as fol4.2 Nonlinear RLC circuit Next, we consider nonlinlows;
ear RLC circuit driven by frequency modulated (FM) signal
iL (t) = IC (t) cos ωc t + IS (t) sin ωc t,(7) shown Fig.6.
vC (t) = VC (t) cos ωc t + VS (t) sin ωc t

- 531 -

R
+ vL -ihave Sine-Cosine circuit as show in Fig.8. So these parameters are given by Fig.8(c)-(f).
+

vq

e(t)

IC

-

i = γφ + δφ3(9.2)
δ = 0.

Assume the waveforms as follow;
q = Qc cos ωt + Q s sin ωt,
φ = Φc cos ωt + Φ s sin ωt

S

.

+

QC

rd

QC
ωQ S
(c) Capacitor cosine element

IS

1[GΩ] 1[F](10.1)

1[µΩ](10.2) VC

S

rd

C

+
S

Q
S
-ωQ C
(d) Capacitor sine element
rs

ΦC

.
+Φ -

Q

1[GΩ] 1[F]

1[H]

rs

S

(b) Sine circuit

.

IC

S

2

2
C

(α+0.75β(Q +Q ))QS(a) Cosine circuit

(9.1)

2

(γ+0.75δ(Φ C+Φ ))Φ

R

+ V S
2

2
C

vq = αq + βq3 ,γ = 0.142857,

2

IS

C

(α+0.75β(Q +Q ))Q C

Let the nonlinear capacitor and inductor be are

β = 7,

2

S

+ V C

E

e(t) = E cos ωt, E = 1[V], R = 1[Ω], ω = 0.01t[rad/sec]
Fig.6 Nonlinear RLC circuit

α = 1,

2
C(γ+0.75δ(Φ +Φ ))Φ

R

1[H]

.

ωΦ S

VS

1[µΩ] + Φ S -ΦS- ωΦ CThe characteristic of nonlinear capacitor is shown in Fig. 7.
(e) Inductor Cosine element

800(f) Inductor sine element

q[C]
400

Fig.8

Sine-Cosine circuit of nonlinear RLC

0

For the capacitor element, equations (13) are realized by Fig.8(c), (d). Where we assume the voltages of two capacitors 1[F]are Qc and Q s . In the same way, the sine-cosine inductor
−800
−4
−2
0
2
4 V[V]
circuits are given by Fig.8 (e), (f). Where we substitute small
dummy resistors r s and large Amy resistors rd to avoid C-E
Fig. 7 Characteristic of nonlinear capacitor
loops and L-J cut sets. Fig.8 is corresponds to the determining
Next, substituting (10.1) into (9.1) and (10.2) into (9.2) re- equation of the HB method.
spectively, we have fundamental harmonic components asfollows,1−400I

vq = (α + 0.75β(Q2c + Q2s ))Qc cos ωt+ (α + 0.75β(Q2c + Q2s ))Q s sin ωt,

0.8(a) Analysis by HB
Stable(11.1)
0.6
Going downi = (γ + 0.75δ(Φ2c + Φ2s ))Φc cos ωt+(γ + 0.75δ(Φ2c + Φ2s ))Φ s sin ωt.

Going up (b) Envelope analysis

Unstable

0.4(11.2)
0.2

neglecting the higher harmonic components. We have from
(10.1) and (10.2).
On the other hand,(c) Envelope analysis(reverse)

0
0

0.5

1

1.5

2
ω[rad/sec]

Fig.9 Comparison between the HB and Envelope analysis.

dq= (Q̇c + ωQ s ) cos ωt + (Q̇ s − ωQc ) sin ωt,
dt(12.1) Thus, we obtain the envelope analysis by solving the circuitwith transient analysis of Spice shown in Fig.9. The curve (a)is the jump phenomena where the dotted line is unstable anddΦvL == (Φ̇c + ωΦ s ) cos ωt + (Φ̇ s − ωΦ s ) sin ωt. (12.2) solid line stable. When the frequency is changed by ω=0.01t,
dt
we have phenomena (b) which has oscillations after jumping
Terefore, sine and cosine components are as follows;
down. When the frequency is changed ω=2-0.01t, we havephenomena(c) which has oscillations after jumping up.
Ic = Q̇c + ωQ s , I s = Q̇ s − ωQc ,(13)
4.3 RF amplifier Consider a simple RF amplifier shown
Vc = Φ̇c + ωΦ s , V s = Φ̇ s − ωΦ sin Fig.10 whose input is assumed by
Qc , Q s , Φc , and Φ s are unknown parameters, therefore it isnecessary to make a circuit to satisfy equation (13). Thus, wev(t) = V1 cos ω1 t + V2 cos ω2 t, for ω2 = ω1 + ∆ω, (14.1)
i=- 532 -

for a small ∆ω. Then, it is written by5

4
3

Envelope analysis

3

2

R

0

-2
-2

-3

-3

-4

Transient analysis
0.004

0.008
0.012
time[msec]

0.016

0.02

-5
0

0.04

Transient analysis
0.08
0.12
time[msec]

0.16

0.2

(a) ∆ω = 0.1,(b) ∆ω = 0.01
Fig.12 Comparisons between the transient and envelope analysis.

L voutThus, the asymptotic behavior of the amplitudes can be calculated by transient analysis of Spice as shown in Fig.12.

-

C0
Q

1

-1

-1

+
C

V[V]

0

-4
0

R1

2

1
V[V]

v(t) = (V1 + V2 cos ∆ωt) cos ω1 t − (V2 sin ∆ωt) sin ω1 t
= VC (t) cos ω1 t + VS (t) sin ω1 t
(14.2)

Envelope analysis

4E

e 1 (t)Table 1 Comparison calcurational time between Envelope andTransient analysis

R2
e (t)
2

∆ω[Mrad/sec]
0.1
0.07
0.05
0.03
0.01

Fig.10 Simple emitter follower amplifier.
C0 = 1[nF], C = 1[nF], L = 10[µF], R1 = R2 = 50[kΩ]
R = 1[kΩ], E = 10[V], ω1 = 107 [rad/sec]
We have used the Ebers-Moll model for the transistor,
id = 10−12 exp[40vd ], αF = 0.99, αR = 0.3,

(15)

and approximate it by Taylor expansion. Then, the Fourier
module can be obtained by MATLAB in the symboric forms
as way of section 2. Thus, the Sine-Cosine circuits of RF
circuit can be given in Fig.10.

Envelope[sec]
3.75
4.83
5.5
5.53
11.12

Transient[sec]
20.45
29.48
40.55
65.05
189.86And calculational time as shown Table 1. Observe that between transient and envelope. When ∆Ω = 0.1, our envelope
analysis is faster about 5 times. However, when ∆ω = 0.01,
faster than about 17 times. Therefore, when ω2 nearly to ω1
namely, ∆ω become more small, envelope analysis is moreequal as shown in Fig.12. And calculational time of our envelope analysis is faster than direct transient analysis shownin Table 1.
5. Conclusions

In this paper, we propose a Spice-oriented envelope analysis which is based on Harmonic Balance method. The coefficients of HB set to slowly varying. Thus, we can obtained efficient asymptotic analysis result. Fourier expansion applies
to nonlinear devices by use MATLAB, and the determining
equation of HB is replaced Sine-Cosine circuit. Therefore,
calculation time of our envelope analysis is faster than direct
Fig.11 Sine-Cosine circuit of the amplifier shown by Fig.10. transient analysis. Also we can analyze the circuit consists
nonlinear inductor and capacitor. So, it is promising to analThe DC circuit shown by Fig.11(a) only has a DC volt- ysis more expansion for nonlinear circuits.
age source, where the input source is removed by the shortReferences
circuit. The Sine-Cosine circuits have slowly varying inputs
δ(k − 1)VC (t), and δ(k − 1)VS (t) to the Cosine and Sine subcircuits, respectively. Observe that all the circuit configura- [1] P.Wambacq and W.Sansen, Distortion Analysis of Analog
Integrated Circuits, Kluwer Academic Pub., 1998.
tions corresponding to the kth harmonic components are thesame topology. Therefore, we can easily formulate the cir- [2] K.K.Clarke and D.T.Hess, Communication Circuits:
cuits by “copy and paste” with Spice simulator. Ih this study,
Analysis and Design, Addison-Wesley Pub. Co., 1971.
we consider maximum harmonic component K=3. Thus, DC
circuit and each three of Sine-Cosine circuit, total number [3] A.Ushida, Y.Yamagami and Y.Nishio, “Frequency responses of nonlinear networks using curve tracing algoof severn circuits are created. All the sub-circuits are courithm,”ISCAS 2002, vol.I, pp.641-644, 2002.
pled with the controlled sourced and the Fourier transistor
modules. The DC circuit is also time-varying. Note that the [4] T. Kinouchi, Y. Yamagami, Y. Nishio, J. Kawata, and A.
resultant Sine-Cosine circuits correspond to the determining
Ushida, “Spice-Oriented Harmonic Balance and Volterra
equations of HB method and it can be solved by the transient
Series Methods,” Proc. of NOLTA’07, pp.513-516, 2007.
analysis of Spice.
kω L I S

+

C

R

L

kω C 0V0S

R1

R1

VC
-

C

R

- kω C 0V0C

E+

δ(k-1) v (t)C(a) DC circuit

Q

C0

-

VS

-

- kω L I C

C

R

L

IS

- kω C VC

C0+ V 0S

R2

R2

R1

IC

kω C VS

Q

+

L

Q

R2

kω

δ(k-1)vS (t)(b) kth Cosine circuit (k=1,2,...,K)(c) kth Sine circuit (k=1,2,...,K)- 533 -