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

### 17074205作业九

P 252    5.1 Exercises

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.

f’(x)=(f(x+h)-f(x-h))/2h=(e^(x+h)-e^(x-h))/2h;

(a)f’(0)=(e^0.1-e^(-0.1))/2*0.1≈1.007

(b)f’(0)=(e^0.01-e^(-0.01))/2*0.01≈1.0000167

(c) f’(0)=(e^0.001-e^(-0.001))/2*0.001≈1.0000002.

P254    5.1 Computer Problems

1. Make a table of the error of the three-point centered-difference formula for f ′ (0), where f (x) = sin x − 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?

1Code

for n=-12:1:-1

h=10^(n)

f=vpa((sin(x+h)-cos(x+h)-sin(x-h)+cos(x-h))/(2*h),11)

error=vpa(1-f,11)

end

>> x=0;

>> Untitled

h =

1.0000e-12

f =

0.99997787828

error =

0.000022121720122

h =

1.0000e-11

f =

1.0000000827

error =

-0.000000082740370999

h =

1.0000e-10

f =

1.0000000827

error =

-0.000000082740370999

h =

1.0000e-09

f =

0.99999997172

error =

0.000000028281931574

h =

1.0000e-08

f =

1.000000005

error =

-0.0000000050247592753

h =

1.0000e-07

f =

0.99999999947

error =

0.0000000005263558478

h =

1.0000e-06

f =

1.0

error =

-0.000000000028755664516

h =

1.0000e-05

f =

0.99999999998

error =

0.000000000015653145446

h =

1.0000e-04

f =

0.99999999833

error =

0.0000000016665548941

h =

1.0000e-03

f =

0.99999983333

error =

0.00000016666667957

h =

0.0100

f =

0.99998333342

error =

0.000016666583335

h =

0.1000

f =

0.99833416647

error =

0.0016658335317

>> x0=[10^(-12) 10^(-11) 10^(-10) 10^(-9) 10^(-8) 10^(-7) 10^(-6) 10^(-5) 10^(-4) 10^(-3) 10^(-2) 10^(-1)];

>> y0=[.2212172e-4 -.827e-7 -.827e-7 .2828e-7 -.50e-8 .53e-9 0 .2e-10 .167e-8 .16667e-6 .1666658e-4 .166583353e-2];

>> c=newtdd(x0,y0,12);

>> x=10^(-12):10^(-1);

>> y=nest(11,c,x,x0);

>> plot(x0,y0,'o',x,y)

2

The minimum error closes to the theoretical expectation.

P263   5.2 Exercises

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)cos x dx  0，π/2

∫(cosx)(0,π/2)=1-0=1

(1)m=1

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

∫(cosx)(0,π/2)=(π/2)/2(y0+y1)=π/4(1+0)=π/4

error=1-π/4≈0.2146

(2)m=2

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

