Docsity
Docsity

Prepare for your exams
Prepare for your exams

Study with the several resources on Docsity


Earn points to download
Earn points to download

Earn points by helping other students or get them with a premium plan


Guidelines and tips
Guidelines and tips

numerical modeling.., Lecture notes of Numerical Methods in Engineering

lecture 10 to be used for the purpose of studying

Typology: Lecture notes

2019/2020

Uploaded on 03/29/2020

jennaxo
jennaxo 🇺🇸

3 documents

1 / 33

Toggle sidebar

This page cannot be seen from the preview

Don't miss anything!

bg1
Numerical Modeling in
Biomedical Systems
BME 125:305
Lecture 10
10/8/19
pf3
pf4
pf5
pf8
pf9
pfa
pfd
pfe
pff
pf12
pf13
pf14
pf15
pf16
pf17
pf18
pf19
pf1a
pf1b
pf1c
pf1d
pf1e
pf1f
pf20
pf21

Partial preview of the text

Download numerical modeling.. and more Lecture notes Numerical Methods in Engineering in PDF only on Docsity!

Numerical Modeling in

Biomedical Systems

BME 125:

Lecture 10

Newton- Cotes

Differentiation

Interpolation

Finite difference methods

Finite Difference Methods, Interpolation, and Integration

Spline interpolation

Newton’s polynomials

Errors in FD methods

Integration

Simpson’s rules

-5 1 2 3 4 5 6 7 8

0

5

10

15

20

x

y

7th order polynomial

x = [1 2 3 4 5 6 7 8]; y = [4 7 5 3 8 0 3 2 ];

Example: 7th^ order polynomial fit to 8 data points

f  x a 1 x^2 a 2 xa 3 Careful – descending powers!

polyfit returns the coefficients p 1 , p 2 ,… in vector “p”. Function polyval takes

these coefficients along with a set of x values (xi) and calculates the polynomial f (x)

xi = [1:0.1:8]; yii = polyval(p,xi);

figure; plot(x,y,'ro',xi,yii); xlabel('x'); ylabel('y'); title('7th order polynomial');

p = polyfit(x,y,7);

20000 2002 2004 2006 2008 2010 2012

20

40

60

80

100

x

f(x)

Shifted fit

20000 2002 2004 2006 2008 2010 2012

20

40

60

80

100

x

f(x)

Raw fit - Large x-values

figure; plot(x,y,'ro',xi,yi); xlim([2000 2012]); ylim([0 100]); xlabel('x'); ylabel('f(x)'); title('Raw fit - Large x-values');

x = [2001 2002 2003 2004 2005 2006 2007 2009 2010 2011]; xi = x-2006; y = [18 28 23 28 26 24 50 28 36 39];

Shifting Data for Polynomial Fitting

(See note on Canvas > Resources > Supplementary Materials for an explanation of why polyfit can experience problems fitting high order polynomials to large x-values)

x = [2001 2002 2003 2004 2005 2006 2007 2009 2010 2011]; y = [18 28 23 28 26 24 50 28 36 39]; p = polyfit(x,y,9); xi = [2001:0.1:2011]; yi = polyval(p,xi); figure; plot(x,y,'ro',xi+2006,yi); xlim([2000 2012]); ylim([0 100]); xlabel('x'); ylabel('f(x)'); title('Shifted fit');

p = polyfit(xi,y,9); xi = [-5:0.1:5]; yi = polyval(p,xi);

Equivalent to interp1(x,Y,xi,’spline’)

-5 0 2 4 6 8

0

5

10

15

20

x

y

Cubic Spline

-5 0 2 4 6 8

0

5

10

15

20

x

y

7th order polynomial

  • Cubic spline interpolation results in reduced fluctuation / overshoot between data points, compared to high-order polynomial fitting

Newton- Cotes

Differentiation

Interpolation

Finite difference methods

Finite Difference Methods, Interpolation, and Integration

Spline interpolation

Newton’s polynomials

Errors in FD methods

Integration

Simpson’s rules

x

f (x)

a (^) b

  • Split the area bounded by the curve f (x) and the x -axis into n vertical segments (strips / rectangles)
  • Calculate the area of each rectangle (height x width) and then sum the area of all rectangles to find the total area under f (x) between a and b.
  • The height of each rectangle equals the value of the function at the mid- point of the rectangle
  • The width of each rectangle is constant, equal to:

n

h ba

1 2 3 … (^) n

h

  • How can we modify this approach to generate a better approximation for the integral?

Numerical Integration of f ( x ):

