IEEE TRANSACTIONS ON COMPUTER-AIDED DESIGN OF INTEGRATED CIRCUITS AND SYSTEMS

1

Simulation of Mutually Coupled Oscillators Using Nonlinear Phase Macromodels Davit Harutyunyan, Joost Rommes, Jan ter Maten, and Wil Schilders

Abstract—Design of integrated RF circuits requires detailed insight in the behavior of the used components. Unintended coupling and perturbation effects need to be accounted for before production, but full simulation of these effects can be expensive or infeasible. In this paper, we present a method to build nonlinear phase macromodels of voltage-controlled oscillators. These models can be used to accurately predict the behavior of individual and mutually coupled oscillators under perturbation at a lower cost than full circuit simulations. The approach is illustrated by numerical experiments with realistic designs. Index Terms—Behavioral modeling, circuit simulation, injection locking, phase noise, pulling, voltage-controlled oscillators (VCOs).

I. I NTRODUCTION

T

HE DESIGN of modern radio-frequency (RF) integrated circuits becomes increasingly more complicated due to the fact that more functionality needs to be integrated on a smaller physical area. In the design process, floor planning, i.e., determining the locations for the functional blocks, is one of the most challenging tasks. Modern RF chips for mobile devices, for instance, typically have an FM radio, Bluetooth, and GPS on one chip. These functionalities are implemented with voltage-controlled oscillators (VCOs) that are designed to oscillate at certain different frequencies. In the ideal case, the oscillators operate independently, i.e., they are not perturbed by each other or any signal other than their input signal. Practically speaking, however, the oscillators are influenced by unintended (parasitic) signals coming from other blocks (such as power amplifiers) or from other oscillators, via for instance (unintended) inductive coupling through the substrate. A possibly undesired consequence of the perturbation is that the oscillators lock to a frequency different than designed for, or show pulling, in which case the oscillators are perturbed from their free-running orbit without locking.

Manuscript received January 19, 2009; revised March 23, 2009 and June 18, 2009. This work was supported by EU Marie-Curie Project O-MOORE-NICE! FP6 MTKI-CT-2006-042477. This paper was recommended by Associate Editor H. E. Graeb. D. Harutyunyan is with CASA group, Department of Mathematics and Computer Science, University of Eindhoven, 5600 MB Eindhoven, The Netherlands (e-mail: [email protected]). J. Rommes, J. ter Maten, and W. Schilders are with NXP Semiconductors, 5656 AE Eindhoven, The Netherlands (e-mail: [email protected]; [email protected]; [email protected]). Color versions of one or more of the figures in this paper are available online at http://ieeexplore.ieee.org. Digital Object Identifier 10.1109/TCAD.2009.2026359

The locking effect was first observed by the Dutch scientist Christian Huygens in the 17th century. He observed that pendulums of two nearby clocks hanging on the same wall after some time moved in unison [1] (in other words, they locked to the same frequency). Similar effects occur also for electrical oscillators. When an oscillator is locked to a different frequency, it physically means that the frequency of the oscillator is changed, and as a result, the oscillator operates at the new frequency. In this case, in the spectrum of the oscillator, we will observe a single peak corresponding to the new frequency of the oscillator. Contrary to the locking case, frequency pulling occurs when the interfering frequency source is not strong enough to cause frequency locking (e.g., weak substrate coupling). In this case, in the spectrum of the pulled oscillator, we will observe several sidebands around the carrier frequency of the oscillator. In Section VII, we will discuss several practical examples of locking and pulling effects. Oscillators appear in many physical systems, and interaction between oscillators has been of interest in many applications. Our main motivation comes from the design of RF systems, where oscillators play an important role [1]–[4] in, for instance, high-frequency phase-locked loops. Oscillators are also used in the modeling of circadian rhythm mechanisms, one of the most fundamental physiological processes [5]. Another application area is the simulation of large-scale biochemical processes [6]. Although the use of oscillators is widely spread over several disciplines, their intrinsic nonlinear behavior is similar, and, moreover, the need for fast and accurate simulation of their dynamics is universal. These dynamics include changes in the frequency spectrum of the oscillator due to small noise signals (an effect known as jitter [2]), which may lead to pulling or locking of the oscillator to a different frequency and may cause the oscillator to malfunction. The main difficulty in simulating these effects is that both phase and amplitude dynamics are strongly nonlinear and spread over separated time scales [7]. Hence, accurate simulation requires very small time steps during time integration, resulting in unacceptable simulation times that block the design flow. Even if computationally feasible, transient simulation only gives limited understanding of the causes and mechanisms of the pulling and locking effects. To some extent, one can describe the relation between the locking range of an oscillator and the amplitude of the injected signal (these terms will be explained in more detail in Section II). Adler [8] shows that this relation is linear, but it is now well known that this is only the case for small injection levels and that the modeling fails for higher injection levels [9]. In addition, other linearized modeling techniques [1] suffer,

