matlab assignment

ENGR 232 Dynamic Engineering

Systems

Lecture 4.2

2nd Order Linear Differential Equations

Dr. Michael Ryan

Agenda • Quick Review

– Homogeneous 2nd Order

• Auxiliary Equations with Complex Roots

• Take away - you will know how to solve 2nd order differential equations with constant coefficients and forcing functions f(t). And… this can be extended to any order constant coefficient D.E. !

Read Schaum’s Differential Equations - Chapters 9

21/27/2016 ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan

Review: Two Independent Solutions

solution, auxiliary equation.

Roots: Discriminant:

If D > 0, 2 distinct real roots, if D < 0, two distinct complex roots, if D = 0 two identical real roots.

• D> 0 Two distinct real roots {r1 r2}

• D = 0 Two identical real roots {r1 = r2}

• D< 0 imaginary roots – next lecture

1/27/2016 ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan Page 3

a ¢¢y + b ¢y + cy = 0, a ¹ 0, y 1 (t), y

2 (t)

y = ert ar2 + br + c = 0

r =

-b ± b2 - 4ac

2a D = b2 - 4ac

𝑦1 𝑡 = 𝑒𝑟1𝑡

𝑦2 𝑡 = 𝑒𝑟2𝑡

𝑦1 𝑡 = 𝑒𝑟1𝑡

𝑦2 𝑡 = 𝒕𝒆𝒓𝟏𝒕

𝑦 𝑡 = 𝐶1𝑒 𝑟1𝑡 + 𝐶2𝑒

𝑟2𝑡

General Solution

𝑦 𝑡 = 𝐶1𝑒 𝑟1𝑡 + 𝐶2𝑡𝑒

𝑟1𝑡

Linear combination

Review: Solve Homogeneous 2nd Order Equation with initial conditions

Equation and initial value:

Auxiliary equation (aka Characteristic equation) and roots:

r1 = 0 r2 = -1

General solution

Apply initial conditions (general notation)

1/27/2016 ENGR 232 Winter 2016 Lecture 4.2 – Dr. M.

Ryan 4

¢¢y + ¢y = 0, y(0) = 2, ¢y (0) = 1

y(t) = c 1 e

r 1 t + c

2 e

r 2 t =

y(0) = c 1 e0 + c

2 e0 = c

1 + c

2 =

¢y (0) = r 1 c

1 e0 + r

2 c

2 e0 =

2

(0)c1 + (-1)c2 = 1

𝐜𝟐 = −𝟏 𝒄𝟏 = 𝟑

𝒚 𝒕 = 𝟑 − 𝒆−𝒕

𝑐1 + 𝑐2𝑒 −𝑡

𝑟2 + 𝑟 = 0 = 𝑟 𝑟 + 1

• The mass is attached to a spring and is free to move horizontally in a single dimension.

• y is the displacement from equilibrium.

• Force exerted by the spring is proportional to y (Hooke’s law).

• Assume the friction force is proportional to the velocity, called viscous friction – other types of friction (static and coulomb) would make the D.E non linear

• Model

• Set b = 0 then is a solution

Review: Mass-Spring Oscillator Model

1/27/2016 ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan Page 5

F

spring = -ky

F

friction = -bv, v =

dy

dt

𝒅𝟐𝒚

𝒅𝒕𝟐 +

𝒃

𝒎

𝐝𝐲

𝐝𝐭 +

𝒌

𝒎 𝐲 = 𝟎

𝑦 𝑡 = cos(𝑎𝑡)

System with friction using numerical values

ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan Page 6

• m = 1, k = 25, b = 6.

assume solution

Does this y(t) satisfy the differential equation?

Substitute into the equation,

1/27/2016

𝑦′′ + 𝑏

𝑚 𝑦′ +

𝑘

𝑚 𝑦 = 0

𝒚′′ + 𝟔𝒚′ + 𝟐𝟓𝒚 = 𝟎

𝒚 𝒕 = 𝒆−𝟑𝒕𝐜𝐨𝐬(𝟒𝒕)

𝒚′ = −𝟑𝒆−𝟑𝒕𝐜𝐨 𝐬 𝟒𝒕 − 𝟒𝒆−𝟑𝒕𝒔𝒊𝒏(𝟒𝒕)

