# Nonlinear Circuit Analysis via Perturbation Methods and Hardware Prototyping

Authors: K. Odame; P. E. Hasler
Journal: VLSI Design
 (2010)

---

## Abstract

Nonlinear signal processing is necessary in many emerging applications where form factor and power are at a premium. In order to make such complex computation feasible under these constraints, it is necessary to implement the signal processors as analog circuits. Since analog circuit design is largely based on a linear systems perspective, new tools are being introduced to circuit designers that allow them to understand and exploit circuit nonlinearity for useful processing. This paper discusses two such tools, which represent nonlinear circuit behavior in a graphical way, making it easy to develop a qualitative appreciation for the circuits under study.


---

## Body

## 1. Introduction

Portable and implantable, always-on electronics stand to benefit from analog signal processing, when only low levels of precision are necessary [1, 2]. To achieve sophisticated signal processing with low power and area overhead, an analog processor can exploit the fundamental nonlinear dynamics that are found in devices and simple circuits. So, the circuit designer must depart from the traditional linear systems paradigm, and learn to analyze and understand circuits from a nonlinear dynamical systems theory perspective.In order to make nonlinear circuit design relevant to an engineer, it must be taught intuitively enough to foster creativity, yet rigorously enough to be of practical benefit. The two most popular tools for studying nonlinear circuits are harmonic balance and Volterra series. Within certain limitations, they are rigorous, but not necessarily intuitive. Since harmonic balance is simulation-based, it can be used to predict nonlinear behavior without ever requiring a deep understanding [3] of the circuit. Volterra series is an analytical tool that quickly leads to high entropy [4] mathematical expressions, from which the causative physical phenomena are hard to discern, much less purposefully manipulate.To bridge the gap between rigor and intuition, we can use visual representation techniques. If the appropriate visualization is formed from rigorous definitions of a circuit's dynamics, then the human vision system, with its pattern recognition ability, will perceive the circuit's qualitative behavior [5].We will present a filter block diagram for analyzing harmonic distortion that is derived from perturbation analysis. Unlike Volterra series kernels, our filter block diagram does not include multidimensional Fourier transforms, and so is accessible to an introductory-level engineering audience. In the second part of this paper, we will discuss the creation and use of phase plane plots of nonlinear circuits. We will describe how to rapidly create the phase plane plots with a reconfigurable hardware platform, instead of with a numerical simulator.This paper is an expansion of the work presented in [6].


## 2. Regular Perturbation

Whenever designers want to get an analytical handle on the sources and causes of distortion, the most commonly-used tool is Volterra series analysis. If a problem is tractable using Volterra series, then it can also be solved with perturbation theory, which will yield asymptotically-identical results [7].There are certain problems for which Volterra series are ill-suited—multiple-time-scale behavior and multiple steady states, for instance [8]—that can be solved with perturbation theory. Despite the power of perturbation theory, it is still a relatively obscure concept in discussions about nonlinearity and distortion in analog circuits.We therefore find it worthwhile to give a basic treatment of regular perturbation—the simplest perturbation method—as applied to distortion analysis of first-order analog circuits. In addition, a filter block diagram representation of the circuit will naturally evolve from our analysis, making it visually clear how the distortion terms are manifested, and how well-known tenets of low-distortion design, such as feedback, come about.Consider the initial value problem(1)ẋ=f(t,x,ϵ);x(t0)=x0(ϵ),
where ϵ is a small perturbation parameter such that ϵ=0 yields an analytically-soluble equation. If f is sufficiently smooth (the specific smoothness requirements of f are discussed in [9]), then the problem has a unique solution x(t,ϵ). As the solution for ϵ≠0 may not be analytical, it can be approximated as a power series in ϵ to an accuracy of O(ϵn+1). That is, we can write the solution as (2)x(t,ϵ)=∑i=0n(xi(t)ϵi)︸x̂(t,ϵ)+O(ϵn+1),
where x̂(t,ϵ) is the approximate solution. To conduct regular perturbation, we apply the substitution x(t,ϵ)≈x̂(t,ϵ) to  (1). The resulting system is then solved by equating like powers of ϵ. The following sections will illuminate this idea.


## 3. The Basic First-Order Circuit