0278-0070/$26.00 © 2009 IEEE

2

IEEE TRANSACTIONS ON COMPUTER-AIDED DESIGN OF INTEGRATED CIRCUITS AND SYSTEMS

despite their simplicity, from the fact that they cannot model nonlinear effects such as injection locking [9], [10]. In this paper, we use the nonlinear phase macromodel introduced in [2] and further developed and analyzed in [7], [9]–[13]. Contrary to linear macromodels, the nonlinear phase macromodel is able to capture nonlinear effects such as injection locking. Moreover, since the macromodel replaces the original oscillator system by a single scalar equation, simulation times are decreased while the nonlinear oscillator effects can still be studied without loss of accuracy. One of the contributions of this paper is that we show how such macromodels can be used in industrial practice to predict the behavior of inductively coupled oscillators. Returning to our motivation, during floor planning, it is of crucial importance that the blocks are located in such a way that the effects of any perturbing signals are minimized. A practical difficulty here is that transient simulation of the full system is very expensive and usually unfeasible during the early design stages. One way to get insight in the effects of inductive coupling and injected perturbation signals is to apply the phase shift analysis [2]. In this paper, we will explain how this technique can be used to estimate the effects for perturbed individual and coupled oscillators, and how this can be of help during floor planning. We will consider perturbations caused by oscillators and by other components such as balanced/ unbalanced transformers (baluns). This paper is organized as follows. In Section II, we summarize the phase noise theory. A practical oscillator model and an example application are described in Section III. Inductively coupled oscillators are discussed in detail in Section IV. In Section V, we give an overview of existing methods to model injection locking of individual and resistively/capacitively coupled oscillators. In Section VI, we show how the phase noise theory can be used to analyze oscillator-balun coupling. Numerical results are presented in Section VII, and the conclusions are drawn in Section VIII. II. P HASE N OISE A NALYSIS OF O SCILLATOR

Fig. 1. VCO: current of the nonlinear resistor is given by f (v) = S tanh((Gn /S)v(t)).

where b(t) ∈ Rn are perturbations to the free-running oscillator. For small perturbations b(t), it can be shown [2] that the solution of (2) can be approximated by xp (t) = xpss (t + α(t)) + y(t)

(3)

where y(t) is the orbital deviation and α(t) ∈ R is the phase shift, which satisfies the following scalar nonlinear differential equation: α(t) ˙ = V T (t + α(t)) · b(t)

α(0) = 0

(4a) (4b)

where V (t) ∈ Rn is called perturbation projection vector (PPV) of (2). It is a special projection vector of the perturbations and is computed based on Floquet theory [2], [13], [19]. The PPV is a periodic function with the same period as the oscillator and can efficiently be computed directly from the PPS solution (see, for example, [20]). Using this simple and numerically cheap method, one can do many kinds of analysis for oscillators, e.g., injection locking, pulling, a priori estimate of the locking range [2], [9]. For small perturbations, the orbital deviation y(t) can be ignored [2], and the response of the perturbed oscillator is computed by xp (t) = xpss (t + α(t)) .

(5)

A general free-running oscillator can be expressed as an autonomous system of differential (algebraic) equations dq(x) + j(x) = 0 dt x(0) = x(T )

III. LC O SCILLATOR (1a) (1b)

where x(t) ∈ Rn are the state variables, T is the period of the free-running oscillator, which is, in general, unknown, q, j : Rn → Rn are (nonlinear) functions describing the oscillator’s behavior and n is the system size. The solution of (1) is called periodic steady state (PSS) and is denoted by xpss . Although finding the PSS solution can be a challenging task in itself, we will not discuss this in this paper and refer the interested reader to, for example, [11], [14]–[18]. A general oscillator under perturbation can be expressed as a system of differential equations dq(x) + j(x) = b(t) dt

(2)

For many applications, oscillators can be modeled as an LC tank with a nonlinear resistor, as shown in Fig. 1. This circuit is governed by the following differential equations for the unknowns (v, i): ! " Gn dv(t) v(t) + + i(t) + S tanh v(t) = b(t) (6a) C dt R S L

di(t) − v(t) = 0 dt

(6b)

where C, L, and R are the capacitance, inductance, and resistance, respectively. The nodal voltage is denoted by v, and the branch current of the inductor is denoted by i. The voltagecontrolled nonlinear resistor is defined by S and Gn parameters, where S has influence on the oscillation amplitude and Gn is the gain [9].

HARUTYUNYAN et al.: SIMULATION OF COUPLED OSCILLATORS USING NONLINEAR PHASE MACROMODELS

3

Fig. 4. Two inductively coupled LC oscillators.

Fig. 2. Side band level of the voltage response versus the injected current amplitude for different offset frequencies.

Fig. 3. Floor plan with relocation option that was considered after nonlinear phase noise analysis showed an intolerable pulling due to unintended coupling. Additionally, shielding was used to limit coupling effects even further.