Numerical Integration / Quadrature

  • Using rectangular strips is the simplest form of a more general strategy:  Replace the integrand f (x) with an approximate function whose area is easier to calculate:

     

b

a

n

b

a

I f x f x f (^) n  x  a 0 a 1 xa 2 x^2 anx^ n

 Replace the original integrand with an nth^ order polynomial function:

x

f (x)

a (^) b Zero-order polynomial (constant)

x

f (x)

a (^) b 1st-order polynomial (linear)

x

f (x)

a (^) b 2nd-order polynomial (quadratic)

The Trapezoid Rule

  ^    

b

a

b

a

I f x f 1 x f 1  x a 0 a 1 x

  • Replaces the exact integrand f (x) with a 1st order polynomial function f 1 (x) :

Recall Lecture 9:

x 0 x x 1

f (x)

f (x 1 )

f 1 (x)

f (x 0 )

1 0

1 0 (^1 0) x x x x

f x f x f x f x  

  

The 1st^ order polynomial that passes through [x 0 , f (x 0 )] and [x 1 , f (x 1 )] is:

x

f (x)

a (^) b 1st-order (linear)

     ^    x a

b a

f x f a f b f a  

   1

We replace the true f ( x ) (black curve) with a 1st^ order polynomial f 1 ( x ) (red line) as an approximation for the purpose of integration

x

f (x)

a b 1st-order (linear)

      

b

a

b

a

I f x f 1 x

 x a dx

b a

f b f a I f a

b

a

 

 (^)  

  

     ^    x a

b a

f b f a f x f a  

 1  

   ^  

2

I  ba f a  f^ b

  • Multiply out the second term in the square brackets, collect terms with x - dependence, integrate, rearrange:

 The trapezoid rule:

f (a)

f (b)

Example

  • Numerically integrate f  x  0. 2  25 x 200 x^2  675 x^3  900 x^4  400 x^5

from a = 0 to b = 0.

(^00) 0.2 0.4 0.6 0.8 1

1

2

3

4

x

y = f(x)

(^00) 0.2 0.4 0.6 0.

1

2

3

4

x

y = f(x)

1.640533 (exact)

… a true relative error of 89.5%

2

  1. 80.^20.^232 2

I  h f a  f b   

(^00) 0.2 0.4 0.6 0.8 1

1

2

3

4

x

y = f(x)

Using the Trapezoid rule gives

Trapezoid rule:

 3  2

Oh f a f b I h   The error associated with using the Trapezoid rule is O(h^3 ): 

(^00) 0.2 0.4 0.6 0.8 1

1

2

3

4

x

y = f(x) h = 0.

(^00) 0.2 0.4 0.6 0.8 1

1

2

3

4

x

y = f(x)

h = 0.

We can break down the interval from a to b into multiple smaller sections, each with smaller h, applying the Trapezoid rule to each section and summing the individual areas:

(^00) 0.2 0.4 0.6 0.8 1

1

2

3

4

x

y = f(x)^ h^ = 0.^

 This is the “multiple application” or Composite Trapezoid Rule

  • The composite trapezoid rule is fine for evaluating well-behaved functions to typical levels of accuracy
  • If higher accuracy is required, or if the underlying function is time consuming to evaluate, the composite trapezoid rule can become inefficient

 Can we develop more efficient numerical integration algorithms with errors that reduce more quickly as the number of segments (n) increases?

  ^    

b

a

n

b

a

I f x f x f (^) n  x  a 0 a 1 xa 2 x^2 anx^ n Use higher order polynomial functions to approximate f (x):

x

f (x)

a b

1st-order (linear)

Trapezoid rule

x

f (x)

a (^) b

2nd-order (quadratic)

Simpson’s 1/3 rule

x

f (x)

a (^) b

3rd-order (cubic)

Simpson’s 3/8 rule

Simpson’s 1/3 rule

      

b

a

b

a

I f x f 2 x f 2  x a 0 a 1 xa 2 x^2

  1. Use a 2nd order polynomial function to approximate f (x):

 ^0 ^4 ^1 ^ ^2  3

f x f x f x h I    2

b a h

This can be integrated and simplified to: (^) where 

“1/3 rule”

 If we know the value of f (x) at x = a, b, and midway between a and b, ( data points) we can use Simpson’s 1/3 rule

x

f (x)

a = x 0 b = x 2

  1. Approximate the integral by determining the area under this polynomial:

2 2

0 0

2 0 1 0 2 0 1

x x

x x

I  (^)  f x  (^) b  b x  x  b x  x x x dx

x 1