𝒚′′ = −𝟕𝒆−𝟑𝒕𝐜𝐨 𝐬 𝟒𝒕 + 𝟐𝟒𝒆−𝟑𝒕𝒔𝒊𝒏(𝟒𝒕)

Some calculus yields

𝒚′′ + 𝟔𝒚′ + 𝟐𝟓𝒚 = −𝟕𝒆−𝟑𝒕𝐜𝐨 𝐬 𝟒𝒕 + 𝟐𝟒𝒆−𝟑𝒕𝒔𝒊𝒏 𝟒𝒕 +

𝟔(−𝟑𝒆−𝟑𝒕𝐜𝐨 𝐬 𝟒𝒕 − 𝟒𝒆−𝟑𝒕𝒔𝒊𝒏 𝟒𝒕 ) + 𝟐𝟓 𝒆−𝟑𝒕𝐜𝐨 𝐬 𝟒𝒕 = 𝟎

y(t) satisfies D. E. coefficients of sin() and cos() =-0

Another set of numerical values

ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan Page 7

• m = 1, k = 25, b = 10.

assume solution

Does this y(t) satisfy the differential equation?

Substitute into the equation,

1/27/2016

𝑦′′ + 𝑏

𝑚 𝑦′ +

𝑘

𝑚 𝑦 = 0

𝒚′′ + 𝟏𝟎𝒚′ + 𝟐𝟓𝒚 = 𝟎

𝐲′ 𝐭 = −𝟓𝐞−𝟓𝐭

𝐲′′ 𝐭 = 𝟐𝟓𝐞−𝟓𝐭

𝒚′′ + 𝟏𝟎𝒚′ + 𝟐𝟓𝒚 = 𝟎 = 𝟐𝟓𝐞−𝟓𝐭 + 𝟏𝟎 −𝟓𝐞−𝟓𝐭 + 𝟐𝟓 𝐞−𝟓𝐭

y 𝐭 = 𝐞−𝟓𝐭

Yes – so it is a solution….There are infinitely more

Summary: Three Cases

ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan Page 8

¢¢y + 25y = 0

y(t) = e-3t cos4ty(t) = cos5t

¢¢y + 6 ¢y + 25y = 0 ¢¢y +10 ¢y + 25y = 0

y(t) = e-5t

1/27/2016

𝑎𝑦′′ + 𝑏𝑦′ + 𝑐𝑦 = 0

general solution is the linear combination of two terms

First 2: 𝐜𝐨𝐬 𝝑𝒕 = 𝒆𝒊𝜽𝒕+𝒆−𝒊𝜽𝒕

𝟐 , the last one has a single term

𝒚 𝒕 = 𝑪𝟏𝒆 𝒓𝟏𝒕 + 𝑪𝟐𝒆

𝒓𝟐𝒕

Brief Review of Complex Numbers

1/27/2016 ENGR 232 Winter 2016 Lecture 4.2 – Dr. M.

Ryan 9

𝑟 = 𝑥2 + 𝑦2

𝜑 = tan−1 𝑦

𝑥

Any complex number can be thought of as a two dimensional vector, where the first dimension represents the real component and the second dimension represents the imaginary component. For example, the complex number 𝑧 = 𝑥 + 𝑖𝑦 can be plotted as:

Example - D.E. having complex roots in auxiliary equation • Homogeneous 2nd order IVP • Auxiliary Equation

• General solution and derivative

• Apply initial conditions

• Set up matrix equation and use linear algebra to solve for coefficients

1 1 −3 + 𝑖4 −3 − 𝑖4

𝑐1 𝑐2

= 2 2

1/27/2016 ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 10

¢¢y + 6 ¢y + 25y = 0, y(0) = 2, ¢y (0) = 2

y = c 1 e

r 1 t + c

2 e

r 2 t =

¢y = r 1 c

1 e

r 1 t + r

2 c

2 e

r 2 t =

r2 + 6r + 25 = 0 𝑟 = −6± 36−4(1)(25)

2(1) = −𝟑 ± 𝒊𝟒

𝑦 0 = 𝑐1 + 𝑐2 = 2