Most common first-order analog ciruits (simple amplifiers, buffers, switches, etc.) are of the form depicted in  Figure1. The governing equation is(3)ẋ=g(u-x)+h(x),
where u is the a.c. input signal, x is the a.c. output signal and g(·) and h(·) are nonlinear functions. The dependence of the system on the output, other than through feedback to the input, is modeled by h(x). In practice, h(x) is typically some nonideality such as finite output resistance.Figure 1
General block diagram form of a first-order circuit. The primary processing block isg(·), which is a nonlinear function of the input u and of x via feedback. The nonlinearity h(x) models such nonidealities as finite output impedance.In order to apply perturbation analysis to  (3), we begin by assuming that the input signal has a small amplitude. This is expressed as u=ϵv, where ϵ is a small perturbation parameter and v is a suitably-scaled version of the input signal. Note that with the definition of u, (3) is solvable via separation of variables for the special case ϵ=0.With the introduction of the perturbation parameterϵ, we can approximate the solution to  (3) with the power series(4)x(t)≈∑i=1nϵixi(t).
Note the ϵ0 term of  (4) is set to 0. This corresponds to analyzing a circuit about its d.c. bias point, where the d.c. bias point is shifted to the origin. For ease of notation, define z=u-x. The approximation of z is defined similarly to  (4), with z1=v-x1 and zi=-xi,for all i>1.Ifϵ is sufficiently small, then the functions g(z) and h(x) can be approximated by their truncated Taylor series as(5)g(z)≈g1z+gn-1zn-1+gnzn,h(x)≈h1x+hn-1xn-1+hnxn.
Functions g and h are assumed to be dominantly (n-1)th-order nonlinearities, with gi=g(i)(0)/i! and hi=h(i)(0)/i!.  Equation (5) assumes g(0)=h(0)=0, which, again, corresponds to analyzing a circuit about its d.c. bias point.Substituting  (4) and (5) into  (3) and collecting powers of ϵ, we get the following set of first-order linear equations(6)ẋ1+(g1-h1)x1=g1v,⋮ẋk+(g1-h1)xk=0∀k<n-1⋮ẋn-1+(g1-h1)xn-1=gn-1z1n-1+hn-1x1n-1ẋn+(g1-h1)xn=gnz1n-ngn-1z1n-1x2+hnx1n-nhn-1x1n-1x2.Theϵ1 equation is the linearized portion of  (3) with input v. Taking the Laplace transform of this equation, we write(7)X1(s)=g1H(s)V(s),
where H(s)=1/(s+g1-h1).Theϵk equations (k<(n-1)) are filters with 0 input. As such, the steady state solutions of these equations is 0.


## 4. Harmonic Distortion Terms

The inputs of theϵn-1 equation are terms of z1n-1 and x1n-1. To understand the implications of these terms to harmonic distortion, assume a single-tone input, v=cos(ωt). This elicits the signals(8)x1=g1|H(jω)|cos(ωt+ϕ(jω)),z1=|1-g1H(jω)|cos(ωt+ϕz1(jw)),=|(s-h1)H(jω)︸Hz1(jω)|cos(ωt+ϕz1(jw)),
Here we have defined Hz1(s)=(1-g1H(s)). The phases ϕ(s) and ϕz1(s) are the arguments of H(s) and Hz1(s), respectively. The signals x1 and z1 are single tones of frequency ω as well, since they are merely linearly-filtered versions of v.Raisingz1 and x1 each to the (n-1)th power produces harmonics as follows. If (n-1) is odd(even), then odd(even) harmonics up to the (n-1)th harmonic are generated. The amplitude of the mω frequency term in x1n-1 is(9)(n-1)!g1((n+m-1)/2)!((n-m-1)/2)!2n-2|H(jω)|,
while that of the mω frequency term in z1n-1 is(10)(n-1)!((n+m-1)/2)!((n-m-1)/2)!2n-2|Hz1(jω)|.
After filtering in the ϵn-1 equation, the amplitudes of these terms will be, respectively,(11)(n-1)!hn-1g1((n+m-1)/2)!((n-m-1)/2)!2n-2|H(jω)||H(jmω)|,(12)(n-1)!gn-1((n+m-1)/2)!((n-m-1)/2)!2n-2|Hz1(jω)||H(jmω)|.Analogous to that of theϵn-1 equation, the input to the ϵn equation has terms in z1n and x1n. In general, the x2 terms are identically zero, except for the special case n=3.


## 5. Feedback and Distortion

We now make some observations about the harmonic distortion results that were discussed in the previous section.In theϵn-1 equation, the amplitude of the mth harmonic that the z1n-1 term contributes is given by  (12). We plot this amplitude expression, along with that of  (11), as a function of frequency in  Figure 2 for the third-order harmonic generated by a dominantly-third order nonlinearity. That is, n=4 and m=3. Also, we chose h1=1, h3=1/3, g1=G, g3=G/3, where G was varied from 10 to 1000.Magnitude-frequency plots of the third harmonic. The “gain”,G of the g(z) function is varied from 10 to 1000. This causes the band-pass shape of the z13-contributed harmonic to shift to the right, while that contributed by x13 falls in magnitude.
(a)(b)Notice from the figure that ifg1≫h1, then, for a given frequency, the amplitude of the z1n-1-contributed harmonic is greatly reduced. In fact, if we ensure gi≫hiforalli, then the harmonic contribution of the x1n-1 terms is negligible. This would mean that the distortion is effectively due only to z1, whose associated harmonics are band-pass filtered. This in turn means that the distortion can be kept small if the circuit is operated well below the corner frequency.These two notions—that frequency and feedback gain can be sacrificed for higher linearity—conform with the traditional rules-of-thumb for low-distortion design.


## 6. Illustrative Examples

### 6.1. Source Follower Amplifier

