当前位置:文档之家› 数值分析作业

数值分析作业

数值分析作业
数值分析作业

数值分析作业及答案

Chap1

1、写出下列语句的运行结果。在MA TLAB 上执行它们以验证所得解答。

a=[1 2 3 ;4 5 6 ]’ b=[9;7;5;3;1] c=b(2:4) d=b(4:-1:1) e=sort(b) f=[3,b ’]

解:a=6352

4

1 b=13579

c=357

d=9

7

5

3 e=9

7

53

1

F=[3 9 7 5 3 1] 3、给定一向量:

a=[4 -1 2 -8 4 5 -3 -1 6 -7]

写一段程序计算a 中正数的和。运行程序并显示结果。 解:a=[4 -1 2 -8 4 5 -3 -1 6 -7]; s=0;

for i=1:length(a) if a(i)>0

s=s+a(i); end end s

6、编写一个函数M 文件fun_es(x),计算如下函数:

23

0.5sin x y e x x =-

其中参数可以为标量,也可以为向量。在MA TLAB 里键入如下命令检验此函数:

fun_es(3) fun_es([1 2 3])

解:function y=fun_es(x) y=0.5*exp(x/3)-x.^2.*sin(x);

chap2

1、设0x >,x 的相对误差为δ,求L nx 的误差。

解:Lnx-Lnx*=dLnx=dx/x=δ

2、设x 的相对误差为2%,求2x 的相对误差。 解:dLnf(x)=xf ’(x)/f(x)dLnx=4%

5、计算球体积要使相对误差限为1%,问度量半径R 时允许的相对误差限是多少?

解:dLnf(x)=xf ’(x)/f(x)dLnx=3dLnx=1% dLnx=0.33%

9、正方形的边长大约为100cm ,应怎样测量才能使其面积误差不超过1cm 2? 解:s=x 2

s-s*=2x(x-x*)=1

x-x*=1/(2x)=1/200=0.5*10-2 即测量边的误差不超过0.005cm 10、设2

12

S gt =

,假定g 是准确的,而对t 的测量有±0.1秒的误差,证明当t 增加时S 的绝对误差增加,而相对误差却减少。

解:s 的绝对误差为: gtdt 当t 增加时,gtdt 增大。

S 的相对误差为:tgt/(gt 2/2)dLnt=2/tdt 当t 增加时,2/tdt 变小。

Chap3

1.当的二次插值多项式求时,)(,4,3,0)(2,1,1x f x f x -=-=

解:Φ(x)=l0(x)y0+l1(x)y1+l2(x)y2

=0+(x-1)(x-2)/((-1-1)(-1-2))*(-3)+(x+1)(x-1)/((2-1)(2+1)*4 =(5x 2

+9x-14)/6

2.求不大于3次的插值多项式,此多项式在7

,6,3,13210

====x x x x 处与的值

x

x f 2

)(=

相同

解:y0=1 y1=9 y2=36 y3=49

Φ(x)=l0(x)y0+l1(x)y1+l2(x)y2+l3(x)y3

=(x-3)(x-6)(x-7)/(-60)*1+(x-1)(x-6)(x-7)/24*9+(x-1)(x-3)(x-7)/(-15)*36+ (x-1)(x-3)(x-6)/24*49

试构造Hermite 插值多项式)(3x H 并计算)5.1(f 的值

解:a0(x)=(1+2*(x-1)/(2-1))((x-2)/(1-2))2=(2x-1)(x-2)2

a1(x)=(1+2*(x-2)/(1-2))((x-1)/(2-1))2=(2x-1)(x-2)2=(5-2x)(x-1)2 b0(x)=(x-1)((x-2)/(1-2))2=(x-1)(x-2)2 b1(x)=(x-2)((x-1)/(2-1))2=(x-2)(x-1)2 H3(x)=2a0(x)+3a1(x)-b1(x)

Chap4

3.用复合梯形求积法计算下列积分,取n=2,4,8,16:

?

??+1

01

/4

21

)()exp( )()tan( )(dx x

c dx x b dx

x a π

解:(a ) I2=h/2(f1+2f2+f3)=pi/16(tan(0)+2tan(pi/8)+tan(pi/4))=0.35901 I4=h/2(f1+2f2+2f3+2f4+f5)=0.34976

I8=h/2(f1+2(f2+f3+f4+f5+f6+f7+f8)+f9)=0.34737

I16=h/2(f1+2(f2+f3+f4+f5+f6+f7+f8+f9+f10+f11+f12+f13+f14+f15+f16)+f17)=0.34677 (b) I2=h/2(f1+2f2+f3)=1/4(exp(0)+2*exp(0.5)+exp(1))=1.75393 I4=h/2(f1+2f2+2f3+2f4+f5)=1.72722

I8=h/2(f1+2(f2+f3+f4+f5+f6+f7+f8)+f9)=1.72052

I16=h/2(f1+2(f2+f3+f4+f5+f6+f7+f8+f9+f10+f11+f12+f13+f14+f15+f16)+f17)=1.71884

(c) I2=h/2(f1+2f2+f3)=1/4(exp(0)+2*exp(0.5)+exp(1))=0.40883 I4=h/2(f1+2f2+2f3+2f4+f5)=0.40619

I8=h/2(f1+2(f2+f3+f4+f5+f6+f7+f8)+f9)= 0.40565

I16=h/2(f1+2(f2+f3+f4+f5+f6+f7+f8+f9+f10+f11+f12+f13+f14+f15+f16)+f17)= 0.40551

4.用复合梯形求积法计算如下积分,取n=2,4,8,25,100:

?

2

/0

)sin(πdx x

并给出数值积分结果与精确结果1之间的误差。

解:clear a=0;b=pi/2; n=1;

for k=1:5,n=2*n; if k==4,n=25; end

if k==5,n=100; end

h=(b-a)/n;i=1:n+1; x=a+(i-1)*h;f=sin(x);

I=(h/2)*(f(1)+2*sum(f(2:n))+f(n+1));

fprintf(' %3.0f %10.5f %10.5f\n',n,I,I-1); end 结果如下:

2 0.94806 -0.05194 4 0.98712 -0.01288 8 0.99679 -0.00321 25 0.99967 -0.0003

3 100 0.99998 -0.00002

5. 将区间[1,9]8等分,试用复化梯形公式求积分x

x d 5691

?-的近似值,计算过程中保留3

位小数.

解:clear a=1;b=9; n=8;

h=(b-a)/n;i=1:n+1; x=a+(i-1)*h;f=sqrt(6*x-5);

I=(h/2)*(f(1)+2*sum(f(2:n))+f(n+1));

fprintf(' %3.0f %10.5f \n',n,I); 运行结果如下: 8 37.81817

8.已知勒让德多项式n n

n n

n dx x d n x P )1(!21)(2-=,证明两点(n=2)的高斯-勒让德求积公式

为)3

1(

)3

1(d )(1

1f f x x f +-

≈?-

证明:p 2(x)=3/2x 2-1/2=0 a1=a2=1

X1=3

1

x2=-

3

1

)3

1

(

)31

()2(2)1(1d )(1

1

f f x f a x f a x x f +-

=+≈?

-

10.用6个节点的高斯求积公式计算下列积分:

?????

++2

5

.01

1

210)log( )()log()sin( )()2exp( )()1log( )()cos(21

)(dx

x e dx

x x d dx

x x c dx

x x b dx

x a π

解:(a)令x=pi/2(t+1)

?

cos 21dx x =?

-++1

1)])