𝑟1 = −3 + 𝑖4 𝑟2 = −3 − 𝑖4

𝑦′ 0 = 𝑟1𝑐1 + 𝑟2𝑐2 = 2 = −3 + 𝑖4 𝑐1 + −3 − 𝑖4 𝑐2 = 2

Use MATLAB to do the work of complex number computation

ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 11

%% Second order differential equation example.

clear, clc, close all

% The roots

r1 = -3+4i; r2 = -3-4i;

% Matrix for matching initial conditions

M = [1 1; r1 r2]

% The initial conditions

IC = [2; 2]

% c1 and c2, the coefficients

c1c2 = inv(M)*IC

c1 = c1c2(1); c2 = c1c2(2);

%% Plot the solution

t = 0:0.01:3;

y = c1*exp(r1*t) + c2*exp(r2*t);

plot(t, y)

grid on

1 1 −3 + 𝑖4 −3 − 𝑖4

𝑐1 𝑐2

= 2 2

How would you check to see if the above plot of solution is “basically correct”?

c1c2 = 1.0000 - 1.0000i 1.0000 + 1.0000i

1/27/2016

Note we obtained a real answer from y(t) = f (complex numbers)

(1) IC (2) Final value (3) shape based on roots of aux eq.

Euler’s Formula Dealing with Complex numbers

If r = (α + βi) what does ert = ???

The first term is fine but what is eiβt ?

Euler’s (Oiler’s) formula: eiθ = cosθ + isinθ (it might be good to memorize this formula – it will be used forever)

Applying Euler’s equation to the second term

we can express a complex exponential as an expression with sines and cosines albeit with an “i”

1/27/2016 ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 12

ert = e(a+ib )t = eateibt = eat cosbt + isinbt( )

e rt = e(a+ib )t = eateibt

1/27/2016

Applying this to the solution of D.E.

Simplify the first term,

c1er1t = (1-i)*e-3t(cos4t + i sin4t)

Recall how to multiply complex numbers: (a + ib)(c + id) = (ac – bd) + i(ad + bc)

In a similar way,

ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 13

¢¢y + 6 ¢y + 25y = 0, y(0) = 2, ¢y (0) = 2

y(t) = c 1 e

r 1 t + c

2 e

r 2 t , c

1 = 1- i, r

1 = -3+ 4i, c

2 = 1+ i, r

2 = -3- 4i

c

1 e

r 1 t = e-3t cos4t + sin4t + i sin4t - cos4t( )é

ë ù û

c

2 e

r 2 t = e-3t cos4t + sin4t - i sin4t - cos4t( )é

ë ù û

y(t) = c 1 e

r 1 t + c

2 e

r 2 t =

𝑦 𝑡 = 𝑐1𝑒 𝑟1𝑡 + 𝑐2𝑒

𝑟2𝑡 = 𝑒𝛼𝑡 𝑑1𝑐𝑜𝑠 𝛽𝑡 + 𝑑2𝑠𝑖𝑛 𝛽𝑡

Note that when we add c1er1t and c2er2t

the complex terms cancel leaving us with only the real terms

𝟐𝒆−𝟑𝒕 𝒄𝒐𝒔𝟒𝒕 + 𝒔𝒊𝒏𝟒𝒕

In general if r = α +/- iβ

Solution

Summary: Solution for D.E. with complex roots of auxiliary equation

Solution method:

(1) Find r1, r2, roots of auxiliary equation (complex numbers)

(2) Set up initial condition equation

(3) Solve for c1 and c2 (complex numbers)

(4) Solution is in terms of complex numbers

(5) Simplify if desired to get form

1/27/2016 ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 14

a ¢¢y + b ¢y + cy = 0, y(0) = Y 0 , ¢y (0) = Y

1

a ¹ 0; a, b, c, Y 0 , Y

1 real; D = b2 - 4ac < 0

c 1 + c

2 = Y

0

r 1 c

1 + r

2 c

2 = Y

1

y(t) = c 1 e

r 1 t + c

2 e

r 2 t

This process is the same as the case where r1 and r2 are real roots

𝑦 𝑡 = 𝑐1𝑒 𝑟1𝑡 + 𝑐2𝑒