∫(cosx)(0,π/2)=(π/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

∫(cosx)(0,π/2)=(π/8)/2(y0+y4+2(y1+y2+y3))

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

≈0.9781

error=1-0.9781=0.0129.

2.Apply the Composite Midpoint Rule with m=1,2,and 4 panels to approximate the integrals in exercise 1,and report the errors.

∫(cosx)(0,π/2)=1-0=1

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

∫(cosx)(0,π/2)=π/2*f(π/2)=0

error=|1-0|=1

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

∫(cosx)(0,π/2)=(π/4)*(f(π/8)+ f(3π/8))= 1.0262

error=|1-1.0262|=0.0262

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

∫(cosx)(0,π/2)=(π/8)*(f(π/16)+f(3π/16)+f(5π/16)+f(7π/16))≈1.0065

error=|1-1.0065|=0.0065

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.

∫(cosx)(0,π/2)=1-0=1

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

∫(cosx)(0,π/2)=(π/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

∫(cosx)(0,π/2)=(π/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

∫(cosx)(0,π/2)=(π/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∫dx/(1+x^2)dx(0,1)

∫dx/(1+x^2)dx(0,1)= 1.557407724654902

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

∫dx/(1+x^2)dx(0,1)=1/6(y0+y2+4y1)=1/6(tan0+tan1+4tan1/2)= 0.623769614005011

error=|1.557407724654902-0.623769614005011|= 0.933638110649891

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

∫dx/(1+x^2)dx(0,1)= 1/12(y0+y4+4(y1+y3)+2y2)= 0.616480519083576

error= |1.557407724654902- 0.616480519083576|=0.940927205571326

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

∫dx/(1+x^2)dx(0,1)= 1/24(y0+y8+4(y1+y3+y5+y7)+2(y2+y4+y6))= 0.615693358232365

error= |1.557407724654902-0.615693358232365 |= 0.941714366422537

p264   5.2 Computer Problems

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.

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

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 is 0.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.

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

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

Sol：代码：

function simpson= simpson(a,b,m,f)

simpson=0;

d=0;

c=0;

for i=1:m-1

if (mod(i,2)==0)

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

else

d=d+f(a+(2*i-1)*(b-a)/2*m)

end

end

simpson=(b-a)/(6*m)*(f(a)+f(b)+2*c+4*d);

end

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

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

ans =

4.5194e+05

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

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

ans =

9.5482e+06

h>> f=@(x)x/((x^4+1)^(1/2));

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

ans =

0.0172

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

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

ans =

0.0065

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

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

The definite integral is 1.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

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

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

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

Sol:代码：

function simpson= simpson(a,b,m,f)

simpson=0;

d=0;

c=0;

for i=1:m-1

if (mod(i,2)==0)

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

else

d=d+f(a+(2*i-1)*(b-a)/2*m)

end

end

simpson=(b-a)/(6*m)*(f(a)+f(b)+2*c+4*d);

end

a

>> f=@(x)exp(x^2);

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

ans =

Inf

>> f=@(x)exp(x^2);

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

ans =

Inf

f>>  f=@(x)cos(exp(x));

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

d =

-0.2292

d =

NaN

ans =

NaN

>>  f=@(x)cos(exp(x));

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

d =

2.9098

d =

NaN

ans =

NaN

P268     5.3Exercises

1. Apply Romberg Integration to find R 33 for the integrals.

(c) ∫ e^x dx   (0,1)

R11=h1/2(f(0)+f(1))=(e+1)/2

R21=1/2R11+h2f((a+b)/2)=(e+1+2e^0.5)/4

R22=(4R21-R11)/3=(e+4e^0.5+1)/6

R32=1/2R21+0.25(f(0.25)+f(0.75))=(e+1+2e^0.5+4e^0.25+4e^0.75)/4

R33=(23e+23+44e^0.5+96e^0.25+96e^0.75)/90

P291 6.1Exercises

3. Use separation of variables to find solutions of the IVP given by y(0) = 1 and the following  differential equations:

(a) y = t         (b) y = t^ 2 y        (c) y = 2(t + 1)y

(d) y = 5t^ 4 y    (e) y = 1/y^ 2        (f) y = t^ 3 /y^2

(a)

dy/dt=t,dy=tdt,y=1/2t^2+C;

y(0)=0+C=1;

C=1;

y=1/2t^2+1.

(b)

(1/y)dy=t^2dt;

ln|y|=1/3t^3+C;

|y|=e^(1/3t^3+C);

y(0)=e^C=1,C=0;

y=e^1/3t^3

(c)

1/ydy=2(t+1)dt=(t+1)^2+C

|y|=e^((t+1)^2+C)

y(0)=e^(1+C)=1,C=-1

y= e^((t+1)^2-1)

dy/dt=5(t^4)y,(1/y)dy=(5t^4)dt;

ln|y|=t^5+C;

|y|=e^(t^5+C);

y(0)=e^C=1,C=0;

So,y=e^(t^5).

(e)

y^2dy=dt

1/3y^3=t+c

y^3=3t+3c

y=(3t+3c)^(1/3)

y(0)=1

c=1

y=(3t+3)^1/3

(f)

y^2dy=t^3dt

1/3y^3=1/4y^4+C

y=(3/4t^4+3c)^1/3

c=1/3

y=(3/4t^4+1)^1/3

5. Apply Euler’s Method with step size h = 1/4 to the IVPs in Exercise 3 on the interval [0,1]. List the w i ,i = 0,...,4, and find the error at t = 1 by comparing with the correct solution.

(a) y = t         (b) y = t^ 2 y        (c) y = 2(t + 1)y

(d) y = 5t^ 4 y    (e) y = 1/y^ 2        (f) y = t^ 3 /y^2

y(0) = 1

(a)

y=1/2t^2+1.

w0=y0=1,y0=1,                             e0=0

w1=w0+hf(t0,w0)=1+1/4f(0,1)=1,y1= 1.03125,    e1=0.03125

w2=w1+1/4f(t1,w1)=1.0625,y2=1.125,           e2=0.0625

w3=w2+1/4f(t2,w2)=1.1875,y3=1.28125,         e3=0.09375

w4=w3+1/4f(t3,w3)= 1.375,y4=1.5,              e4=0.125

Sol:Euler’s Method
w0=y0=1
wi+1=wi+hf(ti,wi)
(a) w=[1.0000,1.0000,1.0625,1.1875,1.3750],error=0.1250
(b) w=[1.0000,1.0000,1.0156,1.0791,1.2309],error=0.1648
(c) w=[1.0000,1.5000,2.4375,4.2656,7.9980]error=12.0875
(d) w=[1.0000,1.0000,1.0049,1.0834,1.5119],error=1.2064
(e) w=[1.0000,1.2500,1.4100,1.5357,1.6417],error=0.0543
(f) w=[1.0000,1.0000,1.0039,1.0349,1.1334],error=0.0717

P301   6.2   Exercises

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

approximation w 0 ,...,w 4 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 = t 2 y      (c) y = 2(t + 1)y

(d) y = 5t 4 y        (e) y = 1/y 2     (f) y = t 3 /y 2

Sol:w0=y0=1
wi+1=wi+h/2(f(ti,wi)+f(ti+h,hf(ti,wi)))
(a) w=[1.0000,1.0313,1.1250,1.2813,1.5000],error=0
(b) w=[1.0000,1.0078,1.0477,1.1587,1.4054],error=0.0097
(c) w=[1.0000,1.7188,3.3032,7.0710,16.7935],error=3.2920
(d) w=[1.0000,1.0024,1.0442,1.3077,2.7068],error=0.0115
(e) w=[1.0000,1.2050,1.3570,1.4810,1.5817],error=0.0003
(f) w=[1.0000,1.0020,1.0193,1.0823,1.2182],error=0.0132