Many work [1], [9] has been done for the simulation of this type of oscillators. Here, we will give an example that can be of practical use for designers. During the design process, early insight in the behavior of system components is of crucial importance. In particular, for perturbed oscillators, it is very convenient to have a direct relationship between the injection amplitude and the side band level. For the given RLC circuit with the following parameters L = 930 · 10−12 H, C = 1.145 · 10−12 F, R = 1000 Ω, S = 1/R, Gn = −1.1/R, and injected signal b(t) = Ainj sin(2πf ), we plot the side band level of the voltage response versus the amplitude Ainj of the injected signal for different offset frequencies (see Fig. 2). The results in Fig. 2 can be seen as a simplified representation of Arnol’d tongues [21], which is helpful in engineering practice. We see, for instance, that the oscillator locks to a perturbation signal with an offset of 10 MHz if the corresponding amplitude is larger than ∼10−4 A (when the signal is locked the side band level becomes 0 dB). This information is useful when designing the floor plan of a chip, since it may put additional requirements on the placement (and shielding) of components that generate, or are sensitive to, perturbing signals. As an example, consider the floor plan in Fig. 3. The analysis described above and in Fig. 2 first helped to identify and

quantify the unintended pulling and locking effects due to the coupling of the inductors (note that the potential causes (inductors) of pulling and locking effects first have to be identified; in practice, designers usually have an idea of potential coupling issues, for instance when there are multiple oscillators in a design). The outcome of this analysis indicated that there were unintended pulling effects in the original floorplan, and hence, some components were relocated (and shielded) to reduce unintended pulling effects. Finally, the same macromodels, but with different coupling factors due to the relocation of components, were used to verify the improved floorplan. Although the LC tank model is relatively simple, it can be of high value particularly in the early stages of the design process (schematic level), since it can be used to estimate the effects of perturbation and (unintended) coupling on the behavior of oscillators. As explained before, this may be of help during floor planning. In later stages, one typically validates the design via layout simulations, which can be much more complex due to the inclusion of parasitic elements. In general, one has to deal with larger dynamical systems when parasitics are included, but the phase noise theory still applies. Therefore, in this paper we do not consider extracted parasitics. However, the values for L, C, R, and coupling factors are typically based on measurement data and layout simulations of real designs. IV. M UTUAL I NDUCTIVE C OUPLING Next, we consider the two mutually coupled LC oscillators shown in Fig. 4. The inductive coupling between these two oscillators can be modeled as di2 (t) di1 (t) +M = v1 (t) dt dt di1 (t) di2 (t) L2 +M = v2 (t) dt dt

L1

(7a) (7b)

√ where M = k L1 L2 is the mutual inductance and |k| < 1 is the coupling factor. This makes the matrix !

L1 M

M L2

"

positive definite, which ensures that the problem is well posed. In this section, all the parameters with a subindex refer to the parameters of the oscillator with the same subindex. If we combine the mathematical model (6) of each oscillator with (7),

4

IEEE TRANSACTIONS ON COMPUTER-AIDED DESIGN OF INTEGRATED CIRCUITS AND SYSTEMS

then the two inductively coupled oscillators can be described by the following differential equations: ! " Gn dv1 (t) v1 (t) + v1 (t) = 0 (8a) + i1 (t) + S tanh C1 dt R1 S di1 (t) di2 (t) − v1 (t) = −M dt dt ! " Gn dv2 (t) v2 (t) C2 + v2 (t) = 0 + i2 (t) + S tanh dt R2 S L1

L2

di2 (t) di1 (t) − v2 (t) = −M . dt dt

(8b) (8c)

Fig. 5. Two resistively coupled LC oscillators.

(8d)

For small values of the coupling factor k, the right-hand side of (8b) and (8d) can be considered as a small perturbation to the corresponding oscillator, and we can apply the phase shift theory described in Section II. Then, we obtain the following simple nonlinear equations for the phase shift of each oscillator: ! " 0 (t + α (t)) · (9a) α˙1 (t) = V T di (t) 1 2 1 −M dt ! " 0 (t + α (t)) · (9b) α˙2 (t) = V T 2 1 (t) 2 −M didt where the currents and voltages are evaluated by using (5)

Fig. 6. Two capacitively coupled LC oscillators.

[15] and references therein for more details on time integration of electric circuits. V. R ESISTIVE AND C APACITIVE C OUPLING

[v1 (t), i1 (t)]T = x1pss (t + α1 (t))

(9c)

[v2 (t), i2 (t)]T = x2pss (t + α2 (t)) .

(9d)

For completeness in this section, we describe how the phase noise theory applies to two oscillators coupled by a resistor or a capacitor.

Small parameter variations have also been studied in the literature by Volterra analysis (see, e.g., [22], [23]).

A. Resistive Coupling

A. Time Discretization