𝑟2𝑡 = 𝑒𝛼𝑡 𝑑1𝑐𝑜𝑠 𝛽𝑡 + 𝑑2𝑠𝑖𝑛 𝛽𝑡

Alternate Method (you may prefer this)

(1) Find r1, r2, roots of auxiliary equation (complex numbers)

(2) Choose one root: r1=α + iβ,

(3) Real Solution has the form

(4) Identify form of equations for initial conditions

(5) Write out solution to IVP

1/27/2016 ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 15

a ¢¢y + b ¢y + cy = 0, y(0) = Y 0 , ¢y (0) = Y

1

a ¹ 0; a, b, c, Y 0 , Y

1 real; D = b2 - 4ac < 0

y(t) = eat d

1 cosbt + d

2 sinbtéë ùû

y(0) : d 1 = Y

0

¢y (0) : ad 1 + bd

2 = Y

1

y(t) = eat d

1 cosbt + d

2 sinbtéë ùû

As a check verify that final solution satisfies the IC’s

Solve for d1 d2

Example: Spring-Mass-Damper System

• Spring-mass-friction system,

• m = 36, b = 12, k = 37, y(0) = 0.7, y’(0) = 0.1

(1) Auxiliary equation (2) Roots

(3) General solution:

(4) Initial value equation:

(5) Solution:

1/27/2016 ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 16

m ¢¢y + b ¢y + ky = 0

r = -

1

6 ± i

y(0) : d 1 = 0.7

¢y (0) : - 1

6 d

1 + d

2 = 0.1

y(t) = e-t/6(0.7cost + 0.2167sint)

𝟑𝟔𝒓𝟐 + 𝟏𝟐𝒓 + 𝟑𝟕 = 𝟎

𝛼 = − 1

6 𝛽 = 1

y(t) = eat d

1 cosbt + d

2 sinbtéë ùû

y(0) : d 1 = Y

0

¢y (0) : ad 1 + bd

2 = Y

1

Plotting solution

• Solution

ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 17

y(t) = e-t/6(0.7cost + 0.2167sint)

Is the IC correct?

What is the final value?

We see the characteristic shape of a damped sinusoidal

Response – this will always be the case if the real part of the root is negative

1/27/2016

Linear Homogeneous equation roots negative real unequal

a = 1,

D > 0, real roots

2 negative roots

ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 18

real or complex solutions to auxiliary equation

r1 = -3 r2 = -1

1/27/2016

We see the sum of two exponentials with different time constants – the time constant is time when

argument of e is -1. The two time constants are { 𝟏

𝟑 , 𝟏}

𝑦 𝑡 = 2𝑒−3𝑡 − 𝑒−𝑡

𝑎𝑠 𝑡 → ∞ 𝑦 𝑡 → 0 𝑎𝑠 𝑡 → −∞ 𝑦 𝑡 → +∞

𝑎𝑦′′ + 𝑏𝑦′ + 𝑐𝑦 = 0, 𝐷 = 𝑏2 − 4𝑐

Linear Homogeneous equation roots positive real unequal

a = 1,

D > 0, real roots

2 positive roots

ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 19

r1 = +3 r2 = +1

1/27/2016

Grows without bound as t increases

Time constant is not defined for this form

𝑦 𝑡 = 2𝑒3𝑡 − 𝑒𝑡

𝑎𝑠 𝑡 → ∞ 𝑦 𝑡 → ∞ 𝑎𝑠 𝑡 → −∞ 𝑦 𝑡 → 0

𝑎𝑦′′ + 𝑏𝑦′ + 𝑐𝑦 = 0, 𝐷 = 𝑏2 − 4𝑐

Linear Homogeneous equation one positive real root one negative real root

a = 1,

D > 0, two real roots

one negative

one positive

ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 20

r1 = +3 r2 = -1

1/27/2016

One term goes to zero and the other grows without bound at t approaches infinity

𝑦 𝑡 = 2𝑒3𝑡 − 𝑒−𝑡

𝑎𝑠 𝑡 → ∞ 𝑦 𝑡 → ∞ 𝑎𝑠 𝑡 → −∞ 𝑦 𝑡 → 0 𝑎𝑠 𝑡 → −∞ |𝑦 𝑡 | → ∞

