Download Ordinary Differential Equations: Solutions to Problems and more Slides Numerical Methods in Engineering in PDF only on Docsity!
CE 341/441 - Lecture 20 - Fall 2004
p. 20.
LECTURE 20SOLUTION TO SINGLE 1ST ORDER INITIAL VALUE PROBLEMS (IVP’s) • Solve
i.c.
- Consider two classes of methods:
- Runge-Kutta type formulae
- single step methods• very simple to program• self starting (only need i.c.’s)
- Multi-step formulae
- Multi-step methods are much more efficient than single step methods (for the
same accuracy)
- Multi-step methods are not self starting
→
use single step method to start up and
then go over to multi-step
dy -----dt
-^
f^
y t,(
y t
o (^
)^
y^ o
CE 341/441 - Lecture 20 - Fall 2004
p. 20.
Runge-Kutta type formulas
→
Single Step Methods
is obtained in terms of
,^
and
evaluated for various values
of
between
and
⇒
self starting
- Self starting since solution involves only information between
⇒
therefore
all information required is available at the 1st step (i.e. the response function of theprevious step only).
INSERT FIGURE NO. 91 • Various orders of accuracy are available:
- 1st order - Euler• 2nd order - Improved Euler, Modified Euler• 4th order - Runge-Kutta
y^
j^
1 +^
y^
j^
f^
y^
j^
t^ ,j
(^
)^
f^
y t,(
y^
t^ j
t^ j
1
t^ j
t^
t^ j
1
x x x x x x
t
solution knownto here
solutiondesired here
y
yj
yj+ tj+
tj
t
CE 341/441 - Lecture 20 - Fall 2004
p. 20.
Runge-Kutta type methods • Solve• Recursive relationship for all Runge-Kutta methods:
where
(^
is therefore the slope as per the definition of
‘s and
‘s
- Select these coefficients such that you minimize errors.• Compare Taylor Series expansion of
and select the recursive relationship coef-
ficients such that you eliminate the appropriate error terms.
dy -----dt
-^
f^
t y,(
y^
)^
y^ o
y^
j^
1 +^
y^
j^
tΦ
y
j^
t^ j
t
,^
(^
a
g 1 1
a
g 2 2
a
g 3
3
a^ n
g n
g
1
f^
t y,(
≅^
g^1
f
g
2
f^
t^
p^1
t y,
p^2
tg
1
(^
g
3
f^
t^
p^3
t y
p^4
tg
2
(^
a^
i^
p^ i
y^
j^
1
CE 341/441 - Lecture 20 - Fall 2004
p. 20.
Euler Method • 1st order method Derivation 1 • Use forward difference approximation for
⇒
⇒
- Simply “march” forward in time from^ INSERT FIGURE NO. 93
dy -----dt
dy -----dt
-^
f^
y t,(
y^
j^
1 +^
y^
j
t
-^
f^
y^
j^
t^ ,j
(^
y^
j^
1 +^
y^
j^
t f
y^
j^
t^ , j
(^
t^
tj+
tj
t j
t^2
t^1
t^00
j^
j+
CE 341/441 - Lecture 20 - Fall 2004
p. 20.
Derivation 2 • Cast into generic Runge-Kutta form with
expanded to only 1 term
where
the local truncation error per time step
(1)
- Develop Taylor Series expansion for
about
y^
j^
1 +^
y^
j^
t a
1
g
1
E
L
E
L^
y^
j^
1 +^
y^
j^
ta
1
f^
t^ j
y^
j , (^
)^
E
L
y^
j^
1 +^
t^ j
y^
j^
1 +^
y^
j^
dyt -----dt
j
t (^
d-
2 y^2 dt
j
O
t (^
dy -----dt
j
f^
t^ j
y j , (^
d
2 y^2 dt
j
˙f t
j^
y^
j , (^
CE 341/441 - Lecture 20 - Fall 2004
p. 20.
(2)
- Now compare Equations (1) and (2)• Comparing terms we note that• Thus the Euler Method (substituting for
into the formula)
(3)
local
truncation error (i.e. per time step),
, is second order!
However this error builds up as we time step and will in fact become first order.
y^
j^
1 +^
y^
j^
t f
t^ j
y j , (^
)^
t (^
-^
˙f t
j^
y^
j , (^
)^
O
t (^
y^
j^
t a
1
f^
t^ j
y^
j , (^
)^
E
L
y^
j^
t f
t^ j
y^
j , (^
)^
t (^
-^
˙f t
j^
y^
j , (^
)^
O
t (^
a
1
E
L
t (^
–^
˙f t
j^
y^
j , (^
a
1
y^
j^
1 +^
y^
j^
t f t
j^
y^
j , (^
E
L
CE 341/441 - Lecture 20 - Fall 2004
p. 20.
(5)
- Taking the difference between Equations (4) and (5) defines the local truncation error
(6)
-^
= the truncation error of the Euler formula
per time step
- This is consistent with the
error in the forward difference approximation used
to evaluate
- This equation assumes that
is exact
- The total solution error at
is due to the truncation error at every time step (local
error)
plus
the error that has accumulated in all previous steps.
- Only for the first step when using the i.c.• For other steps the solution
carries an accumulated error!
Y
j^
1 +^
Y
j^
t f
t^ j
Y
j , (^
)^
t (^
-^
˙f t
j^
Y^
j , (^
)^
t (^
-^
˙˙f
t^
j^
Y
j , (^
)^
E
j^
1 +^
y^
j^
1 +^
Y
j^
1
t (^
-^
˙f t
j^
Y
j , (^
–^
O
t (^
E
j^
1 +^
O
t (^
O
t (^
dy -----dt
Y^
j
t^ j
1
y^
j^
Y^
j
=
y^
j
CE 341/441 - Lecture 20 - Fall 2004
p. 20.
- In general the total solution has an error which can be defined by examining the differ-
ence between the Euler formula, Equation (3) and our Taylor Series expansion, Equa-tion (5).
(7)
- The total solution error:
(8) (9)
- Also it can be shown that (by one of the mean value theorems):
(10)
y^
j^
1 +^
Y
j^
1
-^
y^
j^
Y
j
-^
t^
f^
t^ j
y j , (^
)^
f^
t^ j
Y
j , (^
(^
)^
t (^
–^
˙f t
j^
Y
j , (^
)^
O
t (^
ε^
j^
1 +^
y^
j^
1 +^
Y
j^
1
ε^
j^
y^
j^
Y
j
≡ f^
t^ j
y j , (^
)^
f^
t^ j
Y
j , (^
y^
j^
Y
j
-^
f ∂
y
t
j^
ξ^
j , (^
y^
j^
ξ^
j^
Y
j
CE 341/441 - Lecture 20 - Fall 2004
p. 20.
- Let’s apply Equation (11) to a simplified scenario
and
change every time step
- Estimate
by assuming that
and
are constants over the interval of interest
- Starting with the i.c.’s at
⇒
⇒
p^
j^
E
j^
1
ε^
j^
1 +^
p^
E
ε^
j^
1 +^
ε^
j^
t p
(^
)^
E
=^ j^
εo
ε^1
E
ε^2
E
t p
(^
)^
E
ε^2
E
t p
(^
ε^3
E
t p
(^
)^
t p
(^
)^
E
ε^3
E
t p
t (^
p
2
(^
CE 341/441 - Lecture 20 - Fall 2004
p. 20.
⇒
term > all other terms
ε^4
E
t p
t (^
p
2
(^
)^
tp
(^
)^
E
ε^4
E
p^
t E
p
2
t (^
E
t (^
p
3 E
E
O
t (^
ε^4
×^
O
t (^
pO
t (^
p
2 O
t (^
(^3) p
O
t (^
O
t (^
×
ε^4
E
O
t (^
ε^
j^
1 +^
j^
(^
)E
O
t (^
CE 341/441 - Lecture 20 - Fall 2004
p. 20.
INSERT FIGURE NO.INSERT FIGURE NO.INSERT FIGURE NO.INSERT FIGURE NO.INSERT FIGURE NO.INSERT FIGURE NO.INSERT FIGURE NO.INSERT FIGURE NO.INSERT FIGURE NO.INSERT FIGURE NO.INSERT FIGURE NO. 95
y y^1
slope nolonger exactat j=
slopeexact at j=
j=
j=
t
slope
f(t
,yj
)j y(t)
CE 341/441 - Lecture 20 - Fall 2004
p. 20.
- Convergence: A numerical method is convergent if (assuming no round off) the numer-
ical solution approaches the exact solution as
∆
t^
↓
- Stability: Deals with the artificial amplification of components of the numerical solu-
tion.
- Under certain circumstances, components of the discrete solution (often the short
wavelength components) experience artificial (i.e. not physical) sustained growthfrom time step to time step which ultimately leads to numerical overflow (i.e. thecomputer can not hold the numbers anymore)
→
unstable solution.
- Stability is a property of both the differential equation and the numerical method
→
i.e. the difference equations determine stable/unstable behavior.
- Discrete solution can be
- unconditionally unstable• conditionally stable
→
restrictions on time step
always
inaccurate. A stable scheme may be inaccurate (depends
on the truncation error).
t