Resistive coupling is modeled by connecting two oscillators by a single resistor (see Fig. 5). The current iR0 flowing through the resistor R0 satisfies the following relation:

The system (9) is solved by using implicit backward Euler for the time discretization and the Newton method is applied for the solution of the resulting 2-D nonlinear equations (10a) and (10b), i.e., # m+1 $ α1m+1 = α1m + τ V T + α1m+1 1 t ! " 0 m+1 m · (10a) −M i2 (t τ)−i2 (t ) # m+1 $ α2m+1 = α2m + τ V T + α2m+1 2 t ! " 0 m+1 m · (10b) −M i1 (t τ)−i1 (t ) % # $ &T v1 (tm+1 ), i1 (tm+1 ) = x1pss tm+1 + α1m+1 (10c) % # $ & T v2 (tm+1 ), i2 (tm+1 ) = x2pss tm+1 + α2m+1 (10d) α11

= 0;

α21

= 0;

m = 1, . . .

where τ = tm+1 − tm denotes the time step. For the Newton iterations in (10a) and (10b), we take (α1m , α2m ) as initial guess on the time level (m + 1). This provides very fast convergence (in our applications within around four Newton iterations). See

iR0 =

v1 − v 2 R0

(11)

where R0 is the coupling resistance. Then, the phase macromodel is given by ! " (v1 − v2 )/R0 (t + α (t)) · α˙1 (qt) = V T (12a) 1 1 0 ! " −(v1 − v2 )/R0 T (12b) α˙2 (t) = V 2 (t + α2 (t)) · 0 where the voltages are updated by using (5). More details on resistively coupled oscillators can be found in [7]. B. Capacitive Coupling When two oscillators are coupled via a single capacitor with a capacitance C0 (see Fig. 6), then the current iC0 through the capacitor C0 satisfies iC0 = C0

d(v1 − v2 ) . dt

(13)

HARUTYUNYAN et al.: SIMULATION OF COUPLED OSCILLATORS USING NONLINEAR PHASE MACROMODELS

5

' where Mij = kij Li Lj , i, j = 1, 2, 3, i < j is the mutual inductance and kij is the coupling factor. The parameters of the nonlinear resistor are S = 1/R1 and Gn = −1.1/R1 and the current injection in the primary balun is denoted by I(t). For small coupling factors, we can consider M12 (di2 (t)/ dt) + M13 (di3 (t)/dt) in (15b) as a small perturbation to the oscillator. Then, similar to (9), we can apply the phase shift macromodel to (15a) and (15b). The reduced model corresponding to (15a) and (15b) is ! " dα(t) 0 = V T (t + α(t)) · . 2 (t) 3 (t) −M12 didt − M13 didt dt (16) The balun is described by a linear circuit (15c)–(15f) which can be written in a more compact form Fig. 7.

Oscillator coupled with a balun.

In this case, the phase macromodel is given by ! d(v −v ) " C0 1dt 2 α˙1 (t) = V T 1 (t + α1 (t)) · 0 ! " −C0 d(v1dt−v2 ) (t + α (t)) · α˙2 (t) = V T 2 2 0

E

(14b)

VI. O SCILLATOR C OUPLING W ITH B ALUN In this section, we analyze inductive coupling effects between an oscillator and a balun. A balun is an electrical transformer that can transform balanced signals to unbalanced signals and vice versa, and they are typically used to change impedance (applications in (RF) radio). The (unintended) coupling between an oscillator and a balun typically occurs on chips that integrate several oscillators for, for instance, FM radio, Bluetooth, and GPS, and hence it is important to understand possible coupling effects during the design. In Fig. 7, a schematic view is given of an oscillator which is coupled with a balun via mutual inductors. The following mathematical model is used for oscillator and balun coupling (see Fig. 7): ! " Gn dv1 (t) v1 (t) + v1 (t) = 0 (15a) + i1 (t) + S tanh C1 dt R1 S di1 (t) di2 (t) di3 (t) + M12 + M13 − v1 (t) = 0 (15b) dt dt dt dv2 (t) v2 (t) + C2 + i2 (t) + I(t) = 0 (15c) dt R2

L2

L3

(17)

where (14a)

where the voltages are updated by using (5). Time discretization of (12) and (14) is done according to (10).

L1

di1 (t) dx(t) + Ax(t) + B +C =0 dt dt

di2 (t) di1 (t) d3 (t) + M12 + M23 − v2 (t) = 0 (15d) dt dt dt dv3 (t) v3 (t) + C3 + i3 (t) = 0 (15e) dt R3

di3 (t) di1 (t) di2 (t) + M13 + M23 − v3 (t) = 0 (15f) dt dt dt



 C2 0 0 0 L2 0 M23   0 E=  0 0 C3 0 L3 0 M23 0   0 0 1/R2 1 0 0  −1 0 A=  0 0 1/R3 0 0 0 −1 0 B T = ( 0 M12 0 M13 ) C T = ( I(t) 0 0 0 ) xT = ( v2 (t) i2 (t) v3 (t) i3 (t) ) .

