﻿ 17074187作业9 - 计算数学达人 - 专，学者，数值代数，微分方程数值解

### 17074187作业9

2.Use the three-point centered-difference formula to approximate f'(0),where f(x)=e^x,for

(a)h=0.1 (b) h=0.01 (c) h=0.001.

(a) h=0.1,x=0,f(x)=e^x

f'(x)=[f(x+h)-f(x-h)]/2h=(e^0.1-e^-0.1)/0.2=0.2003/0.2=1.0016675

(b) h=0.01,x=0,f(x)=e^x

f'(x)=[f(x+h)-f(x-h)]/2h=(e^0.01-e^-0.01)/0.02=1.0000167

(c) h=0.001,x=0,f(x)=e^x

f'(x)=[f(x+h)-f(x-h)]/2h=(e^0.001-e^-0.001)/0.002=1.00000015

Computer probleam

1.Make a table of the error of the three-point centered-difference formula for f'(0),where f(x)=sinx-cosx,

with h=10^-1,...,10^-12,as in the table in Section 5.1.2.Draw a plot of the results.Does the minimum error

correspond to the theoretical expectation?

Code:

function E=tcd(f,x)

for i=1:12

h(i)=1/(10.^i);

end

for i=1:12

H(i)=(f(x+h(i))-f(x-h(i)))/(2*h(i));

end

for i=1:12

e(i)=1-H(i);

end

E(1:12,1)=h(1:12);

E(1:12,2)=H(1:12);

E(1:12,3)=e(1:12);

E=vpa(E,14);

x0=h(1:12);

y0=e(1:12);

plot(x0,y0)

Run

>> f=@(x) sin(x)-cos(x);

>> E=tcd(f,0)

Result

E =

[            0.1, 0.99833416646828,         0.0016658335317179]

[           0.01, 0.99998333341667,       0.000016666583334768]

[          0.001, 0.99999983333338,      0.0000001666666240574]

[         0.0001, 0.99999999833289,   0.0000000016671100055987]

[        0.00001, 0.99999999998435, 0.000000000015653367491097]

[       0.000001, 0.99999999997324, 0.000000000026755486715047]

[      0.0000001, 0.99999999947364,  0.00000000052635584779637]

[     0.00000001, 0.99999999947364,  0.00000000052635584779637]

[    0.000000001,  1.0000000272292,   -0.000000027229219767833]

[   0.0000000001,  1.0000000827404,    -0.00000008274037099909]

[  0.00000000001,  1.0000000827404,    -0.00000008274037099909]

[ 0.000000000001,  1.0000333894311,      -0.000033389431109754]

1. Apply the composite Trapezoid Rule with m=1,2,and 4 panels to approximate the integral.Compute the error by comparing with the exact value from calculus.

b)∫（0  π/2）cosx dx

(1)m=1,h=(b-a)/m=(π/2-0)/1=π/2

=(π/2)/2*(y0+y1)=π/4(1+0)=π/4

error=1-π/40.2146

(2)m=2, h=(b-a)/m=(π/2-0)/2=π/4

=(π/4)/2*(y0+y2+2y1)=π/8(1+0+2^(1/2))0.9481

error=1-0.9481=0.0519

(3)m=4, h=(b-a)/m=(π/2-0)/4=π/8

=(π/8)/2*(y0+y4+2(y1+y2+y3))

=π/16*[1+0+2(cos(π/8)+ cos(2π/8)+ cos(3π/8))]

0.9871

error=1-0.9781=0.0129.

3. Apply the composite Simpson’s Rule with m=1,2,and 4 panels to approximate the integrals in exercise 1,and report the errors.

(1)m=1,h=(b-a)/2m=(π/2-0)/2=π/4

=(π/4)/3*(y0+y2+4y1)=π/12(1+0+4*(2)^(1/2)/2)1.0022799

error=|1-1.0022799|=0.0022799

(2)m=2, h=(b-a)/2m=(π/2-0)/4=π/8

=(π/8)/3*(y0+y4+4(y1+y3)+2y2)≈1.0001346

error=|1-1.0001346|=0.0001346.

(3)m=4, h=(b-a)/2m=(π/2-0)/8=π/16

=(π/16)/3*(y0+y8+4(y1+y3+y5+y7)+2(y2+y4+y6))1.0000083

error=|1-1.0000083|=0.0000083.

4. Apply the composite Simpson’s Rule with m = 1,2, and 4 panels to the integrals, and report

the errors.(b)f(1  0)dx/(1+x^2)dx

1. Use the composite Trapezoid Rule with m=16 and 32 panels to approximate the definite integral.Compare with the correct integral and report the two errors.

2. Apply the composite Simpson’s Rule to the integrals in Computer Problem 1. Use m = 16 and

(d)∫(1  3) x^2ln x dx    (h)∫(0 1)(x dx)/(x^4+1)^(1/2)

(d)(1)Answer:The definite integral is 7.0098 and 7.0014.

(2)Code:function trapezoid = trapezoid(a,b,m,f)

trapezoid=0;

for i=1:m-1

trapezoid=trapezoid+f(a+i*(b-a)/m);

end