1(2/cos[2(2dt t ππ

=A1f(t1)+A2f(t2)+A3f(t3)+A4f(t4)+A5f(t5)+A6f(t6)

=0.17132449*(pi/(4+2*cos(pi/2*(-0.93246951+1))))+0.36076157*(pi/(4+2*cos(pi/2*(-0.6612093

9+1))))+ 0.46791393*(pi/(4+2*cos(pi/2*(-0.23861919+1))))+ 0.46791393*(pi/(4+2*cos(pi/2*(0.23861919+1))))+ 0.36076157*(pi/(4+2*cos(pi/2*(0.66120939+1))))+ 0.17132449*(pi/(4+2*cos(pi/2*(0.93246951+1)))) =1.8138

(b)令x=3/2+t/2

?

+21

)

1log(dx x x =?

-++11

3)

2/2/5log(dt t t

=A1f(t1)+A2f(t2)+A3f(t3)+A4f(t4)+A5f(t5)+A6f(t6)

=0.17132449*(log(2.5-0.93246951/2)/(3-0.93246951))+0.36076157*(log(2.5-0.66120939/2)/(3-0.66120939))+ 0.46791393*(log(2.5-0.23861919/2)/(3-0.23861919)) + 0.46791393*(log(2.5+0.23861919/2)/(3+0.23861919)) + 0.36076157*(log(2.5+0.66120939/2)/(3+0.66120939)) + 0.17132449*(log(2.5+0.93246951/2)/(3+0.93246951)) =0.6143

(c)令x=1/2+t/2

?

1

2dx xe x

=

?

-++1

114

1)1(dt e t t

=A1f(t1)+A2f(t2)+A3f(t3)+A4f(t4)+A5f(t5)+A6f(t6)=2.0973 (d)令x=1/2+t/2

?1

)log(sin dx x x =?

-++1

121

212

1)log()sin(dt

t t

=A1f(t1)+A2f(t2)+A3f(t3)+A4f(t4)+A5f(t5)+A6f(t6)=-0.24

(e)令x=5/4+3*t/4

?

2

5

.0)log(dx x =?-+1

145

34

3)log(dt t

=A1f(t1)+A2f(t2)+A3f(t3)+A4f(t4)+A5f(t5)+A6f(t6)=0.2329

19.用三点公式计算2

1+1

=

)

()(x x f 在x =1.0,1.1,1.2处的导数值.已知函数值 f (1.0)=0.250 000, f (1.1)=0.226757, f (1.2)=0.206 612 解:h=0.1,f k ’=(f k+1-f k-1)/2h

f(1.0)’=(f(1.1)-f(0.9))/0.2=-0.2513 f(1.1)’=(f(1.2)-f(1.0))/0.2=-0.2163 f(1.2)’=(f(1.3)-f(1.1))/0.2=-0.1886

chap6

2. 求方程在=1.5附近的一个根,将方程改写成下列等价形式,并建立相应迭代公式.

(1) ,迭代公式.

(2) ,迭代公式.

(3),迭代公式.

试分析每种迭代公式的收敛性,并选取一种收敛最快的方法求具有4位有效数字的近似根. 解:(1)1+1/x 2的一阶导数绝对值<1,故迭代收敛 X1=1.5, xk+1=1+1/xk^2, 得到x=1.4656 (2) (1+x 2)^1/3的一阶导数绝对值>1,故迭代发散 (3) 1/(x-1)^1/2的一阶导数绝对值>1,故迭代发散 14.用牛顿下山法求方程01)(3

=--=

x x f x

的根,取,6.00=x 计算精确到10

5

-。

解:迭代公式为:

1

31)

()

(123---+-=-=k k k k k x x x k

k x df x f k k k p x p

x x

X(1)=0.6 p(1)=1 p(k)=p(k-1)/2 x=1.3395

Chap7

解:正规方程为: 5a0+10a1=20.1 10a0+22.5a1=44.95

解方程得a0=0.22 a1=1.9 拟合直线为s(x)=0.22+1.9x

正规方程为:

6a0+8.9a1+14.17a2=14.58

8.9a0+14.17a1+24.023a2=22.808

14.17a0+24.023a1+42.8629a2=38.0962

解方程得a0=0.5966 a1=1.2533 a2=-0.0109

S2(x)=0.5966+1.2533x-0.0109 x2

设s3(x)=a0+a1x+a2x2+a3x3

正规方程为:

6a0+8.9a1+14.17a2+24.023a3=14.58

8.9a0+14.17a1+24.023a2+42.8629a3=22.808

14.17a0+24.023a1+42.8629a2+79.5192a3=38.0962

24.023a0+42.8629a1+79.5192a2+151.801a3=67.1883

解方程得a0=0.629 a1=1.185 a2=-0.0353 a3=-0.01

S3(x)= 0.629+1.185x-0.0353 x2-0.01 x3

5.用正交二次和三次多项式求解上题。

解:设s2(x)=a1t1(x)+a2t2(x)+a3t3(x)

t1(x)=1

t2(x)=x-sum(x)/6=x-1.4833

t3(x)=x2-sum(x.^2.*t2(x))/sum(t2(x.^2))*t2(x)-sum(x.^2)/sum(t1(x.^2))*t1(x) =x2-3.103x+2.241

a1=(s,t1)/(t1,t1)=sum(y*t1)/6=2.43

a2=(s,t2)/(t2,t2)=sum(y*t2)/sum(t2*t2)=1.2201

a3=(s,t3)/(t3,t3)=sum(y*t3)/sum(t3,t3)=-0.0416

于是,正交二次多项式为:

s2(x)=a1t1(x)+a2t2(x)+a3t3(x)

=2.43+1.2201(x-1.4833)-0.0416(x2-3.103x+2.241)

=0.5270+1.3492x-0.0416x2

设s3(x)=a1t1(x)+a2t2(x)+a3t3(x)+a4t4(x)

t1(x)=1

t2(x)=x-sum(x)/6=x-1.4833

t3(x)=x2-sum(x.^2.*t2(x))/sum(t2(x.^2))*t2(x)-sum(x.^2)/sum(t1(x.^2))*t1(x) =x2-3.103x+2.241

t4(x)=x3-(x4,t3)/(t3,t3)*t3(x)-(x4,t2)/(t2,t2)*t2(x)-(x4,t1)/(t1,t1)*t1(x) =x3-14.2273(x2-3.103x+2.241)-16.4619(x-1.4833)-7.1438

=x3-14.2273x2+60.6092x-14.6093

a1=(s,t1)/(t1,t1)=sum(y*t1)/6=2.43

a2=(s,t2)/(t2,t2)=sum(y*t2)/sum(t2*t2)=1.2201

a3=(s,t3)/(t3,t3)=sum(y*t3)/sum(t3,t3)=-0.0416

a4=(s,t4)/(t4,t4)=0.0531

于是,正交三次多项式为:

S3(x)=a1t1(x)+a2t2(x)+a3t3(x)+a4t4(x)

=2.43+1.2201(x-1.4833)-0.0416(x2-3.103x+2.241)

+0.0531(x3-14.2273x2+60.6092x-14.6093)

=0.0531x3-0.7971 x2+4.5675x-0.2488

18.求函数x x x f 2sin 3

1

cos 21)(+=

在区间]1,0[的一次最佳平方逼近多项式 解:设s(x)=a+bx , 正规方程为:

a+1/2b=int(‘f(x)’,’x ’,0,1)=0.6568

1/2a+1/3b=int(‘x*f(x)’,’x ’,0,1)=0.3360 S(x)=a+bx=0.6112+0.0912

19.求函数3)(x x f =在区间]2,0[的正交二次最佳平方逼近多项式 解:设s(x)=a1t1(x)+a2t2(x)+a3t3(x) t1(x)=1

t2(x)=x-(x,t1)/(t1,t1)*t1(x)=x-int(‘x’,’x ’,0,2)/2=x-1 t3(x)=x 2-(x 2,t2)/(t2,t2)*t2(x)-(x 2,t1)/(t1,t1)*t1(x)

=x 2-4(x-1)+4/3=x 2-4x+16/3 于是得到正规方程:

(t1,t1)a1+(t1,t2)a2+(t1,t3)a3=int(‘t1*f(x)’,’x ’,0,2) (t2,t1)a1+(t2,t2)a2+(t2,t3)a3=int(‘t2*f(x)’,’x ’,0,2) (t3,t1)a1+(t3,t2)a2+(t3,t3)a3=i nt(‘t3*f(x)’,’x ’,0,2)

2a1+0a2+16/3a3=4 0a1+2/3a2-4/3a3=12/5

16/3a1-4/3a2+17.0667a3=6.4

a1=-5.9985 a2=9.5989 a3=2.9994

s(x)=-5.9985+9.5989(x-1)+2.9994(x 2-4x+16/3)

20.求函数x

e x

f =)(在区间]1,1[-上的一次最佳一致逼近多项式 解:设s(x)=ax+b

R(x)=f(x)-s(x)=e x -ax-b 交错点组为x1,x2,x3 F’’(x)在[-1,1]内恒正, x1=-1, x3=1 R(x1)=-r(x2)=r(x3) , r ’(x2)=0 e -1+a-b=e 1-a-b

e -1+a-b=-e x2+ax2+b e x2-a=0

解方程组得:a=(e 2-1)/2e =1.1752

x2=log(e2-1)-log(2e)=0.1614 b=(2a-ax2+1/e)/2=1.2643 一次最佳一致逼近为s(x)=1.1752x+1.2643

Chap8

1. 用高斯消去法解线性方程组

12341234

124123423642710642624525

x x x x x x x x x x x x x x x +++=??++-=??

++=??-+++=-?

解:一、二式不能用gauss 消去法,改用列元素消去法:

6 4 0 2 6 4 2

7 -1 10 2 1 3 1 6

-2 4 5 2 -5

1.0000 0.6667 0 0.3333 1

0 -0.6667 7.0000 -2.3333 6 0 -0.3333 3.0000 0.3333 4 0 5.3333 5.0000 2.6667 -3

1.0000 0.6667 0 0.3333 1

0 5.3333 5.0000 2.6667 -3 0 -0.3333 3.0000 0.3333 4 0 -0.6667 7.0000 -2.3333 6

1.0000 0.6667 0 0.3333 1

0 1.0000 0.9375 0.5000 -0.5625 0 0 3.3125 0.5000 3.8125 0 0 7.6250 -2.0000 5.6250

1.0000 0.6667 0 0.3333 1

0 1.0000 0.9375 0.5000 -0.5625 0 0 7.6250 -2.0000 5.6250 0 0 3.3125 0.5000 3.8125

1.0000 0.6667 0 0.3333 1

0 1.0000 0.9375 0.5000 -0.5625 0 0 1.0000 -0.2623 0.7377 0 0 0 1.3689 1.3689

1.0000 0.6667 0 0.3333 1

0 1.0000 0.9375 0.5000 -0.5625 0 0 1.0000 -0.2623 0.7377 0 0 0 1.0000 1 X=( 2 -2 1 1)

2.用高斯列主元消去法解线性方程组

1231215342421042x x x -????????????=????????????-??????

解:结果如下:

3 4 2 4

1 -

2 1 5 2 10 4 -2

1.0000 1.3333 0.6667 1.3333 0 -3.3333 0.3333 3.6667 0 7.3333

2.6667 -4.6667

1.0000 1.3333 0.6667 1.3333 0 7.3333

2.6667 -4.6667 0 -

3.3333 0.3333 3.6667

1.0000 1.3333 0.6667 1.3333 0 1.0000 0.3636 -0.6364 0 0 1.5455 1.5455

1.0000 1.3333 0.6667 1.3333 0 1.0000 0.3636 -0.6364 0 0 1.5455 1.5455

1.0000 1.3333 0.6667 1.3333 0 1.0000 0.3636 -0.6364 0 0 1.0000 1.0000 X=(

2.0000 -1.0000 1.0000)’ 4.应用直接三角分解法解方程组

???

??=-+-=-+-=+-1

2322534321

321321x x x x x x x x x 解:L=

1.0000 0 0 -0.2500 1.0000 0

0.5000 2.0000 1.0000 U=

4.0000 -3.0000 1.0000 0 1.2500 -1.7500 0 0 2.0000

Ly=b y=( 5.0000 -1.7500 2.0000) Ux=y x=(1 0 1) 10.讨论方程组

123121221361313x x x -??????

??????=????????????-??????

的Jacobi 迭代法的收敛性

解:B=I-D -1A=0

3

1

302

1

20

---- ,它的特征值分别为:-3.2 -0.9 4.1。 最大特征值的模=4.1>1,Jacobi 迭代法发散。 11.应用Jacobi 迭代法解方程组

121232

31091027106

x x x x x x x -=??

-+-=??-+=? 解:X (k+1)=BX (k)+F= (I-D -1A)X (k)+D -1b x (0)=[0 0 0]; x (1)=[ 0.9000 0.7000 0.6000]; x (2)=[ 0.9700 0.9100 0.6700]; x (3)=[ 0.9910 0.9310 0.6910]; x 4)=[ 0.9931 0.9373 0.6931]; x (5)=[ 0.9938 0.9381 0.6938];

12. 应用Gauss-Seidel 迭代法解线性方程组

121232

31091027106

x x x x x x x -=??

-+-=??-+=? 解:x (0)=[0 0 0];

x (1)=[ 0.9000 0.7900 0.6790];

x (2)=[ 0.9934 0.9380 0.6938]; x (3)=[ 0.9938 0.9381 0.6938]; x (4)=[ 0.9938 0.9381 0.6938]; x (5)=[ 0.9938 0.9381 0.6938];

Chap9

1. 用幂法计算下列矩阵的主特征值及对应的特征向量:

(a) 1732341213A -?? ?=- ? ?--??,(b )2343463331A -?? ?

=- ? ???

当特征值有3位小数稳定时迭代终止。 解:(a )[i,j]=power_method(a1,x,1e-6,30) i=9.6056;

j=[-1,-0.6056,0.3944]’

(b) [i,j]=power_method(a2,x,1e-6,30) i=8.8699;

j=[-0.6044,1.0,0.1508]’

2、利用反幂法求矩阵

621231111?? ? ? ???

的最接近于6的特征值及对应的特征向量。 解:a=[6 2 1;2 3 1;1 1 1] a =

6 2 1 2 3 1 1 1 1 x=[1 1 1]' x =

1 1 1

[i,j]=INV_shift(a,x,1e-6,20)

Chap10

5、对初值问题

2

2

12 , 011(0)0

y y x x y ?'=-≤≤?+?

?=? 试用Euler 法取步长h=0.1和h=0.2计算其近似解,用Matlab 实现求解。

解:dyfun=inline(‘1/(1+x^2)-2*y^2’)

[x,y]=naeuler(dyfun,[0,1],0,0.1);[x,y] ans =

0 0 0.1000 0.1000 0.2000 0.1970 0.3000 0.2854 0.4000 0.3609 0.5000 0.4210 0.6000 0.4656 0.7000 0.4957 0.8000 0.5137

0.9000 0.5219

1.0000 0.5227

[x,y]=naeuler(dyfun,[0,1],0,0.2);[x,y]

ans =

0 0

0.2000 0.2000

0.4000 0.3763

0.6000 0.4921

0.8000 0.5423

1.0000 0.5466

12、用变形的Euler方法、二阶Heun方法解初值问题

'=-++≤≤=

1,01,(0)1

y y t t y

并Matlab实现。

解:变形的Euler方法:

dyfun=inline('-y+x+1')

[ x,y1]=naeuler1(dyfun,[0,1],1,0.1);[x,y1]

ans =

0 1.0000

0.1000 1.0050

0.2000 1.0190

0.3000 1.0412

0.4000 1.0708

0.5000 1.1071

0.6000 1.1494

0.7000 1.1972

0.8000 1.2500

0.9000 1.3072

1.0000 1.3685

二阶Heun方法:

[x,y3]=naeuler3(dyfun,[0,1],1,0.1);[x,y3]

ans =

0 1.0000

0.1000 1.0050 0.2000 1.0190 0.3000 1.0412 0.4000 1.0708 0.5000 1.1071 0.6000 1.1494 0.7000 1.1972 0.8000 1.2500 0.9000 1.3072 1.0000 1.3685

13、应用经典的四阶Runge-Kutta 方法解初值问题

21

(),13,(1)2y y y t y t

'=+≤≤=-

取步长为h=0.5

解:dyfun=inline('(y^2+y)/x')

x,y4]=naeuler4(dyfun,[1,3],-2,0.5);[x,y4] ans =

1.0000 -

2.0000 1.5000 -1.4954 2.0000 -1.3306 2.5000 -1.2480

3.0000 -1.1985

Chap11

2.将下列线性规划问题化为标准型 (1)min z =3 x 1+5x 2-4x 3+2x 4

???????≥=--+-≥-+-≤+-+0

,,95341322318362:.421432143214321x x x x x x x x x x x x x x x t s

解:化为标准型如下:

min z =3 x 1+5x 2-4)(87x x -+2x 4

???????≥=---+-=---+-=++--+0

8,7,6,5,,,95)(34132)(23183)(62:.42148721648721548721x x x x x x x x x x x x x x x x x x x x x x x x t s

(3)max 1234f 342x x x x =+++

??????

?≤≥=+-+-≥-+--≤--+0

,0,152342722351232:.42143214

3214321x x x x x x x x x x x x x x x t s 解:化为标准型如下:

min –f=-3x 1-x 2-4(x 5-x 6)+2x 7

??????

?≥=-+-=-++--=++-+0

,,,,,,152-)x -(x 3427)x -(x 223512)x -(x 32:.9876521765219

76521876521x x x x x x x x x x x x x x x x x x t s 3. 用单纯型法求解以下线性规划问题

(1)max 232z x x =- 12323123

3412.:212,,0x x x s t x x x x x ++=??

-≤??≥?

解:LP 化为标准形式如下: Min -z=0x 1-x 2+2x 3

X 1+3x 2+4x 3=12 2x 1-x 3+x 4=12 X 1,x 2,x 3,x 4≥0 单纯形表如下:

0 -1 2 0

X1 x2 x3 x4 b θ 0 x1 1 3 4 0 12 4 0 x4 0 2 -1 1 12 6 Zj 0 0 0 0 0 λj 0 -1 2 0

-1 x2 1/3 1 4/3 0 4 0 x4 -2/3 0 -11/3 1 4 Zj -1/3 -1 -4/3 0 -4 λj 1/3 0 13/3 0 所以x2=4,x1=x3=0,min=-4, max=4

(3)in m 1234f 235x x x x =--+-

12341234134123424623412.:4,,,0

x x x x x x x x s t x x x x x x x ++-≤??++-≤??++≤??≥?

w

解:LP 化为标准形式如下:

min f=-2x1-x2+3x3-5x4 s.t.x1+2x2+4x3-x4+x5=6 2x1+3x2+4x3-x4+x6=12 X1+x3+x4+x7=4

X1,x2,x3,x4,x5,x6,x7≥0

单纯形表如下:

-2 -1 3 -5 0 0 0

X1 x2 x3 x4 x5 x6 x7 b θ 0 x5 1 2 4 -1 1 0 0 6 0 x6 2 3 4 -1 0 1 0 12

0 x7 1 0 1 1 0 0 1 4 4 Zj 0 0 0 0 0 0 0 λj -2 -1 3 -5 0 0 0

0 x5 2 2 5 0 1 0 1 10 5 0 x6 3 3 5 0 0 1 1 16 5.3 -5 x4 1 0 1 1 0 0 1 4 Zj -5 0 -5 -5 0 0 -5 λj 3 -1 8 0 0 0 5

-1 x2 1 1 5/2 0 1/2 0 1/2 5 0 x6 0 0 -5/2 0 -3/2 1 -1/2 1 -5 x4 1 0 1 1 0 0 1 4 Zj -6 -1 -15/2 -5 -1/2 0 -11/2 -25 λj 4 0 21/2 0 1/2 0 1/2 所以x2=5,x4=4,x1=x3=0, min=-25 4.写出指派问题的数学模型

拟分配n 人去干n 项工作,每人干且仅干一项工作,若分配第i 人去干第j 项工作,需花费ij c 单位时间,问应如何分配工作才能使工人花费的总时间最少? 解:得到的数学模型为(假设第i 个人干第j 项工作是x ij ):

∑∑==n i n

j ij ij x c 11

min

11=∑=n

j ij

x

s.t. 11

=∑=n

i ij

x

5.运输问题的数学模型

某商品有m 个产地、n 个销地,各产地的产量分别为m a a ,,1 ,各销地的需求量分别为

n b b ,,1 。若该商品由i 产地运到j 销地的单位运价为ij c ,问应该如何调运才能使总运费最

省?

解:得到的数学模型为(假设从i 产地运到j 销地的数量是x ij ):

∑∑==m i n

j ij ij x c 11

min

i n

j ij

a x

=∑=1

s.t.

j m

i ij

b x

=∑=1

数值分析大作业-三、四、五、六、七

大作业 三 1. 给定初值 0x 及容许误差 ,编制牛顿法解方程f (x )=0的通用 程序. 解:Matlab 程序如下: 函数m 文件:fu.m function Fu=fu(x) Fu=x^3/3-x; end 函数m 文件:dfu.m function Fu=dfu(x) Fu=x^2-1; end 用Newton 法求根的通用程序Newton.m clear; x0=input('请输入初值x0:'); ep=input('请输入容许误差:'); flag=1; while flag==1 x1=x0-fu(x0)/dfu(x0); if abs(x1-x0)

while flag1==1 && m<=10^3 x1=x0-fu(x0)/dfu(x0); if abs(x1-x0)=ep flag=0; end end fprintf('最大的sigma 值为:%f\n',sigma); 2.求下列方程的非零根 5130.6651()ln 05130.665114000.0918 x x f x x +??=-= ?-???解:Matlab 程序为: (1)主程序 clear clc format long x0=765; N=100; errorlim=10^(-5); x=x0-f(x0)/subs(df(),x0); n=1; while nerrorlim n=n+1; else break ; end x0=x; end disp(['迭代次数: n=',num2str(n)]) disp(['所求非零根: 正根x1=',num2str(x),' 负根x2=',num2str(-x)]) (2)子函数 非线性函数f function y=f(x) y=log((513+0.6651*x)/(513-0.6651*x))-x/(1400*0.0918); end

数值分析大作业三 四 五 六 七

大作业 三 1. 给定初值 0x 及容许误差 ,编制牛顿法解方程f (x )=0的通用程序. 解:Matlab 程序如下: 函数m 文件:fu.m function Fu=fu(x) Fu=x^3/3-x; end 函数m 文件:dfu.m function Fu=dfu(x) Fu=x^2-1; end 用Newton 法求根的通用程序Newton.m clear; x0=input('请输入初值x0:'); ep=input('请输入容许误差:');

flag=1; while flag==1 x1=x0-fu(x0)/dfu(x0); if abs(x1-x0)

while flag==1 sigma=k*eps; x0=sigma; k=k+1; m=0; flag1=1; while flag1==1 && m<=10^3 x1=x0-fu(x0)/dfu(x0); if abs(x1-x0)=ep flag=0;

end end fprintf('最大的sigma 值为:%f\n',sigma); 2.求下列方程的非零根 5130.6651()ln 05130.665114000.0918 x x f x x +?? =-= ?-???解: Matlab 程序为: (1)主程序 clear clc format long x0=765; N=100; errorlim=10^(-5); x=x0-f(x0)/subs(df(),x0); n=1;

数值分析作业答案

数值分析作业答案 插值法 1、当x=1,-1,2时,f(x)=0,-3,4,求f(x)的二次插值多项式。 (1)用单项式基底。 (2)用Lagrange插值基底。 (3)用Newton基底。 证明三种方法得到的多项式是相同的。 解:(1)用单项式基底 设多项式为: , 所以: 所以f(x)的二次插值多项式为: (2)用Lagrange插值基底 Lagrange插值多项式为: 所以f(x)的二次插值多项式为: (3) 用Newton基底: 均差表如下: xk f(xk) 一阶均差二阶均差 1 0 -1 -3 3/2 2 4 7/ 3 5/6 Newton插值多项式为: 所以f(x)的二次插值多项式为: 由以上计算可知,三种方法得到的多项式是相同的。 6、在上给出的等距节点函数表,若用二次插值求ex的近似值,要使截断误差不超过10-6,问使用函数表的步长h应取多少? 解:以xi-1,xi,xi+1为插值节点多项式的截断误差,则有 式中 令得 插值点个数

是奇数,故实际可采用的函数值表步长 8、,求及。 解:由均差的性质可知,均差与导数有如下关系: 所以有: 15、证明两点三次Hermite插值余项是 并由此求出分段三次Hermite插值的误差限。 证明:利用[xk,xk+1]上两点三次Hermite插值条件 知有二重零点xk和k+1。设 确定函数k(x): 当或xk+1时k(x)取任何有限值均可; 当时,,构造关于变量t的函数 显然有 在[xk,x][x,xk+1]上对g(x)使用Rolle定理,存在及使得 在,,上对使用Rolle定理,存在,和使得 再依次对和使用Rolle定理,知至少存在使得 而,将代入,得到 推导过程表明依赖于及x 综合以上过程有: 确定误差限: 记为f(x)在[a,b]上基于等距节点的分段三次Hermite插值函数。在区间[xk,xk+1]上有 而最值 进而得误差估计: 16、求一个次数不高于4次的多项式,使它满足,,。

数值计算方法大作业

目录 第一章非线性方程求根 (3) 1.1迭代法 (3) 1.2牛顿法 (4) 1.3弦截法 (5) 1.4二分法 (6) 第二章插值 (7) 2.1线性插值 (7) 2.2二次插值 (8) 2.3拉格朗日插值 (9) 2.4分段线性插值 (10) 2.5分段二次插值 (11) 第三章数值积分 (13) 3.1复化矩形积分法 (13) 3.2复化梯形积分法 (14) 3.3辛普森积分法 (15) 3.4变步长梯形积分法 (16) 第四章线性方程组数值法 (17) 4.1约当消去法 (17) 4.2高斯消去法 (18) 4.3三角分解法 (20)

4.4雅可比迭代法 (21) 4.5高斯—赛德尔迭代法 (23) 第五章常积分方程数值法 (25) 5.1显示欧拉公式法 (25) 5.2欧拉公式预测校正法 (26) 5.3改进欧拉公式法 (27) 5.4四阶龙格—库塔法 (28)

数值计算方法 第一章非线性方程求根 1.1迭代法 程序代码: Private Sub Command1_Click() x0 = Val(InputBox("请输入初始值x0")) ep = Val(InputBox(请输入误差限ep)) f = 0 While f = 0 X1 = (Exp(2 * x0) - x0) / 5 If Abs(X1 - x0) < ep Then Print X1 f = 1 Else x0 = X1 End If Wend End Sub 例:求f(x)=e2x-6x=0在x=0.5附近的根(ep=10-10)

1.2牛顿法 程序代码: Private Sub Command1_Click() b = Val(InputBox("请输入被开方数x0")) ep = Val(InputBox(请输入误差限ep)) f = 0 While f = 0 X1 = x0 - (x0 ^ 2 - b) / (2 * b) If Abs(X1 - x0) < ep Then Print X1 f = 1 Else x0 = X1 End If Wend End Sub 例:求56的值。(ep=10-10)

北航数值分析大作业一

《数值分析B》大作业一 SY1103120 朱舜杰 一.算法设计方案: 1.矩阵A的存储与检索 将带状线性矩阵A[501][501]转存为一个矩阵MatrixC[5][501] . 由于C语言中数组角标都是从0开始的,所以在数组MatrixC[5][501]中检索A的带内元素a ij的方法是: A的带内元素a ij=C中的元素c i-j+2,j 2.求解λ1,λ501,λs ①首先分别使用幂法和反幂法迭代求出矩阵按摸最大和最小的特征值λmax和λmin。λmin即为λs; 如果λmax>0,则λ501=λmax;如果λmax<0,则λ1=λmax。 ②使用带原点平移的幂法(mifa()函数),令平移量p=λmax,求 出对应的按摸最大的特征值λ,max, 如果λmax>0,则λ1=λ,max+p;如果λmax<0,则λ501=λ,max+p。 3.求解A的与数μk=λ1+k(λ501-λ1)/40的最接近的特征值λik (k=1,2,…,39)。 使用带原点平移的反幂法,令平移量p=μk,即可求出与μk最接近的特征值λik。 4.求解A的(谱范数)条件数cond(A)2和行列式d etA。 ①cond(A)2=|λ1/λn|,其中λ1和λn分别是矩阵A的模最大和 最小特征值。

②矩阵A的行列式可先对矩阵A进行LU分解后,detA等于U所有对角线上元素的乘积。 二.源程序 #include #include #include #include #include #include #include #define E 1.0e-12 /*定义全局变量相对误差限*/ int max2(int a,int b) /*求两个整型数最大值的子程序*/ { if(a>b) return a; else return b; } int min2(int a,int b) /*求两个整型数最小值的子程序*/ { if(a>b) return b; else return a; } int max3(int a,int b,int c) /*求三整型数最大值的子程序*/ { int t; if(a>b) t=a; else t=b; if(t

数值分析第一次作业及参考答案

数值计算方法第一次作业及参考答案 1. 已测得函数()y f x =的三对数据:(0,1),(-1,5),(2,-1), (1)用Lagrange 插值求二次插值多项式。(2)构造差商表。(3)用Newton 插值求二次插值多项式。 解:(1)Lagrange 插值基函数为 0(1)(2)1 ()(1)(2)(01)(02)2 x x l x x x +-= =-+-+- 同理 1211 ()(2),()(1)36 l x x x l x x x = -=+ 故 2 20 2151 ()()(1)(2)(2)(1) 23631 i i i p x y l x x x x x x x x x =-==-+-+-++=-+∑ (2)令0120,1,2x x x ==-=,则一阶差商、二阶差商为 011215 5(1) [,]4, [,]20(1) 12 f x x f x x ---= =-= =----- 0124(2) [,,]102 f x x x ---= =- 实际演算中可列一张差商表: (3)用对角线上的数据写出插值多项式 2 2()1(4)(0)1*(0)(1)31P x x x x x x =+--+-+=-+ 2. 在44x -≤≤上给出()x f x e =的等距节点函数表,若用二次插值求x e 的近似值,要使 截断误差不超过6 10-,问使用函数表的步长h 应取多少 解: ()40000(), (),[4,4],,,, 1.x k x f x e f x e e x x h x x h x x th t ==≤∈--+=+≤考察点及

(3) 2000 4 43 4 3 () ()[(()]()[()] 3! (1)(1) (1)(1) 3!3! .(4,4). 6 f R x x x h x x x x h t t t e t h th t h e h e ξ ξ =----+ -+ ≤+??-= ≤∈- 则 4 36 ((1)(1) 100.006. t t t h - -+± << Q在点 得 3.求2 () f x x =在[a,b]上的分段线性插值函数() h I x,并估计误差。 解: 22 22 11 1 111 22 11 11 1 () () k k k k h k k k k k k k k k k k k k k k k k k x x x x x x I x x x x x x x x x x x x x x x x x x x x x ++ + +++ ++ ++ + --- =+= --- ?-? -=+- - [] 2 11 22 11 ()()()[()] 11 ()() 44 h h k k k k k k k k R x f x I x x x x x x x x x x x x x h ++ ++ =-=-+- =--≤-= 4.已知单调连续函数() y f x =的如下数据 用插值法计算x约为多少时() 1. f x=(小数点后至少保留4位) 解:作辅助函数()()1, g x f x =-则问题转化为x为多少时,()0. g x=此时可作新 的关于() i g x的函数表。由() f x单调连续知() g x也单调连续,因此可对() g x的数值进行反插。的牛顿型插值多项式为 1()0.110.097345( 2.23)0.451565( 2.23)( 1.10) 0.255894( 2.23)( 1.10)(0.17) x g y y y y y y y - ==-+++++ -++-

数值分析大作业

数值分析报大作业 班级:铁道2班 专业:道路与铁道工程 姓名:蔡敦锦 学号:13011260

一、序言 该数值分析大作业是通过C语言程序编程在Microsoft Visual C++ 6.0编程软件上运行实现的。本来是打算用Matlab软间来计算非线性方程的根的。学习Matlab也差不多有一个多月了,感觉自己编程做题应该没什么问题了;但是当自己真心的去编程、运行时才发现有很多错误,花了一天时间修改、调试程序都没能得到自己满意的结果。所以,我选择了自己比较熟悉的C程序语言来编程解决非线性的求值问题,由于本作业是为了比较几种方法求值问题的收敛速度和精度的差异,选择了一个相对常见的非线性函数来反映其差异,程序运行所得结果我个人比较满意。编写C语言,感觉比较上手,程序出现问题也能比较熟练的解决。最终就决定上交一份C程序语言编程的求值程序了!

二、选题 本作业的目的是为了加深对非线性方程求根方法的二分法、简单迭代法、、牛顿迭代法弦截法等的构造过程的理解;能将各种方法的算法描述正确并且能够改编为程序并在计算机上实现程序的正确合理的运行,能得到自己满意的结果,并且能调试修改程序中可能出现的问题和程序功能的增减修改。本次程序是为了比较各种方法在求解同一非线性方程根时,在收敛情况上的差异。 为了达到上面的条件我选择自己比较熟悉的语言—C语言来编程,所选题目为计算方程f(x)=x3-2x-5=0在区间[2,3]内其最后两近似值的差的绝对值小于等于5 ?的根的几种方法的比较。 110- 本文将二分法、牛顿法、简单迭代法、弦截法及加速收敛法这五种方法在同一个程序中以函数调用的方式来实现,比较简洁明了,所得结果能很好的比较,便于分析;发现问题和得出结论。

北航数值分析大作业第二题精解

目标:使用带双步位移的QR 分解法求矩阵10*10[]ij A a =的全部特征值,并对其中的每一个实特征值求相应的特征向量。已知:sin(0.50.2)() 1.5cos( 1.2)(){i j i j ij i j i j a +≠+== (i,j=1,2, (10) 算法: 以上是程序运作的逻辑,其中具体的函数的算法,大部分都是数值分析课本上的逻辑,在这里特别写出矩阵A 的实特征值对应的一个特征向量的求法: ()[]()() []()[]()111111I 00000 i n n n B A I gause i n Q A I u Bu u λλ-?-?-=-?-?? ?-=????→=??????→= ?? ? 选主元的消元 检查知无重特征值 由于=0i A I λ- ,因此在经过选主元的高斯消元以后,i A I λ- 即B 的最后一行必然为零,左上方变 为n-1阶单位矩阵[]()()11I n n -?-,右上方变为n-1阶向量[]()11n Q ?-,然后令n u 1=-,则 ()1,2,,1j j u Q j n ==???-。

这样即求出所有A所有实特征值对应的一个特征向量。 #include #include #include #define N 10 #define E 1.0e-12 #define MAX 10000 //以下是符号函数 double sgn(double a) { double z; if(a>E) z=1; else z=-1; return z; } //以下是矩阵的拟三角分解 void nishangsanjiaodiv(double A[N][N]) { int i,j,k; int m=0; double d,c,h,t; double u[N],p[N],q[N],w[N]; for(i=0;i

数值分析作业答案part

6.4.设??? ? ? ??=5010010a b b a A ,0det ≠A ,用a ,b 表示解线性方程组f Ax =的雅可比迭代与 高斯—塞德尔迭代收敛的充分必要条件。 解 雅可比迭代法的迭代矩阵 ? ??? ??? ? ??----=???? ? ??----????? ??=-050100100100000001010101 a b b a a b b a B J , ?? ? ?? -=-1003||2ab B I J λλλ,10||3)(ab B J = ρ。 雅可比迭代法收敛的充分必要条件是3 100 ||

上海大学_王培康_数值分析大作业

数值分析大作业(2013年5月) 金洋洋(12721512),机自系 1.下列各数都是经过四舍五入得到的近似值,试分别指出它 们的绝对误差限, 相对误差限和有效数字的位数。 X1 =5.420, x 2 =0.5420, x 3=0.00542, x 4 =6000, x 5=50.610? 解:根据定义:如果*x 的绝对误差限 不超过x 的某个数位的半个单位,则从*x 的首位非零数字到该位都是有效数字。 显然根据四舍五入原则得到的近视值,全部都是有效数字。 因而在这里有:n1=4, n2=4, n3=3, n4=4, n5=1 (n 表示x 有效数字的位数) 对x1:有a1=5, m1=1 (其中a1表示x 的首位非零数字,m1表示x1的整数位数) 所以有绝对误差限 143 11 (1)101022 x ε--≤ ?=? 相对误差限 31() 0.510(1)0.00923%5.4201 r x x x εε-?= == 对x2:有a2=5, m2=0 所以有绝对误差限 044 11 (2)101022 x ε--≤ ?=? 相对误差限 42() 0.510(2)0.00923%0.54202 r x x x εε-?= == 对x3:有a3=5, m3=-2 所以有绝对误差限 235 11 (3)101022 x ε---≤ ?=? 相对误差限 53() 0.510(3)0.0923%0.005423 r x x x εε-?= == 对x4:有a4=0, m4=4 所以有绝对误差限 4411(4)1022 x ε-≤?= 相对误差限 4() 0.5 (4)0.0083%6000 4 r x x x εε= = = 对x5:有a5=6, m5=5 所以有绝对误差限 514 11(5)101022 x ε-≤ ?=? 相对误差限 45() 0.510(5)8.3%600005 r x x x εε?= ==

北航数值分析报告第三次大作业

数值分析第三次大作业 一、算法的设计方案: (一)、总体方案设计: x y当作已知量代入题目给定的非线性方程组,求(1)解非线性方程组。将给定的(,) i i

得与(,)i i x y 相对应的数组t[i][j],u[i][j]。 (2)分片二次代数插值。通过分片二次代数插值运算,得到与数组t[11][21],u[11][21]]对应的数组z[11][21],得到二元函数z=(,)i i f x y 。 (3)曲面拟合。利用x[i],y[j],z[11][21]建立二维函数表,再根据精度的要求选择适当k 值,并得到曲面拟合的系数矩阵C[r][s]。 (4)观察和(,)i i p x y 的逼近效果。观察逼近效果只需要重复上面(1)和(2)的过程,得到与新的插值节点(,)i i x y 对应的(,)i i f x y ,再与对应的(,)i i p x y 比较即可,这里求解 (,)i i p x y 可以直接使用(3)中的C[r][s]和k 。 (二)具体算法设计: (1)解非线性方程组 牛顿法解方程组()0F x =的解* x ,可采用如下算法: 1)在* x 附近选取(0) x D ∈,给定精度水平0ε>和最大迭代次数M 。 2)对于0,1, k M =执行 ① 计算() ()k F x 和()()k F x '。 ② 求解关于() k x ?的线性方程组 () ()()()()k k k F x x F x '?=- ③ 若() () k k x x ε∞∞ ?≤,则取*()k x x ≈,并停止计算;否则转④。 ④ 计算(1) ()()k k k x x x +=+?。 ⑤ 若k M <,则继续,否则,输出M 次迭代不成功的信息,并停止计算。 (2)分片双二次插值 给定已知数表以及需要插值的节点,进行分片二次插值的算法: 设已知数表中的点为: 00(0,1,,) (0,1,,)i j x x ih i n y y j j m τ=+=???=+=?? ,需要插值的节点为(,)x y 。 1) 根据(,)x y 选择插值节点(,)i j x y : 若12h x x ≤+ 或12 n h x x ->-,插值节点对应取1i =或1i n =-,

数值分析第二次大作业

《数值分析》计算实习报告 第二题 院系:机械工程及自动化学院 学号: 姓名: 2017年11月

一、题目要求 试求矩阵A =[a ij ]10×10的全部特征值,并对其中的每一个实特征值求相应的特征向量,已知 a ij ={ sin (0.5i +0.2j ) i ≠j 1.52cos (i +1.2j ) i =j (i,j =1,2, (10) 说明: 1.用带双步位移的QR 方法求矩阵特征值,要求迭代的精度水平为ε=10?12。 2.打印以下内容: (1)全部源程序; (2)矩阵A 经过拟上三角化后所得的矩阵A (n?1); (3)对矩阵A (n?1)实行QR 方法迭代结束后所得的矩阵; (4)矩阵A 的全部特征值λi =(R i ,I i ) (i =1,2,?,10),其中R i =Re(λi ),I i = Im(λi ) 。若λi 是实数,则令I i =0; (5)A 的相应于实特征值的特征向量。 3.采用e 型数输出实型数,并且至少显示12位有效数字。 二、算法设计思路和方案 1. 将矩阵A 拟上三角化得到矩阵A (n?1) 为了减少计算量,一般先利用Householder 矩阵对矩阵A 作相似变换,把A 化为拟上三角矩阵A (n?1),然后用QR 方法计算A (n?1)的全部特征值,而A (n?1)的特征值就是A 的特征值。具体算法如下: 记(1)A A =,()r A 的第r 列至第n 列的元素为(r)(1,2, ,;,1,,)ij a i n j r r n ==+。 对于1,2,,2r n =-执行 (1)若() (2,3,,)r ir a i r r n =++全为零,则令(1)()r r A A +=,转(5);否则转(2)。

北航数值分析大作业第一题幂法与反幂法

《数值分析》计算实习题目 第一题: 1. 算法设计方案 (1)1λ,501λ和s λ的值。 1)首先通过幂法求出按模最大的特征值λt1,然后根据λt1进行原点平移求出另一特征值λt2,比较两值大小,数值小的为所求最小特征值λ1,数值大的为是所求最大特征值λ501。 2)使用反幂法求λs ,其中需要解线性方程组。因为A 为带状线性方程组,此处采用LU 分解法解带状方程组。 (2)与140k λλμλ-5011=+k 最接近的特征值λik 。 通过带有原点平移的反幂法求出与数k μ最接近的特征值 λik 。 (3)2cond(A)和det A 。 1)1=n λλ2cond(A),其中1λ和n λ分别是按模最大和最小特征值。 2)利用步骤(1)中分解矩阵A 得出的LU 矩阵,L 为单位下三角阵,U 为上三角阵,其中U 矩阵的主对角线元素之积即为det A 。 由于A 的元素零元素较多,为节省储存量,将A 的元素存为6×501的数组中,程序中采用get_an_element()函数来从小数组中取出A 中的元素。 2.全部源程序 #include #include void init_a();//初始化A double get_an_element(int,int);//取A 中的元素函数 double powermethod(double);//原点平移的幂法 double inversepowermethod(double);//原点平移的反幂法 int presolve(double);//三角LU 分解 int solve(double [],double []);//解方程组 int max(int,int); int min(int,int); double (*u)[502]=new double[502][502];//上三角U 数组 double (*l)[502]=new double[502][502];//单位下三角L 数组 double a[6][502];//矩阵A int main() { int i,k; double lambdat1,lambdat2,lambda1,lambda501,lambdas,mu[40],det;

