数值稳定性验证实验报告

实验课程：数值计算方法专业：数学与应用数学班级：08070141

学号：37

*名：***

中北大学理学院

实验1 赛德尔迭代法

【实验目的】

熟悉用塞德尔迭代法解线性方程组 【实验内容】

1.了解MATLAB 语言的用法

2.用塞德尔迭代法解下列线性方程组

1234123412341234

541012

581034

x x x x x x x x x x x x x x x x ---=-⎧⎪-+--=⎪⎨

--+-=⎪⎪---+=⎩ 【实验所使用的仪器设备与软件平台】

计算机，MATLAB7.0

【实验方法与步骤】

1.先找出系数矩阵A ，将前面没有算过的x j 分别和矩阵的(,)A i j 相乘,然后将累加的和赋值给sum ，即(),j sum sum A i j x =+⋅.算出

()/(,)i i x b sum A i i =-,依次循环,算出所有的i x 。

2.若i x 前后两次之差的绝对值小于所给的误差限ε，则输出i x .否则重复以上过程，直到满足误差条件为止.

【实验结果】

(A 是系数矩阵，b 是右边向量，x 是迭代初值，ep 是误差限)

function y=seidel(A,b,x,ep) n=length(b); er=1; k=0; while er>=ep

k=k+1;

for i=[1:1:n]

q=x(i);

sum=0;

for j=[1:1:n]

if j~=i

sum=sum+A(i,j)*x(j);

end

end

x(i)=(b(i)-sum)/A(i,i);

er=abs(q-x(i));

end

end

fprintf('迭代次数k=%d\n',k)

disp(x')

【结果分析与讨论】

>> A=[5 -1 -1 -1;-1 10 -1 -1;-1 -1 5 -1;-1 -1 -1 10]; b=[-4 12 8 34];

seidel(A,b,[0 0 0 0],1e-3)

迭代次数k=6

0.99897849430002

1.99958456867649

2.99953139743435

3.99980944604109

实验课程：数值计算方法专业：数学与应用数学班级：08070141

学号：37

姓名：汪鹏飞

中北大学理学院

实验2 最小二乘法的拟合

【实验目的】

熟悉最小二乘法的拟合方法 【实验内容】

1.了解MATLAB 语言的用法

2.给出数据

希望用一次，二次多项式利用最小二乘法拟合这些数据，试写出正规方程组，并求出最小平方逼近多项式。

【实验所使用的仪器设备与软件平台】

计算机，MATLAB7.0，

【实验方法与步骤】

1.先算出 1

n j j i i l x ==∑，再算出 (1)1

n

j j i i i b x y -==∑

2.在正规方程组

20121111

231

01211111

120121

1111n n n n

m

i i m i i

i i i i n n n n n

m i i i m i i i

i i i i i n n n n n

m m m m m m i i i m i i i

i i i i i a n a x a x a x y a x a x a x a x x y a x a x a x a x x y ====+=====+++=====⎧++++=⎪

⎪

⎪++++=⎪⎨⎪⎪

⎪++++=⎪⎩∑∑∑∑∑∑∑∑∑∑∑∑∑∑ 中令1

n j j i i l x ==∑，(1)1

n

j j i i i b x y -==∑

3.令正规方程组的系数矩阵为A ，当j m ≤时，将j l 的值赋给(1,1)A j +， 当j>m 时，将j l 的值赋给(1,1)A j m m +-+，再将矩阵A 的其他元素写出来，于是正规矩阵可写成Aa b =,最后用1a A b -=即可算出向量a,向量a 的

元素依次是常数项,一次项的系数,二次项的系数m次项的系数.

4.由系数即可写出拟合的多项式曲线.

【实验结果】

（x，y为给定的对应值，m是要求的拟合次数)

function a=zxecf(x,y,m)

n=length(x);

A(1,1)=n;

for j=[1:1:2*m]

l(j)=0;

for i=[1:1:n]

l(j)=l(j)+x(i)^j;

end

if j<=m

A(1,1+j)=l(j);

else A(j+1-m,m+1)=l(j);

end

end

for j=[1:1:m+1]

b(j)=0;

for i=[1:1:n]

b(j)=b(j)+y(i)*x(i)^(j-1);

end

end

for i=[2:1:m+1]

k=-1;

for j=[m+1:-1:1]

k=k+1;

A(i,j)=l(m+i-1-k);