(18a)

(18b) (18c) (18d) (18e)

With these notations, (16) and (17) can be written in the following form: ! " dα(t) 0 T = V (t + α(t)) · (19) −B T dx(t) dt dt E

di1 (t) dx(t) + Ax(t) + B +C =0 dt dt

(20)

where i1 (t) is computed by using (5). This system can be solved by using a finite difference method. VII. N UMERICAL E XPERIMENTS It is known that a perturbed oscillator either locks to the injected signal or is pulled, in which case side band frequencies all fall on one side of the injected signal (see, e.g., [9]). We will see that contrary to the single oscillator case, where side band frequencies all fall on one side of the injected signal, for (weakly) coupled oscillators, a double-sided spectrum is formed. In Section VII-A–C, we consider two LC oscillators with different kinds of coupling and injection. The inductance and resistance in both oscillators are L1 = L2 = 0.64 nH and R1 = R2 = 50 Ω, respectively. The first oscillator is designed to have a free-running frequency f1 = 4.8 GHz with capacitance C1 = 1/(4L1 π 2 f12 ) = 1.7178 pF. Then, the inductor current in

6

IEEE TRANSACTIONS ON COMPUTER-AIDED DESIGN OF INTEGRATED CIRCUITS AND SYSTEMS

Fig. 8. Inductive coupling. Comparison of the output spectrum of the first oscillator obtained by the phase macromodel and by the full simulation for a different coupling factor k. (a) k = 0.0005. (b) k = 0.001. (c) k = 0.005. (d) k = 0.01.

the first oscillator is A1 = 0.0303 A, and the capacitor voltage is V1 = 0.5844 V. In a similar way, the second oscillator is designed to have a free-running frequency f2 = 4.6 GHz with the inductor current A2 = 0.0316 A, and the capacitor voltage V2 = 0.5844 V. For both oscillators, we choose Si = 1/Ri , Gn = −1.1/Ri with i = 1, 2. In Section VII-D, we describe experiments for an oscillator coupled to a balun. The values for L, C, R, and (mutual) coupling factors are based on measurement data and layout simulations of real designs. In all the numerical experiments, the simulations are run until Tfinal = 6 · 10−7 s with the fixed time step τ = 10−11 . Simulation results, with the phase shift macromodel, are compared with simulations of the full circuit using the CHORAL [24], [25] one-step time integration algorithm, hereafter referred to as full simulation. All experiments have been carried out in Matlab 7.3. We would like to remark that in all experiments simulations with the macromodels were typically ten times faster than the full circuit simulations. In all experiments, for a given oscillator or balun, we use the response of the nodal voltage to plot the spectrum (spectrum composed of discrete harmonics) of the signal. A. Inductively Coupled Oscillators Numerical simulation results of two inductively coupled oscillators (see Fig. 4) for different coupling factors k are shown in Fig. 8, where the frequency is plotted versus power spectral density (PSD1 ). In Fig. 8, we present results for the first oscillator. Similar results are obtained for the second 1 Matlab

code for plotting the PSD is given in [3].

Fig. 9. Inductive coupling. Phase shift α1 (t) of the first oscillator with k = 0.001.

oscillator around its own carrier frequency. For small values of the coupling factor, we observe a very good approximation with the full simulation results. As the coupling factor grows, small deviations in the frequency occur [see Fig. 8(d)]. Because of the mutual pulling effects between the two oscillators, a double-sided spectrum is formed around each oscillator carrier frequency. The additional sidebands are equally spaced by the frequency difference of the two oscillators. The phase shift α1 (t) of the first oscillator for a certain time interval is shown in Fig. 9. We note that it has a sinusoidal behavior. For a single oscillator under perturbation, a completely different behavior is observed: In locked condition, the phase shift changes linearly, whereas in the unlocked case, the phase shift has a nonlinear behavior different than a sinusoidal (see, for example, [26]).

HARUTYUNYAN et al.: SIMULATION OF COUPLED OSCILLATORS USING NONLINEAR PHASE MACROMODELS

7

Fig. 10. Capacitive coupling. Comparison of the output spectrum of the first oscillator obtained by the phase macromodel and by the full simulation for a different coupling factor k. (a) k = 0.0005. (b) k = 0.001. (c) k = 0.005. (d) k = 0.01.

indication that the frequency of the first oscillator is changed and is locked to a new frequency, which is equal to (1 + a)f1 . The change of the frequency can be explained as follows: As noted in [6], capacitive coupling may change the free-running frequency because this kind of coupling changes the equivalent tank capacitance. From a mathematical point of view, it can be explained in the following way. For the capacitively coupled oscillators, the governing equations can be written as dv1 (t) v1 (t) + + i1 (t) dt R ! " Gn dv2 (t) v1 (t) = C0 + S tanh S dt

