Math245 Computer Lab set #6, Spring 2009
Simulation of the forced 2ndorder ODE
Consider the forced, springmassdamper model described by the linear 2nd order ODE
d
2
y
dt
2
+ 2
λω
dy
dt
+
ω
2
y
=
f
0
sin
ω
f
t.
(1)
λ
,
ω
,
f
0
and
ω
f
are all constant parameters (natural frequency:
ω
=
p
k/m
, damping parameter:
λ
=
c/
(2
√
mk
), forcing amplitude:
f
0
=
F
0
/m
, forcing frequency:
ω
f
).
Example 1
Choose
ω
= 1,
λ
= 1
/
8,
f
0
= 1,
ω
f
= 1 and
y
(0) =
y
0
(0) = 0. Simulate the response of this system.
By using the simulated results, compute the peak amplitude
of the steady response and the phase shift
between the steady response and the forcing function.
•
Step #1:
Let
y
=
y
1
and
y
0
=
y
2
, and express our 2ndorder ODE by a set of 1storder ODEs.
y
0
1
=
y
2
,
(2)
y
0
2
=

2
λωy
2

ω
2
y
1
+
f
0
sin
ω
f
t,
(3)
RHS of these equations are coded in Matlab function ﬁle
fode.m
as follows:
function dy = fode(t,y,lmd,omg,omgf,f0)
dy = [ y(2)
2*lmd*omg*y(2)omg^2*y(1)+f0*sin(omgf*t) ];
where we set
λ
,
ω
,
f
0
and
ω
f
as parameter variables
lmd
,
omg
,
f0
and
omgf
, respectively.
•
Step #2:
Log on to Blackboard. Download Matlab function ﬁle
force.m
from /Course docu
ment/Discussion Sections/ Mﬁles/ to your current working directory. We will use
force
to compute
the peak amplitude and phase shift from the simulated results. To see how to use
force
, type
help
force
after Matlab command line. Following is the summary:
USE:
[ypeak, pshift] = force(t, y, omega_f)
INPUT:
t
: independent variable vector obtained by ode45