Description

1 attachmentsSlide 1 of 1attachment_1attachment_1

Unformatted Attachment Preview

Assignment #6

BME 3420 – Fall 2021

Dec 02, 2021

1- [5 points]

Another approach to solve an ODE using the Explicit Euler formula, is to use the formula to directly

find an equation that predicts S(t

j+1

)

given S(t ). For example, for the linearized pendulum

j

problem, where

0

dS(t)

= [

dt

−

1

] S(t)

g

0

l

we can use the Explicit Euler Formula to obtain

0

S(tj+1 ) = S(tj ) + h [

1

S(tj+1 ) = [

0

1

g

−

0

l

0

0

1

] S(tj ) + h [

1

S(tj+1 ) = [

1

g

−

l

0

] S(tj )

h

] S(tj )

gh

−

] S(tj )

1

l

However, this solution has the same problem as before. Here is a code that implements this solution

for the pendulum problem with with

1

S(0) = [

]

0

and

g

√

= 4

l

In [6]:

import numpy as np

import matplotlib.pyplot as plt

plt.style.use(‘seaborn-poster’)

h=0.1

time = np.arange(0,5+h,h)

s0 = np.array([[1],[0]])

s = np.zeros((len(time),2))

s[0,:] = s0.T

w = 4 #(g/l)^(1/2)

for i in range(1,len(time)-1):

s[i,:] = np.array([[1, h],[-(w**2)*h, 1]])@s[i-1,:]

plt.figure(figsize = (5, 5))

plt.plot(time, np.cos(w*time),label=’Exact Solution’)

plt.plot(time,s[:,0], label = ‘Explicit Euler Formula’)

plt.xlabel(‘Time (s)’)

plt.ylabel(‘$Theta(t)$’)

plt.legend()

plt.show()

An alternative is to implement the Intrinsic Euler Formula, which approximate the state S(t

dS(tj+1 )

S(tj+1 ) = S(tj ) + (tj+1 − tj )

Note that this equation depends on

dS(tj+1 )

dt

,

dt

which is unknonw!

However, we can use some algebra to solve this problem and obtain

dS(tj+1 )

0

= [

dt

1

g

−

Reemplazing in the Implicit Euler formula we obtain

l

0

] S(tj+1 )

j+1 )

by

0

S(tj+1 ) = S(tj ) + h [

[

1

0

0

−

1

[

g

l

0

l

] S(tj+1 ) = S(tj )

−h

h

1

] S(tj+1 ) = S(tj )

S(tj+1 ) = [

for any S(t

]

gh

j

−1

−h

S(tj )

1

l

)

] S(tj+1 ),

1

1

j+1

0

l

g

1

Which allow us to compute S(t

−

0

] S(tj+1 ) − h [

1

g

)

a) [2 points] Implement an algorithm that allows you to solve the pendulum problem using the

Intrisic Euler Formula. Plot the numerical and exact predictions, discuss your plot.

Hint : take a look at the code for the Extrinsic Euler Formula and make the necesary changes

As you observed in point a, the Intrisic Euler Formula does not provide a correct solution either. One

approach overshoots the solution and another undershoots the solution. A reasonable approach

then, is to find the average between the Intrinsic and Extrinsic formulas, this is called the trapezoidal

formula, and is given by

dS(tj )

h

S(tj+1 ) = S(tj ) +

[

dS(tj+1 )

+

dt

2

],

dt

b) [1 point] Demonstrate that for the pendulum problem, where

dS(tj )

0

= [

dt

−

dS(tj+1 )

the trapezoidal formula is given by

l

0

= [

dt

1

g

0

1

g

−

] S(tj )

l

0

] S(tj+1 )

⎡ 1

S(tj+1 ) =

⎣

−

gh

2

⎤

⎦

1

2l

−1

h

⎡

h

1

⎣−

2

gh

⎤

S(tj )

1 ⎦

2l

You should write the mathematical operations to derive this formula using pen and paper.

c) [2 points] Implement an algorithm that allows you to solve the pendulum problem using the

Trapezoidal Formula. Plot the numerical and exact predictions, discuss your plot.

2- [5 points]

A model that describes the relation between position, velocity, and acceleration of a mass in a

simple sprin-mass-damper system is given by

mẍ + bẋ + kx = F (t),

where F (t) is the force applied to the mass.

a) [3 points] Assume that F (t)

= 0

when b

. Create a plot for each condition and discuss your observations.

,

= 0 b = 1

, and b

= 10

For this problem, assume that m

for all t. Solve the ODE to determine the position of the mass

= 1

and k

= 5

. The initial conditions are

x(0) = 1

ẋ(0) = 0

and you are interested in the solution between [0, 20] with h

= 0.01

Hint: Remember that numerical methods can only deal with first order equations, so you need to

reduce the ODE order before solving the problem.

b) [2 points]

⎧0

F (t) = ⎨ 1

⎩

0

if t < 3
if 4 ≤ t ≤ 8
otherwise
Solve the ODE to determine the position of the mass when b
condition and discuss your observations.
= 1
and b
= 10
. Create a plot for each
For this problem, assume that m
= 1
and k
= 5
. The initial conditions are
x(0) = 1
ẋ(0) = 0
and you are interested in the solution between [0, 20] with h
In [ ]:
= 0.01
Purchase answer to see full
attachment
User generated content is uploaded by users for the purposes of learning and should be used following Studypool's honor code & terms of service.