This preview shows page 1. Sign up to view the full content.
Unformatted text preview: MAE107 Introduction to Modeling and Analysis of Dynamic Systems
Prof. R.T. M’Closkey, Mechanical and Aerospace Engineering, UCLA
Subject outline:
1. Intro to ﬁrst order ODEs: ode, ivp, solutions, convolution, transfer function, frequency response, stability
2. Intro to second order ODEs: same subjects
3. Intro to nth order ODEs
4. Poles, zeros, coprime assumptions, block diagram manipulation with transfer functions vs. ODEs
5. Frequency response asymptotic approximations with examples.
6. The unit impulse function (various realizations), unit step function, impulse response and step response, relations between the two, deﬁnitions of 0− and 0+
7. Superposition=convolution, properties of convolution
8. Fourier series and response to periodic inputs; energy spectral density; identiﬁcation with periodic inputs;
discretetime Fourier series for periodic sequences; aliasing
9. Fourier transform; relation between frequency response and FT of impulse response
10. Correlation functions and their transforms for signals with ﬁnite energy; use in identiﬁcation; extension to power
signals including power spectral density functions (also give block diagram interpretation)
11. Bilateral Laplace, ROC, numerical examples
12. Unilateral Laplace and application to solving IVPs
13. Statespace representation theory. Do not distribute these notes 1 Prof. R.T. M’Closkey, UCLA First Order Systems Do not distribute these notes 2 Prof. R.T. M’Closkey, UCLA The Starting Point: First Order Linear Linear Differential Equations The ﬁrst order linear ODE is given by: d
(y (t)) = ay (t) + bu(t)
dt (1) where y = dependent variable that is determined by “solving" the ODE
t = time, the independent variable, speciﬁed on some interval I
u = the “forcing" or “nonhomogeneous" term that is known or speciﬁed for t ∈ I
a, b = coefﬁcients that are assumed to be constant in these notes.
This ODE is linear because y and its time derivative appear linearly in the expression. Another shorthand way of
writing this expression is: y = ay + bu.
˙
In terms of using this model to predict the response of a physical system, the following initial value problem is often
of interest:
Initial Value Problem. The time interval of interest is I = [t0 , t1 ] and an initial condition, denoted y0 , is speciﬁed
at time t0 (in other words, y (t0 ) = y0 is speciﬁed) and the continuous forcing function u is known for all t ∈ I , then,
ﬁnd y for all t ∈ I .
Notes: • the “starting time" t0 is often 0
• t1 may be ∞
• this scenario assumes the system modeled by this ODE is causal (future inputs do not affect the current value
of the dependent variable) Do not distribute these notes 3 Prof. R.T. M’Closkey, UCLA Block diagrams are a useful way to pictorially represent relationships between systems and signals. To derive a
block diagram for (7) note that the unique solution of the IVP also satisﬁes
t y (t) = y0 + (ay (τ ) + bu(τ ))dτ.
t0 A block with an integral symbol is used to denote that the output of the block is the time integral of the signal going
into the block plus a suitable constant. Note that an integrator can “hold" a nonzero output even if the input is
speciﬁed to be zero so the integrator constant is set to the initial condition y0 . This constant is represented in the
block diagram shown below,
summing junction u “integrator” block y
˙ b signal line y0 y a
initial output of
integrator block gain blocks An implicit assumption in block diagrams is that the output of a block (integrator, gain, summer, etc) only depends
on the input(s) to the block. In other words, if the input(s) to a block is some signal, then the block output depends
only on this signal and not what the block output is connected to. Thus, the block is able to provided the necessary
“force" or “current" to drive any block downstream. If this is not the case, then different tools must be used to
analyze the behavior of the connected blocks. This situation is called “loading", and block diagram analysis of
systems in which loading is an issue will give erroneous results. Do not distribute these notes 4 Prof. R.T. M’Closkey, UCLA First Order Examples Lowpass ﬁlter circuit example:
ÿ + R1
ÿ
+
Vin
−
C1 Vout
−
ÿ
ÿ Vin (t) = R1 i(t) + Vout (t)
˙
i(t) = C1 Vout (t) ˙
=⇒ Vout (t) + 1
1
Vout (t) =
Vin (t).
R1 C1
R1 C1 First order linear ODE relating Vin and Vout : 1
1
˙
Vout (t) +
Vin (t)
Vout (t) = −
R1 C1
R1 C1
a (2) b Why is this called a “low pass" ﬁlter? Let Vin = cos(ωt), where ω is the frequency in radians per second. Stable
linear ODE’s driven with a sinusoidal “forcing" function will posses a sinusoidal response after some initial transient.
Thus, we can assume the sinusoidal response of Vout to be Vout (t) = α cos(ωt) + β sin(ωt),
where α and β are to be determined (they are functions of ω ). Substituting this expression into (2) yields ω (−α sin(ωt) + β cos(ωt)) = − 1
1
(α cos(ωt) + β sin(ωt)) +
cos(ωt).
R1 C1
R1 C1 Suppose ω is chosen such that ω <<
Do not distribute these notes 1
,
R1 C1
5 Prof. R.T. M’Closkey, UCLA then the terms on the lefthand side are much smaller than the terms on the righthand side so that we ﬁnd α ≈ 1 and β ≈ 0.
In other words, Vout (t) ≈ Vin (t).
On the other hand, if ω >> 1
,
R1 C1 then α ≈ 0 and β ≈ 1
<< 1,
R1 C1 ω 1
so Vout (t) ≈ R1 C1 ω sin(ωt). Thus, the amplitude of Vout is much smaller than the amplitude of Vin . Furthermore,
there is also a phase difference between the input and output sinusoids: the output lags the input by about one
quarter of a cycle, which can also be expressed as “−90 degrees." When ω ≈ 1/(R1 C1 ) the relationship between α and β and ω is more complicated. This circuit is called a lowpass
ﬁlter precisely because Vout (t) ≈ Vin (t) when ω << 1/(R1 C1 ) and because the amplitude of Vout is much less than
the amplitude of Vin when ω >> 1/(R1 C1 ). This circuit “passes" sinusoids with low frequency and attenuates sinusoids with high frequency, where low frequency versus high frequency is determined by the parameter 1/(R1 C1 ).
The block diagram representation is Vin 1
RC − Do not distribute these notes 6 Vout 1
RC Prof. R.T. M’Closkey, UCLA Highpass ﬁlter circuit example: C1
ÿ +
Vin
−
ÿ
+
R1 Vout
−
ÿ
ÿ Vin (t) = Vcap (t) + Vout (t) 1
˙
˙
Vout (t) = R1 i(t)
=⇒ Vout (t) +
Vout (t) = Vin (t). R1 C1 ˙
C1 Vcap (t) = i(t) As in the lowpass ﬁlter case we can study this circuit’s response to sinusoidal forcing by letting Vin (t) = cos(ωt)
and setting Vout (t) = α cos(ωt) + β sin(ωt),
Substituting into the ODE yields ω (−α sin(ωt) + β cos(ωt)) + 1
(α cos(ωt) + β sin(ωt)) = −ω sin(ωt)
R1 C1 If ω << 1/(R1 C1 ), then α ≈ 0 and β ≈ −R1 C1 ω << 1 so Vout (t) ≈ −R1 C1 ω sin(ωt). Thus, the output amplitude is smaller than the input amplitude and, furthermore, the
output phase leads the input phase by one quarter of a cycle (“+90 degree”).
If ω >> 1/(R1 C1 ), then α ≈ 1 and β ≈ 0 so Vout (t) ≈ Vin (t) = cos(ωt). Thus, this circuit “passes" high frequency sinusoids and attenuates low frequency
sinusoids, where the distinction between high and low frequency is again the parameter 1/(R1 C1 ).
Do not distribute these notes 7 Prof. R.T. M’Closkey, UCLA The ODE is not in our standard y = ay + bu form, though. A change of dependent variable, however, shows how a
˙
highpass ﬁlter can be constructed from a lowpass ﬁlter in series with a “feedthrough" term. Deﬁne V = Vout − Vin
so ˙
˙
˙
V = Vout − Vin
1
˙
˙
=−
Vout + Vin − Vin
RC
1
=−
(V + Vin )
RC
1
1
V−
Vin .
=−
RC
RC Thus, V is a lowpass ﬁltered version of Vin which is also inverted in polarity. The output voltage Vout is reconstructed from Vout = V + Vin . The highpass ﬁlter is expressed with the following block diagram 1 Vin 1
−
RC ˙
V − Do not distribute these notes 8 V Vout 1
RC Prof. R.T. M’Closkey, UCLA Sliding mass example.
Consider a mass sliding with velocity v along a surface with friction and being pulled by a force f :
viscous
friction, c →v
m f
@
R
@
¡¡¡¡¡¡¡¡¡¡ Newton: mv (t) = f (t) − cv (t)
˙
where c is the viscous friction constant. The input variable is assumed to be the force and the output variable is
assumed to be the mass velocity. The ODE can be arranged into standard form v=−
˙ 1
c
v + f.
m
m Block diagram representation: f v
˙ 1
m
− Do not distribute these notes 9 v
c
m Prof. R.T. M’Closkey, UCLA Water tank example.
A water tank with volume ﬂow rate into the tank Qin (the input variable) and water level h (the dependent variable): Pa ↓ Qin 6 h
P1
Variable deﬁnitions: Qin =
Qout =
h=
V=
Pa =
P1 = −→ Qout volume ﬂow rate into tank (m3 /s)
volume ﬂow rate out of tank (m3 /s)
height of water in tank (m)
volume of water in tank (m3 )
atmospheric pressure (N/m2 )
pressure in water at oriﬁce (N/m2 ) Fundamental relation governing water volume in the tank: d
(V (t)) = Qin (t) − Qout (t).
dt
Model of ﬂow through the oriﬁce: Qout (t) = cA 2(P1 (t) − Pa )
,
ρ where c is the oriﬁce coefﬁcient, A is the oriﬁce area, and ρ is the density of water. The pressure at the oriﬁce is
from the hydrostatic pressure law: P1 (t) = Pa + ρgh(t),
Do not distribute these notes 10 Prof. R.T. M’Closkey, UCLA where g = 9.81m/s2 . Combining the last two expressions yields, 2g Qout (t) = cA
The volume as a function of h is V (t) = h(t). π2
D h(t),
4 for a circular tank of diameter D. Combining these expressions yields a nonlinear differential equation for h π 2˙
D h(t) = Qin (t) − cA
4 2g h(t). (3) The block diagram corresponding to this nonlinear ODE is Qin ˙
h 4
π D2 −
√
4cA 2g
π D2 h
√ · nonlinear operation The nonlinear equation can be linearized about an equilibrium point. An equilibrium occurs when Qin (t) =
˙
Qout (t) = Q0 , where Q0 is a positive constant, so h(t) ≡ 0. In this case, the water height is constant and can
be parameterized as a function of Q0 : h(t) = Q2
0
,
(cA)2 2g for all t. The equilibrium point represents a special solution of the differential equation. The signals in the block diagram at
equilibrium are Do not distribute these notes 11 Prof. R.T. M’Closkey, UCLA ˙
h(t) = 0
Qin (t) = Q0 4
π D2 − Q2
0
h(t) =
(cA)2 2g
√ √
4cA 2g
π D2 · Linearization is the process of studying the behavior of small deviations in h relative to small deviations in Qin
relative to their equilibrium values. In this case we let Qin (t) = Q0 + δQ (t), (4) where δQ (t) represents the perturbation in the volume ﬂow rate into the tank. Similarly, we represent the water
height as the nominal water height at its equilibrium value plus a “small" perturbation in height, δh , Q2
0
h(t) =
+ δh (t).
(cA)2 2g (5) The block diagram showing the perturbation variables δh and δQ is Q0 + δQ 4
π D2 −
√
4cA 2g
π D2 Do not distribute these notes 12 Q2
0
+ δh
(cA)2 2g
√ · Prof. R.T. M’Closkey, UCLA Substituting (4) and (5) into (3) yields π 2˙
D δh (t) = Q0 + δQ (t) − cA
4 2g Q2
0
+ δh (t).
(cA)2 2g (6) No approximations has been carried out at this point, only a “translation in coordinates." The approximation comes
Q
in when the square root term is expanded in a Taylor series about the point (cA)0 2g and then discarding terms in δh
2
with powers greater than one: Q0
Q2
0
√
+ δh (t) =
(cA)2 2g
cA 2g
Q0
√
=
cA 2g 1+ (cA)2 2g
δh (t)
Q2
0 (cA)2 2g
(cA)2 2g 2
1+
δh (t) −
δ (t) + · · ·
2Q2
4Q2 h
0
0 Thus, the linear approximation is Q2
Q0
0
√
+ δh (t) ≈
2 2g
(cA)
cA 2g
which is valid when δh (t) <<
motion: Do not distribute these notes Q0
(cA)2 2g . (cA)2 2g
1+
δh (t) ,
2Q2
0 Substituting this approximation into (6) yields the linearized equation of π 2˙
(cA)2 2g
D δh (t) = δQ (t) −
δh (t).
4
2Q0 13 Prof. R.T. M’Closkey, UCLA The block diagram for the linearized equation follows, with the understanding that the input and output variables
are relative to nominal values determined by equilibrium conditions is δQ 4
π D2 ˙
δh
− δh 4(cA)2 g
π Q0 D2 Do not distribute these notes 14 Prof. R.T. M’Closkey, UCLA Solving First Order Linear ODEs
Given the ﬁrst order, linear, constant coefﬁcient ODE y (t) = ay (t) + bu(t),
˙ (7) where u is deﬁned on the interval I = [t0 , t1 ] and a and b are real constants, we label two classes of functions that
are associated with this ODE:
1. Homogeneous solution. Recall that the homogeneous differential equation is the one in which u(t) ≡ 0: y (t) = ay (t).
˙
A nonzero solution to the homogeneous problem, denoted yh , satisﬁes yh (t) = ayh (t).
˙
Note that homogeneous solutions are not unique: a scalar multiple of yh , say αyh for some constant α, is
another homogeneous solution since d
(αyh (t)) = αyh (t) = αayh (t) = a(αyh (t)).
˙
dt
Furthermore, the weighted sum of two homogeneous solutions is another homogeneous solution: let y1 and
y2 be two homogeneous solutions, i.e. y1 = ay1 and y2 = ay2 , then α1 y1 + α2 y2 , where α1 and α2 are any
˙
˙
constants, is also a homogeneous solution since d
d
d
(α1 y1 (t) + α2 y2 (t)) = (α1 y1 (t)) + (α2 y2 (t))
dt
dt
dt
= α1 y1 (t) + α2 y2 (t)
˙
˙
= α1 ay1 (t) + α2 ay2 (t)
= a (α1 y1 (t) + α2 y2 (t)) .
For ﬁrstorder differential equations, a fundamental solution is a homogeneous solution in which yh (t0 ) = 1.
In fact, for ﬁrstorder constant coefﬁcient differential equations all homogeneous solutions are of the form αeat
for any nonzero constant α.
Do not distribute these notes 15 Prof. R.T. M’Closkey, UCLA 2. Particular solution. For the nonhomogeneous ODE, a particular solution, denoted yp , is any function that
satisﬁes (7) (includes the nonhomogeneous term bu(t)). In other words yp (t) = ayp (t) + bu(t), for all t ∈ I .
˙
Note that particular solutions are not unique: adding a homogeneous solution to a particular solution yields
another particular solution.
We can now state the initial value problem for (7):
Initial Value Problem (IVP). Let an initial condition y0 be speciﬁed at time t0 (in other words, y (t0 ) = y0 is speciﬁed)
and the forcing function u is known for time in the interval t ∈ I (note that the beginning of the time interval coincides
with the time when the initial condition is speciﬁed). Then, ﬁnd y over the same time interval.
Solution of the IVP. The solution of the IVP is unique and may be found by several different methods.
1. Method 1: Ad hoc approach. This method is called ad hoc because you must ﬁrst produce any particular
solution, then, add to it a homogeneous solution with “free" parameter α, i.e., y (t) = yp (t) + αyh (t),
where you adjust the constant α such that y (t0 ) matches the initial condition, i.e. y (t0 ) = yp (t0 ) + αyh (t0 ) = y0 .
The difﬁculty is computing the particular solution. There are certain classes of nonhomogeneous terms,
however, for which it is very easy to compute a particular solution. These will be discussed shortly.
2. Method 2: Integrating Factor. The homogeneous problem y = ay has as its fundamental solution, denoted
˙
yf , the following function yf (t) = ea(t−t0 ) , t ≥ t0 . The integrating factor is a function v such that v (t)yf (t) = 1,
Do not distribute these notes 16 t ≥ t0 .
Prof. R.T. M’Closkey, UCLA It is easily seen that the integrating factor is v (t) = ea(−t+t0 ) , t ≥ t0 . Multiplying the ODE by the integrating factor yields d
(y (t))ea(−t+t0 ) − ay (t)ea(−t+t0 ) = ea(−t+t0 ) bu(t)
dt
d
a(−t+t0 )
)
dt (y (t)e
Integrating both sides: t
t0 d
y (τ )ea(−τ +t0 ) dτ =
dτ =⇒ y (t)ea(−t+t0 ) − y (t0 ) =
=⇒ y (t) = ea(t−t0 ) y0 + t t t ea(−τ +t0 ) bu(τ )dτ t0 ea(−τ +t0 ) bu(τ )dτ t0 ea(t−τ ) bu(τ )dτ t0 This way of expressing the solution is convenient since it separately shows the contributions due to the initial
condition versus the contribution due to the “input" u.
Method 3: Variation of Parameters. Variation of parameters uses the fundamental solution as a “basis" for
generating a particular solution that satisﬁes the nonhomogeneous differential equation: yp (t) = c(t)ea(t−t0 ) , t ≥ t0 , where c is to be determined by substitution into the nonhomogeneous ODE: c(t)ea(t−t0 ) + ac(t)ea(t−t0 ) = ac(t)ea(t−t0 ) + bu(t)
˙
=⇒ c(t) = ea(−t+t0 ) bu(t)
˙
t =⇒ c(t) − c(t0 ) =
Do not distribute these notes ea(−τ +t0 ) bu(τ )dτ. t0 17 Prof. R.T. M’Closkey, UCLA The initial value c(t0 ) can be taken to be zero without loss of generality because if it is not zero then it just
generates a homogeneous solution that will be compensated by the addition of another homogeneous solution
with a free parameter in order to match the initial condition of the IVP. Thus,
t c(t) =
t0 and so ea(−τ +t0 ) bu(τ )dτ, t ∈ [t0 , t1 ], yp (t) = c(t)ea(t−t0 )
t =
t0 ea(t−τ ) bu(τ )dτ, t ∈ [t0 , t1 ]. Since yp (t0 ) = 0, the homogeneous term that we need to add to compute the IVP solution is yh (t) = ea(t−t0 ) x0 .
Thus, adding these particular and homogeneous solutions yields the solution to the IVP, i.e. the solution
satisﬁes the ODE and the initial condition, y (t) = yh (t) + yp (t)
a(t−t0 ) =e t y0 +
t0 (8) ea(t−τ ) bu(τ )dτ, t ∈ [t0 , t1 ]. Note that this is exactly the same expression obtained with the integrating factor method. Do not distribute these notes 18 Prof. R.T. M’Closkey, UCLA Terminology “Zerostate response" vs. “Free response”
The solution to the initial value problem (8) is unique and this particular form is useful because it is decomposed
into the sum of the response to the input u when the initial condition is zero with the response of the system due
to its nonzero initial condition when the input is zero. These “pieces" of the solution have special names as shown
below (t0 = 0 in the following expression):
t at y (t) = e y0 +
“free
response" 0 ea(t−τ ) bu(τ )dτ , t ≥ 0. “zerostate response" The term “zerostate response" comes from the fact that this is the system’s response when u is applied to the
system when it is initially in an equilibrium, or rest, state. There are other ways to parse the solution to the IVP
and in doing so we will be able to deﬁne the terms “forced motion", “natural motion", “steadystate response", and
“transient response." Do not distribute these notes 19 Prof. R.T. M’Closkey, UCLA Linearity and Superposition
Analyzing the solution (8) of the IVP gives us the correct insight into the concept of “superposition". Suppose an
initial condition at t0 is speciﬁed as y0,1 and the input u1 is deﬁned on the interval I = [t0 , t1 ]. Then, the solution to
the IVP, denoted y1 , is y1 (t) = e a(t−t0 ) t y0,1 +
t0 ea(t−τ ) bu1 (τ )dτ, t ∈ [t0 , t1 ]. (9) If another initial condition y0,2 at t0 is given, along with another input u2 on the interval I , then the solution to this
initial value problem, denoted y2 , is y2 (t) = e a(t−t0 ) t y0,2 +
t0 ea(t−τ ) bu2 (τ )dτ, t ∈ [t0 , t1 ]. (10) A fair question to ask is
“to what initial value problem is y1 (t) + y2 (t), t ∈ I , the solution?"
By adding (9) and (10) this question is answered: y1 (t) + y2 (t) = e a(t−t0 ) t y0,1 + e a(t−τ ) a(t−t0 ) bu1 (τ )dτ + e t y0,2 + t0 = ea(t−t0 ) (y0,1 + y0,2 ) +
initial
condition ea(t−τ ) bu2 (τ )dτ t0
t ea(t−τ ) b(u1 (τ ) + u2 (τ ))dτ t0 input Thus, the initial condition y0,1 + y0,2 with input u1 + u2 produces the solution y1 + y2 . In other words, for initial value
problems, it is important to sum inputs and sum initial conditions.
Example.
Consider the following lowpass ﬁlter circuit equations when RC = 0.1, ˙
Vout = −10Vout + 10Vin .
Do not distribute these notes 20 Prof. R.T. M’Closkey, UCLA For purposes of standardizing notation we will use the form y = −10y + 10u.
˙
Let t0 = 0 and y0,1 = 1 volt (this initial conditioning is determined by the charge on the capacitor at t = 0). The
input is u1 (t) = sin(10t), t ≥ 0. The solution to the IVP, denoted y1 , is shown below where it has been partitioned
according to (8).
IVP with u1 and x0,1
2
IVP sol.
solution with IC=x0,1, u=0
IC=0, non zero u
input u 1.5 y1 = −10y1 + 10u1
˙ =⇒
y1 (0) = y0,1 = 1 u1 (t) = sin(10t), t ≥ 0 x1 (volts) 1 0.5 0 0.5 1 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 seconds Now consider a second IVP in which t0 = 0, y0,2 = −0.25 volts and u2 (t) = 1, t ≥ 0. This solution is denoted y2 , Do not distribute these notes 21 Prof. R.T. M’Closkey, UCLA IVP with u2 and x0,2
2
IVP sol.
solution with IC=x , u=0 0,2 solutoin with IC=0, non zero u
input u 1.5 y2 = −10y2 + 10u2 ˙ y2 (0) = y0,2 = −0.25 =⇒ u2 (t) = 1, t ≥ 0 x2 (volts) 1 0.5 0 0.5 1 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 seconds It’s clear from these ﬁgures that the IVP that generates the solution y1 + y2 (the sum of the blue traces in the two
preceding ﬁgures) has initial condition y0,1 + y0,2 and input u1 + u2 as shown below, Do not distribute these notes 22 Prof. R.T. M’Closkey, UCLA IVP with u +u and x
1 2 +x 0,1 0,2 2
IVP sol.
solution when IC=x0,1+x0,2, u=0
solution when IC=0, u = u1 + u2
input u 1.5 y = −10y + 10u
˙ =⇒
y (0) = y0,1 + y0,2 = 0.75 u(t) = 1 + sin(10t), t ≥ 0 x (volts) 1 0.5 0 0.5 1 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 seconds ♣ Do not distribute these notes 23 Prof. R.T. M’Closkey, UCLA Timeinvariance property The consequence of timeinvariance (that is, a and b are constant in (7)) is the following: suppose an initial condition
is speciﬁed at t0 and the input u is speciﬁed for t ≥ t0 and let y be the subsequent “output", deﬁned for t ≥ t0 , that
˜
˜
is obtained by solving the IVP; if the same initial condition is shifted to time t0 + T , for any T , and if the new input
u is shifted by the same amount, then the solution to the IVP is just y shifted by T . This is shown more rigorously
˜
using (8).
Let the initial condition be speciﬁed at “starting time" t0 as y0 , and let the input, denoted u, be deﬁned for all t ≥ t0 .
˜
Then, the solution to the IVP, denoted y , is
˜ y (t) = e
˜ a(t−t0 ) t y0 + ea(t−τ ) bu(τ )dτ,
˜ t ≥ t0 . t0 Now let’s shift the starting time to t0 + T , where T is any desired value. The input, denoted u, for this new initial
value problem is just u shifted by T , too:
˜ u(t) = u(t − T ),
˜ t ≥ t0 + T. The solution to this IVP is
a(t−(t0 +T ) y (t) = e t y0 + = ea(t−(t0 +T ) y0 +
= ea((t−T )−t0 ) y0 + t 0 +T
t ea(t−τ ) bu(τ )dτ, ea(t−τ ) bu(τ − T )dτ,
˜ t ≥ t0 + T, ea((t−T )−s) bu(s)ds,
˜ t ≥ t0 + T, t 0 +T
t−T
t0 = y (t − T ),
˜ t ≥ t0 + T, t ≥ t0 + T, which demonstrates that shifting the initial condition and input by T just shifts the solution by T . Do not distribute these notes 24 Prof. R.T. M’Closkey, UCLA For timeinvariant systems we make take t0 = 0 without loss of generality. Thus, t0 = 0 will be the default “starting
time" that we assume from now on.
Example.
Consider the lowpass ﬁlter circuit used in the superposition example, y = −10y + 10u.
˙
Suppose the starting time is t0 = 0.5 seconds and y0 = −0.25 and u(t) = 1, t ≥ 0.5. Then the solution to this IVP
is,
IVP solution
2
IVP sol.
solution when IC= 0.25, u=0
solution when IC=0, u = 1
input u 1.5 y = −10y + 10u ˙ y (0.5) = −0.25
=⇒ u(t) = 1, t ≥ 0.5 x (volts) 1 0.5 0 0.5 1 0 0.2 0.4 0.6 0.8 1 1.2 1.4 1.6 1.8 2 seconds Note that this it the solution of a previously studied IVP that is just shifted to the new “starting time". This property
does not hold if the coefﬁcients in the ODE, i.e. a and b, are functions of time.
♣ Do not distribute these notes 25 Prof. R.T. M’Closkey, UCLA Asymptotic Stability of First Order ODEs
The ﬁrst order ODE y = ay + bu,
˙
is called asymptotically stable if all solutions to the homogeneous IVP converge to zero as t → ∞. Since the
solutions to any homogeneous IVP is given by yh (t) = ea(t−t0 ) y0 ,
then asymptotic stability is equivalent to lim ea(t−t0 ) = 0, t→∞ for which the necessary and sufﬁcient condition is Re(a) < 0. Note that the deﬁnition of asymptotic stability does
not depend on the input –this is a characteristic associated with any linear system. In nonlinear systems there is
domain of attraction associated with a given equilibrium point and in this case it may be possible for a “large" input
to kick the dependent variable out of this domain.
One consequence of asymptotic stability is that any bounded input produces a bounded dependent variable.
In other words, if u(t) < m for all t ≥ t0 , where m > 0 is the bound on the magnitude of u, then a similar bound
can be computed for y if and only if Re(a) < 0, y (t) ≤ e
≤ e a(t−t0 ) t y0  + ea(t−τ ) u(τ )dτ , t0 (Re(a)+j I m(a))(t−t0 )
t ≤ y0  + t0 y0  + t0 ea(t−τ ) u(τ )dτ ea(t−τ ) mdτ
t = y0  + m t t ≥ t0 t0 eRe(a)(t−τ ) ej I m(a)(t−τ ) dτ m
≤ y0  −
.
Re(a)
Do not distribute these notes 26 Prof. R.T. M’Closkey, UCLA This bound is positive for any Re(a) < 0 and any initial condition. Note that as the convergence rate becomes
weaker, i.e. Re(a) → 0, the bound on y becomes larger. Do not distribute these notes 27 Prof. R.T. M’Closkey, UCLA Another View of First Order Systems: Convolution Representation
The integral in the solution we derived for the IVP is called a convolution of the functions eat and u(t),
t at y (t) = e y0 + 0 ea(t−τ ) bu(τ )dτ , t ≥ t0 . (11) convolution The lower limit of integration corresponds to the time at which the initial condition is speciﬁed, i.e. at t = 0. The
form of (11) suggests another representation of the dependent variable as a function of the input: y (t) = ∞
−∞ h(t − τ )u(τ )dτ, (12) where the function h is called the impulse response of the system. The impulse response concept will be studied
in great detail in this course.
In the case of the ﬁrst order, linear, “low pass" form of the ODE, h is given by h(t) = 0 t<0 eat b t ≥ 0 . (13) Since h = 0 for t < 0, the dependent variable y at any given time does not depend on values of u that are in the
future relative to that time. This is an entirely reasonable property that is consistent with our observation of the
physical world. Such systems are called causal. Thus, for causal systems (12) can be modiﬁed to
t y (t) =
−∞ Do not distribute these notes h(t − τ )u(τ )dτ. 28 Prof. R.T. M’Closkey, UCLA This expression can be further split into two intervals of integration to recover (11),
t y (t) = ea(t−τ ) bu(τ )dτ −∞
0 = e
−∞
at a(t−τ ) 0 =e t bu(τ )dτ +
0 a(0−τ ) e ea(t−τ ) bu(τ )dτ
t bu(τ )dτ + −∞ 0 (assume t > t0 ) ea(t−τ ) bu(τ )dτ y0 The interpretation is that the initial condition y (0) = y0 is created by the action of the input u on the interval
(−∞, 0]. With the exception of certain special cases, the idea that an input acting for t < 0 can “set up" a value of
the dependent variable at t = 0 extends to linear ODEs of any order.
The representation (12) of the dependent variable as a function of time has an appealing simplicity and only
requires that we • extend the domain of deﬁnition of u to (−∞, ∞),
• deﬁne h according to (13)
We will see that h plays a fundamental role in the study of linear systems and is the response of the system to a
unit impulse function and, hence, is called the system’s impulse response. Do not distribute these notes 29 Prof. R.T. M’Closkey, UCLA Example. Consider the circuit example again (the ODE is y = −10y + 10u) with the following input deﬁned on the
˙
interval (−∞, ∞), u(t) = −0.3e2t t < 0
1
t≥0 h(t) = 0
t<0
,
10e−10t t ≥ 0 The impulse response of the circuit is and is shown below
Impulse Response
10 8 volts 6 4 2 0
1 Do not distribute these notes 0.5 0
seconds 30 0.5 1 Prof. R.T. M’Closkey, UCLA Computing (12) is done over two time intervals,
t For t < 0, y (t) =
−∞ 10e−10(t−τ ) (−0.3)e2τ dτ = −0.25e2t
0 For t ≥ 0, y (t) = 10e
−∞ −10(t−τ ) −10t = −0.25e t 2τ (−0.3)e dτ +
−10t +1−e 0 10e−10(t−τ ) dτ . The input and output of the circuit are shown below,
Input Output 1.5 1 1
Vout (volts) 2 1.5 Vin (volts) 2 0.5 0.5 0 0 0.5 0.5 1
1 0.5 0
seconds 0.5 1 1
1 0.5 0
seconds 0.5 1 Note that for t ≥ 0, u and y match one of the previous examples we considered when the IC was speciﬁed as
y0 = −0.25, however, in the present example, we did not specify an initial condition but rather expanded the
domain of deﬁnition of the input.
♣ Do not distribute these notes 31 Prof. R.T. M’Closkey, UCLA Yet Another View of First Order Systems: the Transfer Function
Suppose the input has “exponential form": u(t) = est where s may be a complex number (this is denoted by s ∈ C).
Although it may seem confusing as to why we would allow an input to take on complex values and, furthermore,
question its physical signiﬁcance when the engineering quantities we work with are realvalued, it is very convenient
to work with complex algebra when studying linear systems even when we are interested in solutions for problems
in which all of the data are real. A wide variety of functions can be represented by est where s ∈ C: • constant: s = 0 =⇒ Re(est ) = 1 for all t
• sinusoid with frequency ω , phase φ and amplitude α > 0: s = jω , β ∈ C such that β  = α and ∠β = φ =⇒
Re(βest ) = α cos(ωt + φ) for all t
• damped sinusoid with frequency ω and decay rate σ < 0: s = σ + jω =⇒ Re(est ) = eσt cos(ωt) for all t
The notation Re(·) means “take the real part of its argument". Now let’s ﬁnd a particular solution also of exponential
form yp (t) = Hest , where H is a (possibly complex) constant that is determined by substitution into the nonhomogeneous ODE: d
H est = aHest + best
dt =⇒ Hsest = aHest + best =⇒
=⇒ Hs = aH + b
H (s − a) = b. (because est = 0 for all t) If we further assume s = a then H = b/(s − a) and hence a particular solution is yp (t) = Hest = b st
e.
s−a When s = a then the assumed particular solution cannot satisfy the ODE so you must continue your search (a
particular solution of the form teat works). Do not distribute these notes 32 Prof. R.T. M’Closkey, UCLA We can consider H to be a function of s deﬁned by b
,
s−a H (s) = where s ∈ C but s = a. This function is called the transfer function of the system.
Example. Consider the lowpass ﬁlter ODE y = −10y + 10u. The transfer function of this system is
˙ H (s) = 10
.
s + 10 Suppose the input is deﬁned to be u(t) = e−2t cos(5t). This function can be represented as an exponential with
complex constant s = −2 + j 5: u(t) = Re e(−2+j 5)t . If the function e(−2+j 5)t is applied to the system, then the
particular solution that is also a scaled version of this input is (the scaling constant is the transfer function evaluated
at s = −2 + j 5), 10
e(−2+j 5)t .
−2 + j 5 + 10
This particular solution satisﬁes the ODE when the forcing function is given by est , not the actual input u(t) =
Re(est ). The particular solution can be split into its real and imaginary parts, i.e. yp (t) = Re(yp (t)) + I m(yp (t))
yp (t) = H (−2 + j 5)e(−2+j 5)t = and substituted into the ODE: d
(Re(yp (t)) + I m(yp (t))) = −10 (Re(yp (t)) + I m(yp (t))) + 10e(−2+j 5)t
dt
d
=⇒
Re(yp (t)) = −10Re(yp (t)) + 10 Re e(−2+j 5)t
dt
u d
I m(yp (t)) = −10I m(yp (t)) + 10I m e(−2+j 5)t
dt
The middle equation shows that Re(yp ) is a particular solution associated with the actual input u, d
Re(yp (t)) = −10Re(yp (t)) + 10e−2t cos(5t).
dt
Note that the imaginary part of yp also satisﬁes the ODE when the input is I m e(−2+j 5)t . A plot of u and Re(yp )
are shown below. Do not distribute these notes 33 Prof. R.T. M’Closkey, UCLA Exponential form input and corresponding particular solution
6
x p u 5 4 u(t) = e−2t cos(5t) = Re e(−2+j 5)t
=⇒ particular solution is
Re H (s)est , s = −2 + j 5
10
= Re
e(−2+j 5)t
−2 + j 5 + 10 3 2 1 0 1 2 3 4 1 0.8 0.6 0.4 0.2 0 0.2 0.4 0.6 0.8 1 seconds Thus, the transfer function enables us to quickly produce a particular solution when the input can be expressed as
Re(est ) for some s ∈ C. The transfer function has a more general role in describing solutions of ODEs but these
details will have to wait until we introduce the Laplace transform.
♣ Do not distribute these notes 34 Prof. R.T. M’Closkey, UCLA A Special Case of the Transfer Function: The Frequency Response Function When u is sinusoidal it can be expressed as the real part of a complex exponential function u(t) = cos(ωt) =
R ejωt . Thus, the sinusoidal particular solution is √
yp (t) = Re H (jω )ejωt
(j = −1)
= H (jω ) cos(ωt + ∠H (jω )). When setting s = jω in the transfer function, i.e. H (jω ), the resulting expression is of such fundamental importance in the study of linear systems that it is given its own name: H (jω ) is the frequency response function
of the system. In general, H (jω ) is complexvalued, but it can be expressed in terms of two real quantities: the
magnitude, which is denoted H (jω ), and the phase, which is denoted ∠H (jω ). The magnitude and phase can
be plotted versus frequency to generate frequency response plots. These plots are an extremely effective way of
assessing the behavior of the system to sinusoidal inputs. For asymptotically stable systems we will also see that
there is a deep connection between the impulse response of the system and its frequency response function.
Facts pertaining to asymptotically stable systems
If u(t) = α cos(ωt + φ), α, ω, φ ∈ R, where φ is some arbitrary phase, then y (t) = eat y0 − H (jω )α cos(φ + ∠H (jω ) + H (jω )α cos(ωt + φ + ∠H (jω )), t ≥ 0, where y0 is the initial condition and H (jω ) = b/(jω − a). This unique solution is most easily derived using the ad
hoc approach discussed earlier. If you have enough fortitude you can show that the convolution representation, i.e.
t at y (t) = e y0 + 0 ea(t−τ ) α cos(ωτ + φ) dτ,
u(τ ) reduces to y (t) = eat y0 −eat H (jω )α cos(φ + ∠H (jω ) + H (jω )α cos(ωt + φ + ∠H (jω )) .
zero state response Do not distribute these notes 35 Prof. R.T. M’Closkey, UCLA Regrouping these terms reveals what is termed the steadystate response and transient response, y (t) = eat y0 − H (jω )α cos(φ + ∠H (jω ) + H (jω )α cos(ωt + φ + ∠H (jω )),
steadystate response transient response t ≥ 0. If the ﬁrst order system is asymptotically stable then y converges to the steadystate response, denoted yss , which
is given by yss (t) = H (jω )α cos(ωt + φ + ∠H (jω )).
This demonstrates that if the dependent variable is observed until the “transient" term decays to an undetectable
level, the frequency response function at frequency ω can be empirically measured by normalizing the amplitude
of yss by the amplitude of u, which provides and estimate of H (jω ), and by measuring the phase difference
between yss and u, which provides an estimate of ∠H (jω ). Thus, one robust way of measuring the frequency
response of a system is to choose a grid of test frequencies, perform the sinusoidal input experiments, and make
the measurements of the magnitude and phase as previously described.
Facts pertaining to neutrally stable or unstable systems
If the system is neutrally stable or unstable, then yss is never empirically observed because the transient term
does not decay to zero. If the initial condition is y0 = H (jω )α cos(φ + H (jω )), then the transient term is
zero and so y (t) = yss (t) for all t ≥ 0, however, in practice is not possible to specify the initial condition with
arbitrary accuracy and, furthermore, there are likely other disturbance sources which drive the system, too, thereby
producing additional transient terms which are not under control of the engineer. Thus, frequency response testing
of such systems must be done in closedloop (it is possible to determine the openloop frequency response function
from the closedloop data). Do not distribute these notes 36 Prof. R.T. M’Closkey, UCLA Transfer Functions and Block Diagram Algebra
The expression s H (s)est = a H (s)est +b est ,
y (t) u(t) y (t) can be rearranged to yield block diagrams that represent the same underlying differential equation. For example,
writing y (t) = 1
(ay (t) + bu(t)) ,
s yields the block diagram u 1
s b y u
⇐⇒ b y a a
On the other hand, writing y (t) = 1
(sy (t) − bu(t)) ,
a (a = 0), yields the block diagram u b − 1
a y u
⇐⇒ − 1
a y d
dt s
Do not distribute these notes b 37 Prof. R.T. M’Closkey, UCLA Finally, if α ∈ C, then (s + α − α)y (t) = ay (t) + bu(t) =⇒ y (t) = 1
((a + α)y (t) + bu(t)) ,
s+α yields the block diagram α u b 1
s+α y u
⇐⇒ b a+α − y a+α The block diagram algebra is particularly straightforward when representing “dynamic" blocks as transfer functions.
In this case, the input can be assumed to be of exponential form, i.e. u(t) = est , and all signals in the diagram
take the form βest , where β is an appropriate (complex) constant that is determined by writing the relationships
amongst the signals and blocks. The output of a “gain" block is just the product of the input to the block times the
gain. The output of a “dynamic" block is the transfer function associated with the block times the exponential input
to the block.
Example. Consider the last block diagram in the three arrangements just discussed. Assume we don’t know the
transfer function of the system but wish to ﬁnd it from the block diagram. Let u(t) = est and y (t) = βest , where
β is to be determined (and represents the transfer function of the system). The input to the 1/(s + α) block is
best + (a + α)βest . Hence, 1
best + (a + α)βest = βest
s+α
from which we compute β = b/(s − a).
Do not distribute these notes ♣
38 Prof. R.T. M’Closkey, UCLA Integrator or Differentiator?
The three “realizations" derived above represent the same system, i.e. the ﬁrst order differential equation y =
˙
ay + bu, however, in most texts realizations with integrators, not differentiators, are typically used. This choice
is due to the legacy associated with electrical analog computers in which differential equations were analyzed by
building an electric circuit whose voltages represented the dependent variable and its time derivatives. The active
elements of the circuit (operational ampliﬁers) generate noise which would be ampliﬁed when implementing the
differential equation using differentiating elements as compared to a circuit with an equivalent transfer function
employing only integrating elements. Initial conditions are also easier to specify with integrating elements since
the integral operation is achieved by using a feedback capacitor in the op amp: the initial charge on the capacitor
creates a voltage potential that can be scaled to an appropriate initial condition. Do not distribute these notes 39 Prof. R.T. M’Closkey, UCLA Three Sides of the Same Coin
Firstorder, linear, constant coef. ODE y = ay + bu
˙
y (0) = y0 u(t) deﬁned for t ≥ 0 Convolution representation y (t) = ∞ −∞ “impulse response
function” (t ≥ 0) t ea(t−τ ) bu(τ )dτ 0 often a result of ﬁrstprinciples modeling
useful form for numerical analysis including
numerical solutions h(t − τ )u(τ )dτ u deﬁned for t ∈ (−∞, ∞) =⇒ y (t) = e y0 +
at
0
h(t) =
eat b Firstorder,
linear,
timeinvariant
system t<0
t≥0 Transfer Function
Representation “transfer function” s = jω
“frequency response function” H (s) = b
s−a H (j ω ) = b
jω − a u(t) = est =⇒ yp (t) = H (s)est for any s = a the impulse response function is often obtained as the result
of a test on a physical system; the next step is to match an
analytical model to the empirical impulse response function useful for computing particular solutions when the input is of
exponential form; more general interpretation once Laplace
transforms are studied; useful for manipulating block diagrams the convolution representation may also be used to deﬁne
linear systems that cannot be described by ODEs (in other
words, all ODEs “generate” an appropriate impulse response
function, however, not all impulse response functions are
associated with an ODE!) frequency response functions give graphical insight into system
behavior; more general interpretation once Fourier transforms
are studied; frequency response testing is reliable method for
identifying models of physical systems Do not distribute these notes 40 Prof. R.T. M’Closkey, UCLA ...
View
Full
Document
This note was uploaded on 08/08/2011 for the course MAE 107 taught by Professor Tsao during the Spring '06 term at UCLA.
 Spring '06
 TSAO

Click to edit the document details