(C1 + C0 )

Fig. 11. Capacitive coupling. Phase shift of the first oscillator with k = 0.001.

B. Capacitively Coupled Oscillators The coupling capacitance in Fig. 6 is chosen to be C0 = k · Cmean , where Cmean = (C1 + C2 )/2 = 1.794 · 10−12 , and we call k the capacitive coupling factor. Simulation results for the first oscillator for different capacitive coupling factors k are shown in Fig. 10 (similar results are obtained for the second oscillator around its own carrier frequency). For a larger coupling factor k = 0.01, the phase shift macromodel shows small deviations from the full simulation results Fig. 10(d). The phase shift α1 (t) of the first oscillator and a zoomed section for some interval are shown in Fig. 11. In a long run, the phase shift seems to change linearly with a slope of a = −0.00052179. The linear change in the phase shift is a clear

di1 (t) − v1 (t) = 0 dt dv2 (t) v2 (t) + + i2 (t) (C2 + C0 ) dt R ! " Gn dv1 (t) v2 (t) = C0 + S tanh S dt L1

L2

di2 (t) − v2 (t) = 0. dt

(21a) (21b)

(21c) (21d)

This shows that the capacitance in each oscillator is changed by C0 and that the new frequency of each oscillator is f˜i =

2π

'

1 L1 (Ci + C0 )

,

i = 1, 2.

In the zoomed figure within Fig. 11, we note that the phase shift is not exactly linear but that there are small wiggles. By numerical experiments, it can be shown that these small wiggles

8

IEEE TRANSACTIONS ON COMPUTER-AIDED DESIGN OF INTEGRATED CIRCUITS AND SYSTEMS

Fig. 12. Inductive coupling with injection and k = 0.001. (Top) Phase shift. (Bottom) Comparison of the output spectrum obtained by the phase macromodel and by the full simulation with a small current injection. (a) Oscillator 1. (b) Oscillator 2. (c) Oscillator 1. (d) Oscillator 2.

are caused by a small sinusoidal contribution to the linear part of the phase shift. As in case of mutually coupled inductors, the small sinusoidal contributions are caused by mutual pulling of the oscillators [right-hand side terms in (21a) and (21c)].

C. Inductively Coupled Oscillators Under Injection As a next example, let us consider two inductively coupled oscillators where in one of the oscillators an injected current is applied. Let us consider the case where a sinusoidal current of the form I(t) = Ainj sin (2π(f1 − foff )t)

(22)

α˙1 (t) = V

(t + α1 (t)) ·

!

−I(t) 2 (t) −M didt

"

.

Finally, consider an oscillator coupled to a balun, as shown in Fig. 7, with the following parameters values: Oscillator

Primary Balun

Secondary Balun

L1 = 0.64 nH C1 = 1.71 pF R1 = 50 Ω

L2 = 1.10 nH C2 = 4.00 pF R2 = 40 Ω

L3 = 3.60 nH C3 = 1.22 pF R2 = 60 Ω

The coupling factors in (15) are chosen to be k12 = 10−3

k13 = 5.96 ∗ 10−3

k23 = 9.33 ∗ 10−3 .

(24)

The injected current in the primary balun is of the form

is injected in the first oscillator. Then, (9a) is modified to T 1

D. Oscillator Coupled to a Balun

(23)

For a small current injection with Ainj = 10 µA and an offset frequency foff = 20 MHz, the spectra of both oscillators and the phase shift with coupling factor k = 0.001 are shown in Fig. 12. It is clear from Fig. 12(a) and (b) that the phase shift of both oscillators does not change linearly, which implies that the oscillators are not in the steady state. As a result in Fig. 12(c) and (d), we observe spectral widening in the spectra of both oscillators. We note that the phase macromodel simulations are good approximations of the full simulation results.

I(t) = Ainj sin (2π(f0 − foff )t)

(25)

where f0 = 4.8 GHz is the oscillator’s free-running frequency and foff = 20 MHz is the offset frequency. Results of numerical experiments done with the phase macromodel and the full simulations are shown in Fig. 13. We note that for a small current injection (Ainj = 10−4 −10−2 A), both the oscillator and the balun are pulled by each other. When the injected current is not strong (Ainj = 10−4 A), the oscillator is pulled slightly, and in the spectrum of the oscillator [Fig. 13(a)], we observe a spectral widening with two spikes around −60 dB (weak “disturbance” of the oscillator). By gradually increasing the injected current, the oscillator becomes

HARUTYUNYAN et al.: SIMULATION OF COUPLED OSCILLATORS USING NONLINEAR PHASE MACROMODELS

9

Fig. 13. Comparison of the output spectrum of the oscillator coupled to a balun obtained by the phase macromodel and by the full simulation for an increasing injected current amplitude Ainj and an offset frequency foff = 20 MHz. (a) Oscillator. (b) Primary balun. (c) Oscillator. (d) Primary balun. (e) Oscillator. (f) Primary balun. (g) Oscillator. (h) Primary balun.