𝑎𝑦′′ + 𝑏𝑦′ + 𝑐𝑦 = 0, 𝐷 = 𝑏2 − 4𝑐

Linear Homogeneous equation complex conjugate roots negative real part

a = 1,

D < 0, complex roots

real part of roots < 0

ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 21

r = -1+/- 10i

We see a sinusoidal waveform whose amplitude decreases as defined by exp(-at) – envelope

Damped sinusoidal

1/27/2016

𝑒−𝑡y(t) = 2*exp(-t)*cos(10*t)

𝑎𝑦′′ + 𝑏𝑦′ + 𝑐𝑦 = 0, 𝐷 = 𝑏2 − 4𝑐

Linear Homogeneous equation complex conjugate roots positive real part

a = 1,

D < 0, complex roots

real part of roots >0

ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 22

r = +1+/- 10i

1/27/2016

Not a desired response for most applications

y(t) = 2*exp(t)*cos(10*t) 𝑒𝑡

We see a sinusoidal waveform whose amplitude increases as defined by exp(t) – envelope Grows without bound

𝑎𝑦′′ + 𝑏𝑦′ + 𝑐𝑦 = 0, 𝐷 = 𝑏2 − 4𝑐

a = 1,

D < 0, complex roots r = +/- iB (no real component, 𝛼 = 0)

corresponds to b = 0 the damping term in D.E.

Response is purely oscillatory, plug into general solution

Since real part is zero the exponential, 𝑒𝛼𝑡 goes to 1

Purely oscillatory response (aka “undamped” response)

Linear Homogeneous equation complex conjugate roots NO real part

1/27/2016 ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 23

𝑦 𝑡 = 𝑐1𝑒 𝑖𝛽𝑡 + 𝑐2𝑒

−𝑖𝛽𝑡 = 𝑑1𝑐𝑜𝑠 𝛽𝑡 + 𝑑2𝑠𝑖𝑛 𝛽𝑡

𝑦 𝑡 = 𝑐1𝑒 𝑟1𝑡 + 𝑐2𝑒

𝑟2𝑡 = 𝑒𝛼𝑡 𝑑1𝑐𝑜𝑠 𝛽𝑡 + 𝑑2𝑠𝑖𝑛 𝛽𝑡

𝑎𝑦′′ + 𝑏𝑦′ + 𝑐𝑦 = 0, 𝐷 = 𝑏2 − 4𝑐

Relating roots of characteristic equation to response y′′ + by′ + 25y = f(t)

ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 24

Consider the response of a system to a step input such as if you were driving along then jumped a curb to a new constant level – f(t) The auxiliary equation can tell us how the system will respond

Overdamped – both roots negative Critically damped – both roots negative and equal Underdamped – complex conjugate roots with negative real parts

1/27/2016

You may be able to relate to the adjustment of shock absorbers or struts

ℎ𝑜𝑚𝑜𝑔𝑒𝑛𝑒𝑜𝑢𝑠 𝑒𝑞𝑢𝑎𝑡𝑖𝑜𝑛, f t is set to 0

y′′ + by′ + 25y = 0 𝑟2 + 𝑏𝑟 + 25 = 0

Summary • Second order homogeneous linear differential equations

• Auxiliary equation (characteristic equation) – find roots – general solution

– Roots define response • Look at y(t) as t approaches infinity

• General solution (either in terms of complex numbers on left or sines and cosines on right)

• Initial values allow determining of coefficients of general solution • Damping term (overdamped , underdamped or critically damped response – for a step input) value of b in diff

eq. way to classify systems

1/27/2016 ENGR 232 Winter 2016 Lecture 4.2 – Dr. M. Ryan 25

𝑦 𝑡 = 𝑐1𝑒 𝑟1𝑡 + 𝑐2𝑒

𝑟2𝑡 = 𝑒𝛼𝑡 𝑑1𝑐𝑜𝑠 𝛽𝑡 + 𝑑2𝑠𝑖𝑛 𝛽𝑡

𝑦 𝑡 = 𝑐1𝑒 𝑟1𝑡 + 𝑐2𝑒

𝑟2𝑡