According to KCL, the circuit equation of the source follower amplifier in  Figure3(a) is(13)CdVout(t)dt=F(Vin,Vout)-Ibias,
where the function F is defined as(14)F(Vin,Vout)=K2(κVin(t)-Vout(t)-Vth)2,
if M1 is in above-threshold saturation, and(15)F(Vin,Vout)=Ioe(κVin(t)-Vout(t))/UT,
if it is in subthreshold saturation. The parameter K depends on transistor dimensions and doping and Vth is the threshold voltage. Also, κ, Io, and UT have their usual meanings from the EKV MOSFET model [10].Source follower amplifier. (a) Circuit schematic. (b) Block diagram representation of source follower output. The fundamental harmonic is a low-pass filtered version of the input. The second order terms are generated by high-pass filtering the input, squaring and then low pass filtering. The total output is a power series ofϵ terms.
(a)(b)Note thatIbias=F(Vg,VS), where Vg and VS are the d.c. bias-points of the gate and source of M1, respectively. Let us define a characteristic voltage, Vc, as(16)Vc={(κVg-VS-Vth)2,abovethresholdUT,subthreshold.
Now,  (13) can be nondimensionalized [9] by making the substitutions(17)τ=Ibias(CVc)·t;u=κνinVc;x=νoutVc,
where νin and νout are the a.c. portions of Vin and Vout. This gives the state-space equation of the source follower as(18)dxdτ=u-x+(u-x)24,
for above threshold, and(19)dxdτ=u-x+(u-x)22,
for the truncated Taylor expansion in subthreshold. The point is that, regardless of region of operation of M1, the nonlinear equation that describes the source follower has the same functional form. Relating the source follower equations to  (3), we have g(z)~z+z2 and h(x)=0. As such, we expect the harmonic distortion terms to have a band-pass-like dependence on frequency. To show this, we will apply regular perturbation to  (18).First, defineu=ϵv, where the small parameter ϵ is a scaled version of the input amplitude. That is, ϵ=Ain/Vc. Also, taking x=ϵx1+ϵ2x2 and z=u-x and equating like powers of ϵ up to ϵ2, we have
(20)ϵ1:ẋ1=v-x1,(21)ϵ2:ẋ2=z124-x2,
as depicted in  Figure 3(b). Assume a pure-tone input, v=cos(ωt).  Equation (20) is the linear portion of the amplifier.  Equation (21) is a linear filter with input z12/4. The squaring produces a second-harmonic term as well as a d.c. offset. In addition, since z1=v-x1, the second harmonic generated by the squaring is high-pass filtered. The overall effect is that x2 is a band-pass filtered version of a second harmonic of v.  Figure 4 is a plot of experimental data that corroborates our analysis. There is error in the second harmonic measurement due to the small amplitudes involved.Figure 4
Magnitude-frequency response of source follower. Analytical prediction is in bold, and experimental data is plotted as “x”s and “o”s. The fundamental harmonic is a low-pass filtered version of the input. The second harmonic has a bandpass shape, as predicted by perturbation analysis.


### 6.2. Unity Gain Buffer

Consider the unity-gain buffer depicted in  Figure5(a). It is formed by placing an operational transconductance amplifier (OTA) in negative feedback. If we operate the OTA above threshold, the describing equation is(22)CdVoutdt=κβIbias(Vin-Vout)1-κβ(Vin-Vout)24Ibias,
while it is
(23)CdVoutdt=Ibiastanh(κ(Vin-Vout)2UT),
for subthreshold operation. Notice that we have ignored the output conductance term, which is considered very small for OTAs.Unity gain buffer. (a) Circuit schematic. (b) Block diagram representation of output. The fundamental harmonic is a low-pass filtered version of the input. The third-order terms are generated by high-pass filtering the input, cubing and then low pass filtering. The total output is a power series ofϵ terms.
(a)(b)We can define a characteristic voltage,Vc, as(24)Vc={2UTκ,subthresholdIbiasκβ,abovethreshold.
Then, with the following definitions(25)τ=Ibias(CVc)·t;u=νinVc;x=νoutVc,
the nondimensional form of the unity-gain buffer's describing equations (taken to the first few Taylor series terms) is(26)dxdτ={(u-x)-(u-x)34,abovethreshold(u-x)-(u-x)33,subthreshold.
Again, the functional form of the equations is identical, regardless of region of operation.To calculate distortion terms, assumeu=ϵv is a pure-tone signal and proceed as usual. For subthreshold, the separated equations of ϵ are
(27)ϵ1:ẋ1=v-x1,(28)ϵ2:ẋ2=0-x2,(29)ϵ3:ẋ3=z133-x3.
These equations are depicted in the block diagram of  Figure 5(b).  Equation (27) is the linear portion of the amplifier.  Equation (28) is a linear filter with 0 input; it contributes no harmonics at steady state.  Equation (29) is a linear filter with input z13/3. The cubing produces a third-harmonic term as well as a fundamental-frequency term (this fundamental-frequency term will cause gain compression, which is not discussed in this paper). Since z1=v-x1, the overall effect is that x3 is a band-pass filtered version of a third harmonic of v, as shown in  Figure 6.Figure 6
Magnitude-frequency response of unity-gain buffer. Analytical prediction is in bold, and experimental data is plotted as “x”s and “o”s. The fundamental harmonic is a low-pass filtered version of the input. The third harmonic has a bandpass shape, as predicted by perturbation analysis.


### 6.3. Note on Above-Threshold versus Subthreshold Operation

The harmonic behavior of a circuit is similar for above- and subthreshold operation. In absolute numbers, however, above threshold operation yields less distortion. This is because the parameterϵ=Ain/Vc is much smaller for above threshold than for subthreshold. Since the harmonics are multiplied by ϵi, the smaller ϵ seen in above threshold operation translates to lower distortion.


## 6.1. Source Follower Amplifier

According to KCL, the circuit equation of the source follower amplifier in  Figure3(a) is(13)CdVout(t)dt=F(Vin,Vout)-Ibias,
where the function F is defined as(14)F(Vin,Vout)=K2(κVin(t)-Vout(t)-Vth)2,
if M1 is in above-threshold saturation, and(15)F(Vin,Vout)=Ioe(κVin(t)-Vout(t))/UT,
if it is in subthreshold saturation. The parameter K depends on transistor dimensions and doping and Vth is the threshold voltage. Also, κ, Io, and UT have their usual meanings from the EKV MOSFET model [10].Source follower amplifier. (a) Circuit schematic. (b) Block diagram representation of source follower output. The fundamental harmonic is a low-pass filtered version of the input. The second order terms are generated by high-pass filtering the input, squaring and then low pass filtering. The total output is a power series ofϵ terms.
(a)(b)Note thatIbias=F(Vg,VS), where Vg and VS are the d.c. bias-points of the gate and source of M1, respectively. Let us define a characteristic voltage, Vc, as(16)Vc={(κVg-VS-Vth)2,abovethresholdUT,subthreshold.
Now,  (13) can be nondimensionalized [9] by making the substitutions(17)τ=Ibias(CVc)·t;u=κνinVc;x=νoutVc,
where νin and νout are the a.c. portions of Vin and Vout. This gives the state-space equation of the source follower as(18)dxdτ=u-x+(u-x)24,
for above threshold, and(19)dxdτ=u-x+(u-x)22,
for the truncated Taylor expansion in subthreshold. The point is that, regardless of region of operation of M1, the nonlinear equation that describes the source follower has the same functional form. Relating the source follower equations to  (3), we have g(z)~z+z2 and h(x)=0. As such, we expect the harmonic distortion terms to have a band-pass-like dependence on frequency. To show this, we will apply regular perturbation to  (18).First, defineu=ϵv, where the small parameter ϵ is a scaled version of the input amplitude. That is, ϵ=Ain/Vc. Also, taking x=ϵx1+ϵ2x2 and z=u-x and equating like powers of ϵ up to ϵ2, we have
(20)ϵ1:ẋ1=v-x1,(21)ϵ2:ẋ2=z124-x2,
as depicted in  Figure 3(b). Assume a pure-tone input, v=cos(ωt).  Equation (20) is the linear portion of the amplifier.  Equation (21) is a linear filter with input z12/4. The squaring produces a second-harmonic term as well as a d.c. offset. In addition, since z1=v-x1, the second harmonic generated by the squaring is high-pass filtered. The overall effect is that x2 is a band-pass filtered version of a second harmonic of v.  Figure 4 is a plot of experimental data that corroborates our analysis. There is error in the second harmonic measurement due to the small amplitudes involved.Figure 4
Magnitude-frequency response of source follower. Analytical prediction is in bold, and experimental data is plotted as “x”s and “o”s. The fundamental harmonic is a low-pass filtered version of the input. The second harmonic has a bandpass shape, as predicted by perturbation analysis.


## 6.2. Unity Gain Buffer

Consider the unity-gain buffer depicted in  Figure5(a). It is formed by placing an operational transconductance amplifier (OTA) in negative feedback. If we operate the OTA above threshold, the describing equation is(22)CdVoutdt=κβIbias(Vin-Vout)1-κβ(Vin-Vout)24Ibias,
while it is
(23)CdVoutdt=Ibiastanh(κ(Vin-Vout)2UT),
for subthreshold operation. Notice that we have ignored the output conductance term, which is considered very small for OTAs.Unity gain buffer. (a) Circuit schematic. (b) Block diagram representation of output. The fundamental harmonic is a low-pass filtered version of the input. The third-order terms are generated by high-pass filtering the input, cubing and then low pass filtering. The total output is a power series ofϵ terms.
(a)(b)We can define a characteristic voltage,Vc, as(24)Vc={2UTκ,subthresholdIbiasκβ,abovethreshold.
Then, with the following definitions(25)τ=Ibias(CVc)·t;u=νinVc;x=νoutVc,
the nondimensional form of the unity-gain buffer's describing equations (taken to the first few Taylor series terms) is(26)dxdτ={(u-x)-(u-x)34,abovethreshold(u-x)-(u-x)33,subthreshold.
Again, the functional form of the equations is identical, regardless of region of operation.To calculate distortion terms, assumeu=ϵv is a pure-tone signal and proceed as usual. For subthreshold, the separated equations of ϵ are
(27)ϵ1:ẋ1=v-x1,(28)ϵ2:ẋ2=0-x2,(29)ϵ3:ẋ3=z133-x3.
These equations are depicted in the block diagram of  Figure 5(b).  Equation (27) is the linear portion of the amplifier.  Equation (28) is a linear filter with 0 input; it contributes no harmonics at steady state.  Equation (29) is a linear filter with input z13/3. The cubing produces a third-harmonic term as well as a fundamental-frequency term (this fundamental-frequency term will cause gain compression, which is not discussed in this paper). Since z1=v-x1, the overall effect is that x3 is a band-pass filtered version of a third harmonic of v, as shown in  Figure 6.Figure 6
Magnitude-frequency response of unity-gain buffer. Analytical prediction is in bold, and experimental data is plotted as “x”s and “o”s. The fundamental harmonic is a low-pass filtered version of the input. The third harmonic has a bandpass shape, as predicted by perturbation analysis.


## 6.3. Note on Above-Threshold versus Subthreshold Operation

The harmonic behavior of a circuit is similar for above- and subthreshold operation. In absolute numbers, however, above threshold operation yields less distortion. This is because the parameterϵ=Ain/Vc is much smaller for above threshold than for subthreshold. Since the harmonics are multiplied by ϵi, the smaller ϵ seen in above threshold operation translates to lower distortion.


## 7. Field Programmable Analog Array

We have developed a field programmable analog array (FPAA) that can be configured to synthesis and analyze a vast variety of circuits [11].  Figure 7 depicts a simple, second-order filter compiled on our FPAA. More complex circuit configurations are possible, and would involve a larger number of the over 400 components in the FPAA. In this part of the paper we demonstrate the utility and versatility of the FPAA in analyzing the dynamics of a number of fundamental circuit blocks.Figure 7
Compiled Circuit on an FPAA. Component terminals can be connected or disconnected using the switch matrix. In this illustration, connections are indicated by solder dots. The gain of each amplifier is programable over a continuous range of values over a few decades.


## 8. One-Dimensional Systems

### 8.1. Simple Current Mirror

Consider the simple current mirror depicted in  Figure8(a). From Kirchhoff's Current Law (KCL), it obeys the following differential equation:(30)CdVgdt=Ib-f(Vg),
where f(Vg) is the drain current of transistor M1. Assuming M1 and M2 are identical and are both saturated, we have f(Vg)=Iout, which gives(31)CdVgdt=Ib-Iout.Simple current mirror. (a) Circuit that was compiled onto the FPAA. (b) Measured trajectories for different initial conditions. (c) Vector field derived from trajectory measurements. The origin is an unstable equilibrium point, while5 is stable.
(a)(b)(c)For subthreshold operation in saturation, the current through transistors M1 and M2 is [10]  (32)f(Vg)=Iout=Ioe(κVg-VS)/UT,
where Io is a pre-exponential constant dependent on the transistor's size and on doping concentrations. Also, κ is the body-effect coefficient and UT is the thermal voltage. VS is the source voltage, which, for this case, is zero. Setting VS=0 and taking the derivative of  (32) with respect to time, we get
(33)dIoutdt=∂∂Vg(IoeκVg/UT)dVgdt,(34)=κUTIoutdVgdt,
which allows us to rewrite  (31) as(35)CUTκIbdIoutdt=Iout(1-IoutIb),τdIoutdt=Iout(1-IoutIb).
The time constant is identified as τ=(CUT)/(κIb).Equation (35) happens to be the logistic equation, a simple model of population dynamics. It can be solved exactly either by separation of variables followed by partial fractions, or by solving it as Bernoulli's equation. The solution is 		(36)Iout(t)=Ibet/τet/τ-1+Ib/Iout0,
where Iout0 is the initial value of Iout. We are lucky to have an exact solution to  (35), given that it is a nonlinear differential equation. Even so, it is difficult to discern much useful information about Iout's qualitative behavior from  (36). For instance, it is not clear how the behavior of Iout might change with different initial conditions. To answer questions of this sort, it is helpful to do geometric analysis on the system's corresponding vector field.Since the simple current mirror is a one-dimensional system, its vector field is represented as a flow on a line. The direction and speed of the flow are dictated by the right hand side (RHS) of  (35). It is a quadratic, as shown in  Figure 10. The Iout-intercepts are 0 and Ib. There is a maximum at Iout=Ib/2. The vector field is depicted as the arrows on the Iout axis. For positive values of dIout/dt, Iout is increasing, meaning the arrows point to the right. For negative values of dIout/dt, Iout is decreasing, meaning the arrows point to the left. When dIout/dt=0, there is no change in Iout and the circuit is said to be at equilibrium.The vector field provides clear, qualitative information about the behavior ofIout. There are two equilibrium points, namely Iout=0 and Iout=Ib. Note that the vector field flows away from Iout=0. This equilibrium point is unstable, since the system will not recover from slight disturbances away from it. The vector field flows towards Iout=Ib, implying that this is a stable equilibrium point. If the system is initially at Iout=Ib and then experiences a small disturbance, it will tend back to the Iout=Ib point.The vector field in  Figure10 also gives information about the acceleration of Iout as it approaches the Iout=Ib equilibrium point. For 0<Iout<Ib/2, the rate of change of Iout increases until it reaches a peak at Iout=Ib/2. Between Ib/2 and Ib, the system decelerates until the rate of change of Iout eventually becomes zero. For Iout>Ib, the rate of change of Iout steadily decreases until Iout=Ib. It is interesting to note that, for Iout<Ib, the rate of change of Iout is limited to a maximum of Ib/(4τ).The geometric analysis predictions can be checked against experimental measurements of a current mirror that was compiled onto an FPAA.  Figure8(b) depicts various trajectories, or solutions, of the system of  (35) for different initial conditions. Notice that trajectories that start at values lower than Iout=Ib/2 have a sigmoidal shape, with the point of inflection corresponding to the maximum rate of change of current dIout/dt=Ib/(4τ). The parabolic shape of dIout/dt can be extracted from these trajectories, and it is shown in  Figure 8(c).


### 8.2. Simple Peak Detector

Assuming subthreshold operation, the KCL equation for the source follower amplifier of  Figure9(a) is the following.(37)CdVoutdt=Ioe(κVin-Vout)/UT-Ib.
Note that(38)ddteVout/UT=eVout/UTUTdVoutdt,
in which case, the solution to  (37) is(39)Vout=κVin+UTlog(IoIb-(IoIb-e(Vout0-κVin)/UT)e-t/τ),
where τ=CUT/Ib and Vout0 is the initial value of Vout.Source Follower Amplifier Acting as a Simple Peak Detector. (a) Circuit that was compiled onto the FPAA. (b) Measured step responses. (c) Vector field derived from step response measurements. The point(0,0.7) is a stable fixed point.
(a)(b)(c)Figure 10
The vector field of the logistic equation is represented as a flow on theIout axis. For positive values of dIout/dt, Iout is increasing and the flow is to the right. For negative values of dIout/dt, Iout is decreasing and the flow is to the left.The time that it takes forVout to be within 10% of its final value is
(40)t10=τlog|Io/Ib-e(Vout0-κVin)/UTIo/Ib-e0.1κVin/UT|.For a large positive step input,e(Vout0-κVin)≈0, and  (40) is approximately(41)t10=t10+≈τlog|IoIo-Ibe0.1κVin/UT|.For a large negative step input,e(Vout0-κVin)≫Io/Ib, and  (40) becomes(42)t10=t10-≈τlog|Ibe(Vout0-κVin)Io-Ibe0.1κVin/UT|=t10++τ(Vout0UT-κVinUT)log(IbIo).
Equations (41) and (42) indicate that the response of the peak detector is slower for a negative input step than it is for a positive input step. We surmise that if the input is continuously varying at a rate faster than 1/(t10-), then the output will be a reasonable representation of the input's peak values. Explaining the peak detector's behavior with  (41) and (42) is rigorous, but depends on having to manipulate the expression of  (39).One way of avoiding the math is to employ intuitive descriptions of the charging action of the active device (i.e., the transistor) versus the discharging action of the current source [12]. A more rigorous approach is to apply nonlinear geometric analysis to the problem. Consider the plot of dVout/dt versus Vout shown in  Figure 9(c). We constructed it from a number of step response measurements (Figure 9(b)) that we took after compiling the source-follower amplifier onto the FPAA. A large negative input step corresponds to an initial value of Vout0≫Vin. The rate of growth of Vout is bounded by Ib/C. For a large positive input step, however, Vout0≪Vin, and the maximum rate at which Vout approaches Vin can be much greater than Ib/C. The maximum rate of approach in this case is limited only by the initial value, Vout0. As such, there is an asymmetry in the speed of the circuit's response to up-going versus down-going movements on the input. The effect of this asymmetry is that Vout tracks increasing Vin and not decreasing Vin, which is the behavior of a peak detector.


## 8.1. Simple Current Mirror

Consider the simple current mirror depicted in  Figure8(a). From Kirchhoff's Current Law (KCL), it obeys the following differential equation:(30)CdVgdt=Ib-f(Vg),
where f(Vg) is the drain current of transistor M1. Assuming M1 and M2 are identical and are both saturated, we have f(Vg)=Iout, which gives(31)CdVgdt=Ib-Iout.Simple current mirror. (a) Circuit that was compiled onto the FPAA. (b) Measured trajectories for different initial conditions. (c) Vector field derived from trajectory measurements. The origin is an unstable equilibrium point, while5 is stable.
(a)(b)(c)For subthreshold operation in saturation, the current through transistors M1 and M2 is [10]  (32)f(Vg)=Iout=Ioe(κVg-VS)/UT,
where Io is a pre-exponential constant dependent on the transistor's size and on doping concentrations. Also, κ is the body-effect coefficient and UT is the thermal voltage. VS is the source voltage, which, for this case, is zero. Setting VS=0 and taking the derivative of  (32) with respect to time, we get
(33)dIoutdt=∂∂Vg(IoeκVg/UT)dVgdt,(34)=κUTIoutdVgdt,
which allows us to rewrite  (31) as(35)CUTκIbdIoutdt=Iout(1-IoutIb),τdIoutdt=Iout(1-IoutIb).
The time constant is identified as τ=(CUT)/(κIb).Equation (35) happens to be the logistic equation, a simple model of population dynamics. It can be solved exactly either by separation of variables followed by partial fractions, or by solving it as Bernoulli's equation. The solution is 		(36)Iout(t)=Ibet/τet/τ-1+Ib/Iout0,
where Iout0 is the initial value of Iout. We are lucky to have an exact solution to  (35), given that it is a nonlinear differential equation. Even so, it is difficult to discern much useful information about Iout's qualitative behavior from  (36). For instance, it is not clear how the behavior of Iout might change with different initial conditions. To answer questions of this sort, it is helpful to do geometric analysis on the system's corresponding vector field.Since the simple current mirror is a one-dimensional system, its vector field is represented as a flow on a line. The direction and speed of the flow are dictated by the right hand side (RHS) of  (35). It is a quadratic, as shown in  Figure 10. The Iout-intercepts are 0 and Ib. There is a maximum at Iout=Ib/2. The vector field is depicted as the arrows on the Iout axis. For positive values of dIout/dt, Iout is increasing, meaning the arrows point to the right. For negative values of dIout/dt, Iout is decreasing, meaning the arrows point to the left. When dIout/dt=0, there is no change in Iout and the circuit is said to be at equilibrium.The vector field provides clear, qualitative information about the behavior ofIout. There are two equilibrium points, namely Iout=0 and Iout=Ib. Note that the vector field flows away from Iout=0. This equilibrium point is unstable, since the system will not recover from slight disturbances away from it. The vector field flows towards Iout=Ib, implying that this is a stable equilibrium point. If the system is initially at Iout=Ib and then experiences a small disturbance, it will tend back to the Iout=Ib point.The vector field in  Figure10 also gives information about the acceleration of Iout as it approaches the Iout=Ib equilibrium point. For 0<Iout<Ib/2, the rate of change of Iout increases until it reaches a peak at Iout=Ib/2. Between Ib/2 and Ib, the system decelerates until the rate of change of Iout eventually becomes zero. For Iout>Ib, the rate of change of Iout steadily decreases until Iout=Ib. It is interesting to note that, for Iout<Ib, the rate of change of Iout is limited to a maximum of Ib/(4τ).The geometric analysis predictions can be checked against experimental measurements of a current mirror that was compiled onto an FPAA.  Figure8(b) depicts various trajectories, or solutions, of the system of  (35) for different initial conditions. Notice that trajectories that start at values lower than Iout=Ib/2 have a sigmoidal shape, with the point of inflection corresponding to the maximum rate of change of current dIout/dt=Ib/(4τ). The parabolic shape of dIout/dt can be extracted from these trajectories, and it is shown in  Figure 8(c).


## 8.2. Simple Peak Detector

Assuming subthreshold operation, the KCL equation for the source follower amplifier of  Figure9(a) is the following.(37)CdVoutdt=Ioe(κVin-Vout)/UT-Ib.
Note that(38)ddteVout/UT=eVout/UTUTdVoutdt,
in which case, the solution to  (37) is(39)Vout=κVin+UTlog(IoIb-(IoIb-e(Vout0-κVin)/UT)e-t/τ),
where τ=CUT/Ib and Vout0 is the initial value of Vout.Source Follower Amplifier Acting as a Simple Peak Detector. (a) Circuit that was compiled onto the FPAA. (b) Measured step responses. (c) Vector field derived from step response measurements. The point(0,0.7) is a stable fixed point.
(a)(b)(c)Figure 10
The vector field of the logistic equation is represented as a flow on theIout axis. For positive values of dIout/dt, Iout is increasing and the flow is to the right. For negative values of dIout/dt, Iout is decreasing and the flow is to the left.The time that it takes forVout to be within 10% of its final value is
(40)t10=τlog|Io/Ib-e(Vout0-κVin)/UTIo/Ib-e0.1κVin/UT|.For a large positive step input,e(Vout0-κVin)≈0, and  (40) is approximately(41)t10=t10+≈τlog|IoIo-Ibe0.1κVin/UT|.For a large negative step input,e(Vout0-κVin)≫Io/Ib, and  (40) becomes(42)t10=t10-≈τlog|Ibe(Vout0-κVin)Io-Ibe0.1κVin/UT|=t10++τ(Vout0UT-κVinUT)log(IbIo).
Equations (41) and (42) indicate that the response of the peak detector is slower for a negative input step than it is for a positive input step. We surmise that if the input is continuously varying at a rate faster than 1/(t10-), then the output will be a reasonable representation of the input's peak values. Explaining the peak detector's behavior with  (41) and (42) is rigorous, but depends on having to manipulate the expression of  (39).One way of avoiding the math is to employ intuitive descriptions of the charging action of the active device (i.e., the transistor) versus the discharging action of the current source [12]. A more rigorous approach is to apply nonlinear geometric analysis to the problem. Consider the plot of dVout/dt versus Vout shown in  Figure 9(c). We constructed it from a number of step response measurements (Figure 9(b)) that we took after compiling the source-follower amplifier onto the FPAA. A large negative input step corresponds to an initial value of Vout0≫Vin. The rate of growth of Vout is bounded by Ib/C. For a large positive input step, however, Vout0≪Vin, and the maximum rate at which Vout approaches Vin can be much greater than Ib/C. The maximum rate of approach in this case is limited only by the initial value, Vout0. As such, there is an asymmetry in the speed of the circuit's response to up-going versus down-going movements on the input. The effect of this asymmetry is that Vout tracks increasing Vin and not decreasing Vin, which is the behavior of a peak detector.


## 9. A Two-Dimensional System

Figure11 depicts Lyon and Mead's classic second-order section [13]. It is a Gm-C filter with two poles that can be placed anywhere on the real/imaginary plane. We begin our analysis by writing down the governing equations for the circuit, assuming that the OTAs are based on subthreshold MOS transistor differential pairs: (43)C1dV2dt=I2ktanh(κ(V1-V2)2UT),C1dV1dt=I1tanh(κ(Vin-V1)2UT)-I3tanh(κ(V2-V1)2UT),
where I1,2,3 are the bias currents of the OTAs. Also, k is the ratio of the C2 to C1.Figure 11
Second-Order Section. Varying the bias currents of the various amplifiers leads to interesting dynamics.If we define(44)x=κ(V1-Vin)2UT,y=κ(V2-V1)2UT,
then (43) become(45)2UTC1κdxdt=-I1tanh(x)-I3tanh(y),2UTC1κdydt=I1tanh(x)+(I3-I2k)tanh(y).
Further defining(46)I1=Ibias,I2=gkIbias,I3=2rIbias,t=τ·2UTC1κIbias,
where g≥0, we get the following dimensionless equation(47)dxdτ=-tanh(x)-2rtanh(y),dydτ=tanh(x)+(2r-g)tanh(y).


### 9.1. Small Signal Analysis

We can linearize  (47) by replacing the RHS with its Jacobian, giving(48)[ẋẏ]≈[-1-2r12r-g][xy].
The origin is a fixed point. In fact, it is a unique fixed point, since (from  (47)) tanh(x)=-2rtanh(y)⇒x=y=0. The origin is stable for(49)r<1+g2,
and unstable otherwise. It is a spiral for(50)1+g2-g<r<1+g2+g,
and a node otherwise.


### 9.2. Large Signal Analysis

For certain values ofr and g, the nonlinearities of the second-order section causes it to suffer instability. In this region of parameter space, linear analysis accurately predicts that the circuit is small signal stable, but completely fails to recognize that instability would occur for large signals. Mead addresses this issue in [14], but we will present a somewhat more thorough treatment of the problem, using phase-plane analysis and experimental verification with the FPAA.For very large values ofx and y, the tanh functions get saturated, and can each be approximated with a signum function.  Equation (47) becomes(51)dxdτ=-sgn(x)-2rsgn(y),dydτ=sgn(x)+(2r-g)sgn(y).
Figure 12 shows the phase plane (x versus y) that corresponds to  (51). The depicted motion is valid only if 1/2<r<1. The gradient in the first and third quadrants is(52)a=2r-1+g1+2r,
and that in the second and fourth quadrants is(53)b=2r-1-g1-2r.
Observe that, with an initial condition of (1,0),  (51) predicts that the positive x-axis will again be intercepted at (a2/b2,0). If a2/b2>1, then x and y will grow without bound. Stated in terms of the r and g variables, there is large signal instability if(54)r>g+g2+44.
Our analysis of the second-order section can readily be verified experimentally. We compiled the filter on the FPAA, as shown in  Figure 7. The bias currents of all three OTAs are user-programmable, and varying them corresponds to varying the values of r and g. The FPAA thus allows us to explore the parameter space of the filter, and to observe changes in its qualitative behavior. It can effectively be used for bifurcation analysis.Figure 12
Sketch of SOS Phase Plane for Large Signals.Figure13 shows the filter's phase plane plots for various values of r, with g kept fixed. Just as we predicted, there is a unique fixed point, which is initially stable, and gradually changes from a node to a spiral (Figures 13(a) and 13(c)). While linear analysis would predict these three responses as damped, slightly underdamped, and very underdamped, it fails to recognize the possibility of the fourth response, which is large-signal instability. In the fourth panel, r meets the criterion derived from nonlinear analysis,  (54), and we observe oscillation. Further analysis and exploration of parameter space reveals that this second-order section is capable of low-distortion sinusoidal oscillation [15]. Such functionality is valuable in communication systems.SOS Experimental Phase Plane Results for Various Values ofr.
(a)(b)(c)(d)


## 9.1. Small Signal Analysis

We can linearize  (47) by replacing the RHS with its Jacobian, giving(48)[ẋẏ]≈[-1-2r12r-g][xy].
The origin is a fixed point. In fact, it is a unique fixed point, since (from  (47)) tanh(x)=-2rtanh(y)⇒x=y=0. The origin is stable for(49)r<1+g2,
and unstable otherwise. It is a spiral for(50)1+g2-g<r<1+g2+g,
and a node otherwise.


## 9.2. Large Signal Analysis

For certain values ofr and g, the nonlinearities of the second-order section causes it to suffer instability. In this region of parameter space, linear analysis accurately predicts that the circuit is small signal stable, but completely fails to recognize that instability would occur for large signals. Mead addresses this issue in [14], but we will present a somewhat more thorough treatment of the problem, using phase-plane analysis and experimental verification with the FPAA.For very large values ofx and y, the tanh functions get saturated, and can each be approximated with a signum function.  Equation (47) becomes(51)dxdτ=-sgn(x)-2rsgn(y),dydτ=sgn(x)+(2r-g)sgn(y).
Figure 12 shows the phase plane (x versus y) that corresponds to  (51). The depicted motion is valid only if 1/2<r<1. The gradient in the first and third quadrants is(52)a=2r-1+g1+2r,
and that in the second and fourth quadrants is(53)b=2r-1-g1-2r.
Observe that, with an initial condition of (1,0),  (51) predicts that the positive x-axis will again be intercepted at (a2/b2,0). If a2/b2>1, then x and y will grow without bound. Stated in terms of the r and g variables, there is large signal instability if(54)r>g+g2+44.
Our analysis of the second-order section can readily be verified experimentally. We compiled the filter on the FPAA, as shown in  Figure 7. The bias currents of all three OTAs are user-programmable, and varying them corresponds to varying the values of r and g. The FPAA thus allows us to explore the parameter space of the filter, and to observe changes in its qualitative behavior. It can effectively be used for bifurcation analysis.Figure 12
Sketch of SOS Phase Plane for Large Signals.Figure13 shows the filter's phase plane plots for various values of r, with g kept fixed. Just as we predicted, there is a unique fixed point, which is initially stable, and gradually changes from a node to a spiral (Figures 13(a) and 13(c)). While linear analysis would predict these three responses as damped, slightly underdamped, and very underdamped, it fails to recognize the possibility of the fourth response, which is large-signal instability. In the fourth panel, r meets the criterion derived from nonlinear analysis,  (54), and we observe oscillation. Further analysis and exploration of parameter space reveals that this second-order section is capable of low-distortion sinusoidal oscillation [15]. Such functionality is valuable in communication systems.SOS Experimental Phase Plane Results for Various Values ofr.
(a)(b)(c)(d)


## 10. Conclusion

In this paper, we have introduced visual and graphical techniques for analyzing nonlinear circuit dynamics. Our approach to studying harmonic distortion yields information about the various processing flows that are responsible for each harmonic term. The FPAA was used to rapidly create phase plane plots, which concisely encapsulate the nonlinear dynamics of the circuit under study. We have provided various examples of our techniques and have compared our predictions to experimentally-measured data.



---

*Source: 687498-2010-03-18.xml*