first_order_notes - MAE107 Introduction to Modeling and...

Info iconThis preview shows page 1. Sign up to view the full content.

View Full Document Right Arrow Icon
This is the end of the preview. Sign up to access the rest of the document.

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 first 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, definitions of 0− and 0+ 7. Superposition=convolution, properties of convolution 8. Fourier series and response to periodic inputs; energy spectral density; identification with periodic inputs; discrete-time 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 finite energy; use in identification; 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. State-space 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 first 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, specified on some interval I u = the “forcing" or “nonhomogeneous" term that is known or specified for t ∈ I a, b = coefficients 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 specified at time t0 (in other words, y (t0 ) = y0 is specified) and the continuous forcing function u is known for all t ∈ I , then, find 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 satisfies 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 specified 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 Low-pass filter 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" filter? 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 left-hand side are much smaller than the terms on the right-hand side so that we find α ≈ 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 low-pass filter 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 € High-pass filter 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 low-pass filter 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 ˙ high-pass filter can be constructed from a low-pass filter in series with a “feedthrough" term. Define 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 low-pass filtered version of Vin which is also inverted in polarity. The output voltage Vout is reconstructed from Vout = V + Vin . The high-pass filter 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 flow rate into the tank Qin (the input variable) and water level h (the dependent variable): Pa ↓ Qin 6 h P1 Variable definitions: Qin = Qout = h= V= Pa = P1 = −→ Qout volume flow rate into tank (m3 /s) volume flow 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 orifice (N/m2 ) Fundamental relation governing water volume in the tank: d (V (t)) = Qin (t) − Qout (t). dt Model of flow through the orifice: Qout (t) = cA 2(P1 (t) − Pa ) , ρ where c is the orifice coefficient, A is the orifice area, and ρ is the density of water. The pressure at the orifice 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 flow 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 first order, linear, constant coefficient ODE y (t) = ay (t) + bu(t), ˙ (7) where u is defined 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 non-zero solution to the homogeneous problem, denoted yh , satisfies 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 first-order differential equations, a fundamental solution is a homogeneous solution in which yh (t0 ) = 1. In fact, for first-order constant coefficient differential equations all homogeneous solutions are of the form αeat for any non-zero constant α. Do not distribute these notes 15 Prof. R.T. M’Closkey, UCLA 2. Particular solution. For the non-homogeneous ODE, a particular solution, denoted yp , is any function that satisfies (7) (includes the non-homogeneous 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 specified at time t0 (in other words, y (t0 ) = y0 is specified) 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 specified). Then, find 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 first 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 difficulty is computing the particular solution. There are certain classes of non-homogeneous 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 satisfies the non-homogeneous differential equation: yp (t) = c(t)ea(t−t0 ) , t ≥ t0 , where c is to be determined by substitution into the non-homogeneous 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 satisfies 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 “Zero-state 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 non-zero 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. “zero-state response" The term “zero-state 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 define the terms “forced motion", “natural motion", “steady-state 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 specified as y0,1 and the input u1 is defined 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 low-pass filter 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 figures that the IVP that generates the solution y1 + y2 (the sum of the blue traces in the two preceding figures) 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 Time-invariance property The consequence of time-invariance (that is, a and b are constant in (7)) is the following: suppose an initial condition is specified at t0 and the input u is specified for t ≥ t0 and let y be the subsequent “output", defined 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 specified at “starting time" t0 as y0 , and let the input, denoted u, be defined 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 time-invariant 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 low-pass filter 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 coefficients 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 first 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 sufficient condition is Re(a) < 0. Note that the definition of asymptotic stability does not depend on the input –this is a characteristic associated with any linear system. In non-linear 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 specified, 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 first 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 modified 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 definition of u to (−∞, ∞), • define 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 defined 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 specified as y0 = −0.25, however, in the present example, we did not specify an initial condition but rather expanded the domain of definition 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 significance when the engineering quantities we work with are real-valued, 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 find 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 defined 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 low-pass filter ODE y = −10y + 10u. The transfer function of this system is ˙ H (s) = 10 . s + 10 Suppose the input is defined 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 satisfies 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 satisfies 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 complex-valued, 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 steady-state response and transient response, y (t) = eat y0 − |H (jω )|α cos(φ + ∠H (jω ) + |H (jω )|α cos(ωt + φ + ∠H (jω )), steady-state response transient response t ≥ 0. If the first order system is asymptotically stable then y converges to the steady-state 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 closed-loop (it is possible to determine the open-loop frequency response function from the closed-loop 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 find 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 first 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 amplifiers) generate noise which would be amplified 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 First-order, linear, constant coef. ODE y = ay + bu ˙ y (0) = y0 u(t) defined for t ≥ 0 Convolution representation y (t) = ￿ ∞ −∞ “impulse response function” (t ≥ 0) t ea(t−τ ) bu(τ )dτ 0 -often a result of first-principles modeling -useful form for numerical analysis including numerical solutions h(t − τ )u(τ )dτ u defined for t ∈ (−∞, ∞) =⇒ y (t) = e y0 + at ￿ ￿ 0 h(t) = eat b First-order, linear, time-invariant 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 define 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.

Ask a homework question - tutors are online