Lecture 26 Numerical Integration NewtonCotes Formulas Trapezoid rule

Lecture 26: Numerical Integration Newton-Cotes Formulas Trapezoid rule Simpson's 3/8 rule Boole’s rule

trapezoid 1. m function f 1=trapezoid 1(func 1, a, b, n) %Integration func 1 %using composite trapezoid rule at n points between a and b h=(b-a)/n; f 1=func 1(a)/2; for i=1: n-1 f 1=f 1+func 1(a+i*h); end f 1=f 1+func 1(b)/2; f 1=h*f 1; trapezoid 1 test. m

%examples of use of function f 1=trapezoid 1(func 1, a, b, n) %Integration func 1 %using composite trapezoid rule at n points between a and b a=1; b=5; f 1=inline('exp(-x)') disp(['a=', num 2 str(a), ' b=', num 2 str(b)]) %display limits of integration intexact=-exp(-b)+exp(-a); i=1: 10; nvec=2. ^(i); for i=1: length(nvec) f 1 int=trapezoid 1(f 1, a, b, nvec(i)); disp(['n=', num 2 str(nvec(i)), ' Error in composite trapezod rule =', num 2 str(f 1 int-intexact)]); %display value of function f(x) end

>> trapezoid 1 test f 1 = Inline function: f 1(x) = exp(-x) a=1 b=5 n=2 Error in composite trapezod rule =0. 11305 n=4 Error in composite trapezod rule =0. 029605 n=8 Error in composite trapezod rule =0. 0074926 n=16 Error in composite trapezod rule =0. 001879 n=32 Error in composite trapezod rule =0. 00047011 n=64 Error in composite trapezod rule =0. 00011755 n=128 Error in composite trapezod rule =2. 9389 e-005 n=256 Error in composite trapezod rule =7. 3474 e-006 n=512 Error in composite trapezod rule =1. 8369 e-006 n=1024 Error in composite trapezod rule =4. 5922 e-007

a=0; b=2*pi; f 1=inline('exp(sin(x)^2)') disp(['a=', num 2 str(a), ' b=', num 2 str(b)]) %display limits of integration intexact=trapezoid 1(f 1, a, b, 200); i=1: 5; nvec=2. ^(i); for i=1: length(nvec) f 1 int=trapezoid 1(f 1, a, b, nvec(i)); disp(['n=', num 2 str(nvec(i)), ' Error in composite trapezod rule =', num 2 str(f 1 int-intexact)]); %display value of function f(x) end

f 1 = Inline function: f 1(x) = exp(sin(x)^2) a=0 b=6. 2832 n=2 Error in composite trapezod rule =-4. 7337 n=4 Error in composite trapezod rule =0. 66447 n=8 Error in composite trapezod rule =0. 0034145 n=16 Error in composite trapezod rule =7. 8954 e-009 n=32 Error in composite trapezod rule =-3. 5527 e-015 f 1 = Inline function: f 1(x) = exp(sin(x)^2) a=0 b=0. 7854 n=2 Error in composite trapezod rule =0. 021343 n=4 Error in composite trapezod rule =0. 005305 n=8 Error in composite trapezod rule =0. 0013228 n=16 Error in composite trapezod rule =0. 00032898 n=32 Error in composite trapezod rule =8. 0649 e-005

a=0; b=pi/4; f 1=inline('exp(sin(x)^2)') disp(['a=', num 2 str(a), ' b=', num 2 str(b)]) %display limits of integration intexact=trapezoid 1(f 1, a, b, 200); i=1: 5; nvec=2. ^(i); for i=1: length(nvec) f 1 int=trapezoid 1(f 1, a, b, nvec(i)); disp(['n=', num 2 str(nvec(i)), ' Error in composite trapezod rule =', num 2 str(f 1 int-intexact)]); %display value of function f(x) end