北航数值分析大作业第二题

数值分析第二次大作业 史立峰 SY1505327

一、 方案 (1)利用循环结构将sin(0.50.2)() 1.5cos( 1.2)() {i j i j ij i j i j a +≠+==(i,j=1,2,……,10)进行赋值,得到需要变换的 矩阵A ; (2)然后,对矩阵A 利用Householder 矩阵进行相似变换,把A 化为上三角矩阵A (n-1)。 对A 拟上三角化,得到拟上三角矩阵A (n-1),具体算法如下: 记A(1)=A ,并记A(r)的第r 列至第n 列的元素为()n r r j n i a r ij ,,1,;,,2,1) ( +==。 对于2,,2,1-=n r 执行 1. 若 ()n r r i a r ir ,,3,2) ( ++=全为零,则令A(r+1) =A(r),转5;否则转2。 2. 计算 () ∑+== n r i r ir r a d 1 2 )( ()( )r r r r r r r r r r d c a d a c ==-=++则取,0sgn ) (,1)(,1若 )(,12r r r r r r a c c h +-= 3. 令 () n T r nr r r r r r r r r R a a c a u ∈-=++) ()(,2)(,1,,,,0,,0 。 4. 计算 r r T r r h u A p /)(= r r r r h u A q /)(= r r T r r h u p t /= r r r r u t q -=ω T r r T r r r r p u u A A --=+ω)()1( 5. 继续。 (3)使用带双步位移的QR 方法计算矩阵A (n-1)的全部特征值,也是A 的全部特征值,具体算法如下: 1. 给定精度水平0>ε和迭代最大次数L 。 2. 记n n ij n a A A ?-==][) 1()1()1(,令n m k ==,1。

数值分析作业答案(第5章)

5.1.设A 是对称矩阵且011≠a ,经过一步高斯消去法后,A 约化为 ?? ????21 110 A a a T 证明2A 是对称矩阵。 证明 由消元公式及A 的对称性,有 ,,,3,2,,)2(111 11111 )2(n j i a a a a a a a a a a ji i j ji j i ij ij ==-=- = 故2A 对称。 5.2.设n ij a A )(=是对称正定矩阵,经过高斯消去法一步后,A 约化为 ?? ????21 110 A a a T 其中1)2(2)(-=n ij a A 。证明: (1).A 的对角元素;,,2,1,0n i a ii => (2).2A 是对称正定矩阵。 证明 (1).因为A 对称正定,所以 n i e Ae a i i ii ,,2,1,0),( =>=, 其中T i e )0,,0,1,0,,0( =为第i 个单位向量。 (2).由A 的对称性及消元公式,有 ,,,3,2,,)2(111 11111 )2(n j i a a a a a a a a a a ji i j ji j i ij ij ==-=- = 故2A 也对称。 又由A L A a a T 121110=????? ?,其中

??? ?????- =? ????? ? ?????????--=-111 1 11111 21101 1011n n I a a a a a a L , 可见1L 非奇异,因而对任意0≠x ,由A 的正定性,有 ,0),(),(,011111>=≠x AL x L x AL L x x L T T T T 故T AL L 11正定。 由,000110211 111121111 1?? ? ?? ?=????????-??????=-A a I a a A a a AL L n T T T 而011>a ,故知2A 正定

数值分析大作业 三、四、五、六、七

大作业 三 1. 给定初值0x 及容许误差 ,编制牛顿法解方程f (x )=0的通用 程序. 解:Matlab 程序如下: 函数m 文件:fu.m function Fu=fu(x) Fu=x^3/3-x; end 函数m 文件:dfu.m function Fu=dfu(x) Fu=x^2-1; end 用Newton 法求根的通用程序Newton.m clear; x0=input('请输入初值x0:'); ep=input('请输入容许误差:'); flag=1; while flag==1 x1=x0-fu(x0)/dfu(x0); if abs(x1-x0)

while flag1==1 && m<=10^3 x1=x0-fu(x0)/dfu(x0); if abs(x1-x0)=ep flag=0; end end fprintf('最大的sigma 值为:%f\n',sigma); 2.求下列方程的非零根 5130.6651()ln 05130.665114000.0918 x x f x x +?? =- = ?-???解:Matlab 程序为: (1)主程序 clear clc format long x0=765; N=100; errorlim=10^(-5); x=x0-f(x0)/subs(df(),x0); n=1; while nerrorlim n=n+1; else break ; end x0=x; end disp(['迭代次数: n=',num2str(n)]) disp(['所求非零根: 正根x1=',num2str(x),' 负根x2=',num2str(-x)]) (2)子函数 非线性函数f function y=f(x) y=log((513+0.6651*x)/(513-0.6651*x))-x/(1400*0.0918); end

哈尔滨工程大学数值分析大作业2014-附fortran程序

B班大作业要求: 1. 使用统一封皮; 2. 上交大作业内容包含: 一摘要 二数学原理 三程序设计(必须对输入变量、输出变量进行说明;编程无语言要求,但程序要求通过)四结果分析和讨论 五完成题目的体会与收获 3. 提交大作业的时间:本学期最后一次课,或考前答疑;过期不计入成绩; 4. 提交方式:打印版一份;或手写大作业,但必须使用A4纸。 5. 撰写的程序需打印出来作为附录。

课程设计 课程名称: 设计题目: 学号: 姓名: 完成时间:

题目一:非线性方程求根 一 摘要 非线性方程的解析解通常很难给出,因此非线性方程的数值解就尤为重要。本实验通过使用常用的求解方法二分法和Newton 法及改进的Newton 法处理几个题目,分析并总结不同方法处理问题的优缺点。观察迭代次数,收敛速度及初值选取对迭代的影响。 用Newton 法计算下列方程 (1) 310x x --= , 初值分别为01x =,00.45x =,00.65x =; (2) 32943892940x x x +-+= 其三个根分别为1,3,98-。当选择初值02x =时给出结果并分析现 象,当6 510ε-=?,迭代停止。 二 数学原理 对于方程f(x)=0,如果f(x)是线性函数,则它的求根是很容易的。牛顿迭代法实质上是一种线性化方法,其基本思想是将非线性方程f(x)=0逐步归结为某种线性方程来求解。 设已知方程f(x)=0有近似根x k (假定k f'(x )0≠) ,将函数f(x)在点x k 进行泰勒展开,有 k k k f(x)f(x )+f'(x )(x-x )+≈??? 于是方程f(x)=0可近似的表示为 k k k f(x )+f'(x )(x-x )=0 这是个线性方程,记其根为x k+1,则x k+1的计算公式为 k+1k () x =x -'() k k f x f x ,k=0,1,2,… 这就是牛顿迭代法或简称牛顿法。

数值分析作业答案.doc

第2章 插值法 1、当x=1,-1,2时,f(x)=0,-3,4,求f(x)的二次插值多项式。 (1)用单项式基底。 (2)用Lagrange 插值基底。 (3)用Newton 基底。 证明三种方法得到的多项式是相同的。 解:(1)用单项式基底 设多项式为:2 210)(x a x a a x P ++=, 所以:64 211111 1111122 2 211 200 -=-==x x x x x x A 3 76144 211111114241 13110111)() ()(22 221120 022 2 22 11 120 00-=-= ---==x x x x x x x x x f x x x f x x x f a 2 3694211111114411 31101111)(1)(1 )(122 221120 02 2 22112 001=--= --==x x x x x x x x f x x f x x f a 6 5654 2 1 1111114 2 1 3 11011111) (1)(1)(122 2 21120 022 11 00 2=--= ---==x x x x x x x f x x f x x f x a 所以f(x)的二次插值多项式为:26 52337)(x x x P ++-= (2)用Lagrange 插值基底 )21)(11() 2)(1())(())(()(2010210-+-+=----=x x x x x x x x x x x l )21)(11() 2)(1())(())(()(2101201------=----=x x x x x x x x x x x l ) 12)(12() 1)(1())(())(()(1202102+-+-=----= x x x x x x x x x x x l

北航数值分析报告大作业第八题

北京航空航天大学 数值分析大作业八 学院名称自动化 专业方向控制工程 学号 学生姓名许阳 教师孙玉泉 日期2014 年11月26 日

一.题目 关于x , y , t , u , v , w 的方程组(A.3) ???? ?? ?=-+++=-+++=-+++=-+++79 .0sin 5.074.3cos 5.007.1cos sin 5.067.2cos 5.0y w v u t x w v u t y w v u t x w v u t (A.3) 以及关于z , t , u 的二维数表(见表A-1)确定了一个二元函数z =f (x , y )。 表A-1 二维数表 t z u 0 0.4 0.8 1.2 1.6 2 0 -0.5 -0.34 0.14 0.94 2.06 3.5 0.2 -0.42 -0.5 -0.26 0.3 1.18 2.38 0.4 -0.18 -0.5 -0.5 -0.18 0.46 1.42 0.6 0.22 -0.34 -0.58 -0.5 -0.1 0.62 0.8 0.78 -0.02 -0.5 -0.66 -0.5 -0.02 1.0 1.5 0.46 -0.26 -0.66 -0.74 -0.5 1. 试用数值方法求出f (x , y ) 在区域}5.15.0,8.00|), {≤≤≤≤=y x y x D (上的近似表达式 ∑∑===k i k j s r rs y x c y x p 00 ),( 要求p (x , y )以最小的k 值达到以下的精度 ∑∑==-≤-=10020 7210)],(),([i j i i i i y x p y x f σ 其中j y i x i i 05.05.0,08.0+==。 2. 计算),(),,(* ***j i j i y x p y x f (i =1,2,…,8 ; j =1,2,…,5) 的值,以观察p (x , y ) 逼 近f (x , y )的效果,其中j y i x j i 2.05.0,1.0**+==。

相关主题
文本预览
相关文档 最新文档