more disturbed, and in the spectrum, we observe widening with higher side band levels, cf., Fig. 13(c)–(f). When the injected current is strong enough (with Ainj = 10−1 A) to lock the

oscillator to the frequency of the injected signal, we observe a single spike at the new frequency. Similar results are also obtained for the secondary balun.

10

IEEE TRANSACTIONS ON COMPUTER-AIDED DESIGN OF INTEGRATED CIRCUITS AND SYSTEMS

VIII. C ONCLUSION In this paper, we have shown how nonlinear phase macromodels can be used to accurately predict the behavior of individual or mutually coupled VCOs under perturbation, and how they can be used during the design process. Several types of coupling (resistive, capacitive, and inductive) have been described and for small perturbations, the nonlinear phase macromodels produce results with accuracy comparable to full circuit simulations, but at much lower computational costs. Furthermore, we have studied the (unintended) coupling between an oscillator and a balun, a case which typically arises during design and floor planning of RF circuits. ACKNOWLEDGMENT The authors would like to thank J.-P. Frambach (STNWireless) for many helpful discussions about voltage- and digitally controlled oscillators and M. Hanssen (NXP Semiconductors) for providing us measurement data for the inductor and capacitors. The authors would also like to thank M. Striebel of the University of Chemnitz for the CHORAL implementation. R EFERENCES [1] B. Razavi, “A study of injection locking and pulling in oscillators,” IEEE J. Solid-State Circuits, vol. 39, no. 9, pp. 1415–1424, Sep. 2004. [2] A. Demir, A. Mehrotra, and J. Roychowdhury, “Phase noise in oscillators: A unifying theory and numerical methods for characterization,” IEEE Trans. Circuits Syst. I, Fundam. Theory Appl., vol. 47, no. 5, pp. 655–674, May 2000. [3] M. E. Heidari and A. A. Abidi, “Behavioral models of frequency pulling in oscillators,” in Proc. IEEE Int. Behavioral Modeling Simul. Workshop, 2007, pp. 100–104. [4] A. Banai and F. Farzaneh, “Locked and unlocked behaviour of mutually coupled microwave oscillators,” Proc. Inst. Elect. Eng.—Microw. Antennas Propag., vol. 147, no. 1, pp. 13–18, Feb. 2000. [5] S. Agarwal and J. Roychowdhury, “Efficient multiscale simulation of circadian rhythms using automated phase macromodelling techniques,” in Proc. Pacific Symp. Biocomputing, 2008, vol. 13, pp. 402–413. [6] X. Lai and J. Roychowdhury, “Fast simulation of large networks of nanotechnological and biochemical oscillators for investigating selforganization phenomena,” in Proc. IEEE Asia South-Pacific Des. Autom. Conf., 2006, pp. 273–278. [7] X. Lai and J. Roychowdhury, “Fast and accurate simulation of coupled oscillators using nonlinear phase macromodels,” in Proc. IEEE MTT-S Int. Microw. Symp. Dig., 2005, pp. 871–874. [8] R. Adler, “A study of locking phenomena in oscillators,” Proc. IRE Waves Electrons, vol. 34, no. 6, pp. 351–357, Jun. 1946. [9] X. Lai and J. Roychowdhury, “Capturing oscillator injection locking via nonlinear phase-domain macromodels,” IEEE Trans. Microw. Theory Tech., vol. 52, no. 9, pp. 2251–2261, Sep. 2004. [10] Y. Wan, X. Lai, and J. Roychowdhury, “Understanding injection locking in negative resistance LC oscillators intuitively using nonlinear feedback analysis,” in Proc. IEEE Custom Integr. Circuits Conf., 2005, pp. 729–732. [11] M. Günther, U. Feldmann, and J. ter Maten, “Modelling and discretization of circuit problems,” in Handbook of Numerical Analysis, ser. Handb. Numer. Anal., XIII, vol. XIII. Amsterdam, The Netherlands: North Holland, 2005, pp. 523–659. [12] P. Maffezzoni and D. D’Amore, “Evaluating pulling effects in oscillators due to small-signal injection,” IEEE Trans. Comput.-Aided Design Integr. Circuits Syst., vol. 28, no. 1, pp. 22–31, Jan. 2009. [13] P. Maffezzoni, “Unified computation of parameter sensitivity and signalinjection sensitivity in nonlinear oscillators,” IEEE Trans. Comput.-Aided Design Integr. Circuits Syst., vol. 27, no. 5, pp. 781–790, May 2008. [14] S. H. J. M. Houben, “Circuits in motion: The numerical simulation of electrical oscillators,” Ph.D. dissertation, Technische Universiteit Eindhoven, Eindhoven, The Netherlands, 2003.