simpson 1. m function f 1=simpson 1(func 1, a, b, n) %Integration func 1 %using composite simpson rule at n points between a and b h=(b-a)/n; f 1=func 1(a); for i=1: (n/2 -1) f 1=f 1+4*func 1(a+(2*i-1)*h)+2*func 1(a+(2*i)*h); end f 1=f 1+4*func 1(a+(n-1)*h)+func 1(b); f 1=(h/3)*f 1; simpson 1 test. m

%examples of use of function f 1=simpson 1(func 1, a, b, n) %Integration func 1 %using composite Simpson's rule at n points between a and b; n must be even a=1; b=5; f 1=inline('exp(-x)') disp(['a=', num 2 str(a), ' b=', num 2 str(b)]) %display limits of integration intexact=-exp(-b)+exp(-a); i=1: 10; nvec=2. ^(i); for i=1: length(nvec) f 1 int=simpson 1(f 1, a, b, nvec(i)); disp(['n=', num 2 str(nvec(i)), ' Error in composite trapezod rule =', num 2 str(f 1 int-intexact)]); %display value of function f(x) end …

>> simpson 1 test f 1(x) = exp(-x) a=1 b=5 n=2 Error in composite trapezod rule =0. 021369 n=4 Error in composite trapezod rule =0. 0017902 n=8 Error in composite trapezod rule =0. 00012176 n=16 Error in composite trapezod rule =7. 7793 e-006 n=32 Error in composite trapezod rule =4. 8892 e-007 n=64 Error in composite trapezod rule =3. 06 e-008 n=128 Error in composite trapezod rule =1. 9132 e-009 n=256 Error in composite trapezod rule =1. 1958 e-010 n=512 Error in composite trapezod rule =7. 4751 e-012 n=1024 Error in composite trapezod rule =4. 6674 e-013 The same function with composite trapezoid rule: a=1 b=5 n=2 Error in composite trapezod rule =0. 11305 n=4 Error in composite trapezod rule =0. 029605 n=8 Error in composite trapezod rule =0. 0074926 n=16 Error in composite trapezod rule =0. 001879 n=32 Error in composite trapezod rule =0. 00047011 n=64 Error in composite trapezod rule =0. 00011755 n=128 Error in composite trapezod rule =2. 9389 e-005 n=256 Error in composite trapezod rule =7. 3474 e-006 n=512 Error in composite trapezod rule =1. 8369 e-006 n=1024 Error in composite trapezod rule =4. 5922 e-007

Exercise Compute integral of function e^(-x^2) between a=0 and b=10 using composite trapezoid and composite Simpson’s rules for n=2, 4, 8, …, 1024. Compare with exact answer = sqrt(pi)/2 Use simpson 1 test. m as example.

>> inclass 26 f 1 = Inline function: f 1(x) = exp(-x^2) a=0 b=10 n=2 Error in composite trapezod rule =0. 78044 n=4 Error in composite trapezod rule =-0. 046459 n=8 Error in composite trapezod rule =-0. 1186 n=16 Error in composite trapezod rule =-0. 0010671 n=32 Error in composite trapezod rule =-6. 2875 e-012 n=64 Error in composite trapezod rule =0 n=128 Error in composite trapezod rule =2. 2204 e-016 n=256 Error in composite trapezod rule =5. 5511 e-016 n=512 Error in composite trapezod rule =0 n=1024 Error in composite trapezod rule =-1. 3323 e-015 The same function with composite trapezoid rule: a=0 b=10 n=2 Error in composite trapezod rule =1. 6138 n=4 Error in composite trapezod rule =0. 3686 n=8 Error in composite trapezod rule =0. 0032014 n=16 Error in composite trapezod rule =1. 8863 e-011 n=32 Error in composite trapezod rule =0 n=64 Error in composite trapezod rule =-1. 1102 e-016 n=128 Error in composite trapezod rule =-1. 1102 e-016 n=256 Error in composite trapezod rule =-4. 4409 e-016 n=512 Error in composite trapezod rule =-6. 6613 e-016 n=1024 Error in composite trapezod rule =-4. 4409 e-016