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

### 17074190作业9

5.1Exercises

(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.

5.1Computer Problems

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=tcf(f,x)
for i=1:12
h(i)=1/(10.^i);
end
for i=1:12
F(i)=(f(x+h(i))-f(x-h(i)))/(2*h(i));
e(i)=F(i)-1;
end
E(1:12,1)=h(1:12);
E(1:12,2)=F(1:12);
E(1:12,3)=e(1:12);
E=vpa(E,14);
x0=h(1:12);
y0=e(1:12);
plot(x0,y0)

Results:

>> f=@(x) sin(x)-cos(x);
>> E=tcf(f,0)

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]

5.2Exercises

(b)(cosx)(0,π/2)

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

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

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

error=1-π/40.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.

When m=1,h=. h*cos()=;

When m=2,h=. h*(cos()+cos())

When m=4,h= h*(cos()+cos()+cos()+cos())

(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.

5.2Computer Problems

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.

(d)

m=16;a=1;b=3;f=@(x) x^2*log(x);

h=(b-a)/(2*m);

T=f(a)+f(b);

for i=1:2*m-1

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

if i/2==0

T=T+2*f(n);

else T=T+4*f(n);

end

end

T=T*h/3

T =

9.1292

m=32;a=1;b=3;f=@(x) x^2*log(x);

h=(b-a)/(2*m);

T=f(a)+f(b);

for i=1:2*m-1

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

if i/2==0

T=T+2*f(n);

else T=T+4*f(n);

end

end

T=T*h/3

T =

9.2294

(h)

m=16;a=0;b=1;f=@(x) x/(x^4+1)^(1/2);

h=(b-a)/(2*m);

T=f(a)+f(b);

for i=1:2*m-1

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

if i/2==0

T=T+2*f(n);

else T=T+4*f(n);

end

end

T=T*h/3

T =

0.5801

m=32;a=0;b=1;f=@(x) x/(x^4+1)^(1/2);

h=(b-a)/(2*m);

T=f(a)+f(b);

for i=1:2*m-1

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

if i/2==0

T=T+2*f(n);

else T=T+4*f(n);

end

end

T=T*h/3

T =

0.5839

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

(a)

m=16;a=0;b=1;f=@(x) exp(x^2);

h=(b-a)/(2*m);

T=f(a)+f(b);

for i=1:2*m-1

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

if i/2==0

T=T+2*f(n);

else T=T+4*f(n);

end

end

T=T*h/3

T =

1.9121

m=32;a=0;b=1;f=@(x) exp(x^2);

h=(b-a)/(2*m);

T=f(a)+f(b);

for i=1:2*m-1

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

if i/2==0

T=T+2*f(n);

else T=T+4*f(n);

end

end

T=T*h/3

T =

1.9310

(f)

m=16;a=0;b=pi;f=@(x) cos(exp(x));

h=(b-a)/(2*m);

T=f(a)+f(b);

for i=1:2*m-1

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

if i/2==0

T=T+2*f(n);

else T=T+4*f(n);

end

end

T=T*h/3

T =

-0.4800

m=32;a=0;b=pi;f=@(x) cos(exp(x));

h=(b-a)/(2*m);

T=f(a)+f(b);

for i=1:2*m-1

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

if i/2==0

T=T+2*f(n);

else T=T+4*f(n);

end

end

T=T*h/3

T =

-0.4974

5.3 Exercises

(c)

h1=b-a=1;

h2=(1/2)*(b-a)=1/2;

h3=(1/4)*(b-a)=1/4;

R11=h1/[f(a)+f(b)]=e+1;

R21=h2/[f(a)+f(b)+2*f((a+b)/2))]=(e+1)/2+e^(1/2);

R31=(1/2)*R21+h3*=(1/4)*(e+1+2*e^(1/2)+e^(1/4)+e^(3/4));

R22=[(2^2)R21-R11]/3=(1/3)*[e+1+4*e^(1/2)]

R32=[(2^2)R31-R21]/3=(1/2)*[e+1+2*e^(1/2)]+e^(1/4)+e^(3/4);

R33=[(4^2)*R32-R22]/[(4^2)-1]=[23*(e+1)+44*e^(1/2)]/45+[16*(e^(1/4)+e^(3/4)]/15;

6.1 Exercises

3. Use separation of variables to nd 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′=5*(t^4)*y (e) y′=1/(y^2) (f) y′=t^3/y^2

（a） Separate variables, y′=dy/dt=t

dy=tdt

y=(1/2)*t^2+c, c is a constant.

Because of y(0)=1, c=1, to sum up y=(1/2)*t^2+1.

（b） Separate variables, y′=dy/dt= y′=(t^2)*y

dy/y=t^2dt

ln|y|=(1/3)*(t^3)+c, c is a constant.

|y|=e^[(1/3)*(t^3)+c]

Because of y(0)=1, c=0, to sum up y=e^[(1/3)*(t^3)].

（c）  Separate variables, y′=dy/dt= y′=2*(t + 1)*y

dy/y=2*(t+1)dt

ln|y|=t^2+2*t+c, c is a constant.

|y|=e^[t^2+2*t+c]

Because of y(0)=1, c=0, to sum up y=e^[t^2+2*t].

（d） Separate variables, y′=dy/dt= y′=5*(t^4)*y

dy/y=5*(t^4)dt

ln|y|=t^5+c, c is a constant.

|y|=e^[t^5+c]

Because of y(0)=1, c=0, to sum up y=e^[t^5].

e of y(0)=1, c=0, to sum up y=e^[t^2+2*t].

（e） Separate variables, y′=dy/dt= y′=1/(y^2)

y^2dy=dt

(1/3)y^3=t+c, c is a constant.

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

Because of y(0)=1, c=1/3, to sum up y=(3*t+1)^(1/3).

（f） Separate variables, y′=dy/dt= y′=t^3/y^2

y^2dy=t^3dt

(1/3)y^3=(1/4)*(t^4)+c, c is a constant.

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

Because of y(0)=1, c=1/3, to sum up y=[(3/4)*(t^4)+1]^(1/3).

(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=2.3583

(f)     w=1.0000    1.0000    1.0039    1.0349    1.1334

error=0.6166

6.2 Exercises

(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.5871

error=2.4129

(f) w=1.0000    1.0020    1.0193    1.0823    1.2182

error=0.5318