trapezoid=(b-a)/(2*m)*(f(a)+f(b)+2*trapezoid);

end

(3)Result:>> f=@(x)    (x^2)*(log(x));

>>trapezoid(1,3,16,f)

ans =

7.0098

>> f=@(x)    (x^2)*(log(x));

>>trapezoid(1,3,32,f)

ans =

7.0014

m=16,error=9ln3-26/9-7.0098;

m=32,error=9ln3-26/9-7.0014.

The definite integral is0.4404 and 0.4406.

(2)Code:function trapezoid = trapezoid(a,b,m,f)

trapezoid=0;

for i=1:m-1

trapezoid=trapezoid+f(a+i*(b-a)/m);

end

trapezoid=(b-a)/(2*m)*(f(a)+f(b)+2*trapezoid);

end

(3)Result:>> f=@(x)  x/((x^4+1)^(1/2));

>>trapezoid(0,1,16,f)

ans =

0.4404

>>trapezoid(0,1,32,f)

ans =

0.4406

m=16,error=1/2(ln(2^(1/2))+1)-0.4404

m=32,error=1/2(ln(2^(1/2))+1)-0.4406.

3Use the composite Trapezoid Rule with m=16 and 32 panels to approximate the definite integral.

4Apply the composite Simpson’s Rule to the integrals of Computer Problem 3, using m = 16

and 32

(a)∫(0 1) (e^x^2) dx   (f)∫(0 π) cos e^x dx

(a)(1)Answer:The definite integral is1.4644 and 1.4631.

(2)Code:function trapezoid = trapezoid(a,b,m,f)

trapezoid=0;

for i=1:m-1

trapezoid=trapezoid+f(a+i*(b-a)/m);

end

trapezoid=(b-a)/(2*m)*(f(a)+f(b)+2*trapezoid);

end

(3)Result:>> f=@(x)  exp(x^2);

>>trapezoid(0,1,16,f)

ans =

1.4644

>>trapezoid(0,1,32,f)

ans =

1.4631

(f)(1)Answer:The definite integral is-0.2780 and -0.3568.

(2)Code:function trapezoid = trapezoid(a,b,m,f)

trapezoid=0;

for i=1:m-1

trapezoid=trapezoid+f(a+i*(b-a)/m);

end

trapezoid=(b-a)/(2*m)*(f(a)+f(b)+2*trapezoid);

end

(3)Result:>> f=@(x)   cos(exp(x));

>>trapezoid(0,pi,16,f)

ans =

-0.2780

>>trapezoid(0,pi,32,f)

ans=

-0.3568

3.

Use separation of variables to fifind solutions of the IVP given by y(0) = 1 and the following

differential equations:

(a) y′ = t (b) y′ = t2y (c) y′ = 2(t + 1)y(d) y′ = 5t4y (e) y′ = 1/y2(f) y′ = t3/y2

(a)dy/dt=t,dy=t*dt,y=1/2*t^2+c,y(0)=1,c=1
y=1/2t^2+1

(b)dy/y=t^2dt,ln|y|=1/3t^3+c,y(0)=1,c=2/3

Ln|y|=1/3t^3+2/3

(c)ln|y|=t^2+2t+c,y(0)=1,c=-2

Ln|y|=t^2+2t-2
(d)dy/y=5*t^4*dt,ln|y|=t^5+c,|y|=e^(t^5+c),y(0)=1,c=0
y=e^(t^5)

(e)dyy^2=dt,y^3=3t+c,y(0)=1,c=-3

Y^3=3t-3

(f)dyy^2=t^3dt,y^3=3/4t^4+c,y(0)=1,c=-3/4

Y^3=3/4t^4-3/4

5. Apply Euler’s Method with step size h = 1/4 to the IVPs in Exercise 3 on the interval

[0,1]. List the wi, i = 0,... ,4, and fifind the error at t = 1 by comparing with the correct

solution.

w0 = y0

wi+1 = wi + hf (ti,wi).

1. Using initial condition y(0) = 1 and step size h = 1/4, calculate the Trapezoid Method

approximation w0,... ,w4 on the interval [0,1]. Find the error at t = 1 by comparing with the

correct solution found in Exercise 6.1.3.

(a) y′ = t (b) y′ = t2y (c) y′ = 2(t + 1)y

(d) y′ = 5t4y (e) y′ = 1/y2(f) y′ = t3/y2

(a)dy/dt=t,dy=t*dt,y=1/2*t^2+c,y(0)=1,c=1
y=1/2t^2+1

(b)dy/y=t^2dt,ln|y|=1/3t^3+c,y(0)=1,c=2/3

Ln|y|=1/3t^3+2/3

(c)ln|y|=t^2+2t+c,y(0)=1,c=-2

Ln|y|=t^2+2t-2
(d)dy/y=5*t^4*dt,ln|y|=t^5+c,|y|=e^(t^5+c),y(0)=1,c=0
y=e^(t^5)

(e)dyy^2=dt,y^3=3t+c,y(0)=1,c=-3

Y^3=3t-3

(f)dyy^2=t^3dt,y^3=3/4t^4+c,y(0)=1,c=-3/4

Y^3=3/4t^4-3/4