[15] W. H. A. Schilders and E. J. W. ter Maten, Eds., Numerical Methods in Electromagnetics, ser. Handbook of Numerical Analysis, vol. 13. Amsterdam, The Netherlands: Elsevier, 2005. [16] T. A. M. Kevenaar, “Periodic steady state analysis using shooting and wave-form-Newton,” Int. J. Circuit Theory Appl., vol. 22, no. 1, pp. 51–60, Jan./Feb. 1994. [17] K. Kundert, J. White, and A. Sangiovanni-Vincentelli, “An envelopefollowing method for the efficient transient simulation of switching power and filter circuits,” in Proc. IEEE ICCAD. Dig. Tech. Papers, Nov. 1988, pp. 446–449. [18] A. Semlyen and A. Medina, “Computation of the periodic steady state in systems with nonlinear components using a hybrid time and frequency domain method,” IEEE Trans. Power Syst., vol. 10, no. 3, pp. 1498–1504, Aug. 1995. [19] H. G. Brachtendorf, Theorie und Analyse von autonomen und quasiperiodisch angeregten elektrischen Netzwerken. Bremen, Germany: Habilitationsschrift, Universität Bremen, 2001. [20] A. Demir, D. Long, and J. Roychowdhury, “Computing phase noise eigenfunctions directly from steady-state Jacobian matrices,” in Proc. IEEE/ACM ICCAD, 2000, pp. 283–288. [21] P. L. Boyland, “Bifurcations of circle maps: Arnol’d tongues, bistability and rotation intervals,” Commun. Math. Phys., vol. 106, no. 3, pp. 353–381, Sep. 1986. [22] J. P. Aikio, “Frequency domain model fitting and Volterra analysis implemented on top of harmonic balance simulation,” Ph.D. dissertation, Faculty Technol. Univ. Oulu, Oulu, Finland, 2007. [23] J. P. Aikio, M. Makitalo, and T. Rahkonen, “Harmonic load-pull technique based on Volterra analysis,” in Proc. Eur. Microw. Conf., 2009, submitted for publication. [24] M. Günther, “Simulating digital circuits numerically—A charge-oriented ROW approach,” Numer. Math., vol. 79, no. 2, pp. 203–212, Apr. 1998. [25] P. Rentrop, M. Günther, M. Hoschek, and U. Feldmann, “CHORAL—A charge-oriented algorithm for the numerical integration of electrical circuits,” in Mathematics—Key Technology for the Future. Berlin, Germany: Springer-Verlag, 2003, pp. 429–438. [26] X. Lai and J. Roychowdhury, “Automated oscillator macromodelling techniques for capturing amplitude variations and injection locking,” in Proc. IEEE/ACM ICCAD, Nov. 2004, pp. 687–694.

Davit Harutyunyan received the M.Sc. degree in mathematics from the University of Groningen, Groningen, The Netherlands, in 2002 and the Ph.D. degree in mathematics from the University of Twente, Enschede, The Netherlands, in 2007. He is currently a Postdoc Researcher with the University of Eindhoven, Eindhoven, The Netherlands.

Joost Rommes received the M.Sc. degree in computational science, the M.Sc. degree in computer science, and the Ph.D. degree in mathematics from Utrecht University, Utrecht, The Netherlands, in 2002, 2003, and 2007, respectively. He is currently a Researcher with NXP Semiconductors, Eindhoven, The Netherlands, and works on model reduction.

HARUTYUNYAN et al.: SIMULATION OF COUPLED OSCILLATORS USING NONLINEAR PHASE MACROMODELS

Jan ter Maten received the M.Sc. degree and Ph.D. degree in mathematics from Utrecht University, Utrecht, The Netherlands, in 1976 and 1984, respectively. In 1983, he joined Philips Electronics and later became a Researcher with Philips Research. Research topics covered simulation of high-frequency problems, multirate time integration, model order reduction, statistics, optimization, and sensitivity analysis used in simulation techniques in circuit simulation. Since 2006, he has been with NXP Semiconductors, Eindhoven, The Netherlands. He is active within scientific computing in electrical engineering and energetically stimulates Europe-wide research cooperation in the aforementioned research areas.

11

Wil Schilders received the M.Sc. degree (cum laude) in mathematics from Radboud University, Nijmegen, in 1978 and the Ph.D. degree from Trinity College Dublin, Dublin, Ireland, in 1980. In 1980, he joined Philips Electronics and has done extensive work on simulation methods for semiconductor devices, electronic circuits, electromagnetics, and other problems related to the electronics industry. In 1999, he was a Part-Time Professor of numerical mathematics for industry with Eindhoven University of Technology. In 2006, the semiconductor activities of Philips were transferred to the new company NXP Semiconductors, Eindhoven, The Netherlands, where he is currently a Mathematics Group Leader. Dr. Schilders served as Vice-President of the European Consortium for Mathematics in Industry and Editor-in-Chief for the Dutch Journal for Mathematicians.