实验3 - 非线性方程AX=0的解法

《数值计算方法》实验报告 start 21 [N,N]=size(A);X=zeros(N,1);C=zeros(1,N+1);Aug=[A B]; p=1 N=length(B); k=1 k=k+1 k>max1 Y N j=1 j=j+1 j>N Y N X(1)=(B(1)-A(1,2)*P(2))/A(1,1) Y j==1 N X(2)=(B(2)-A(2,1)*X(1)'-A(2,3:4)*P(3:4))/A(2,2) Y j==2 N Y j==N-1 X(N-1)=(B(N-1)-A(N-1,N-3:N-2)*X(N-3:N-2)' Y X(N)=(B(N)-A(N,N-1)*(X(N-1))')/A(N,N) N j==N X(j)=(B(j)-A(j,j-1)*X(j-1)'-A(j,j+1)*P(j+1))/A(j,j) err=abs(norm(X'-P)); relerr=err/(norm(X)+eps);P=X'; N err

表4.通过高斯赛德尔迭代得到的x的解

x 1~10 0.46379552381655 0.5000001084251992 0.5000000000017429 0.499999999992555 0.4999978569134675 0.5372846051999656 0.5000002015766873 0.5000000000009169 0.499999999966017 0.4999947932669753 0.5090229246013306 0.500000022610945 0.5 0.4982216344361741 0.4999999862385722 0.4999999999999205 0.5000000004390391 0.5000887238901357 0.4989418602397619 0.499999995873979 0.5 0.4999853512481308 0.5000000005302938 0.5 0.5000887238901357 0.5000000004390388 0.4999999999999212 0.4999999862385723 0.4982216344361741 0.500015318846052 0.5000000000239392 0.5 0.4999947932669753 0.499999999966016 0.5000000000009176 0.5000002015766873 0.5372846051999657 0.4999978569134675 0.4999999999925573 0.5000000000017445 0.500000108425199 0.4637955238165501 11~20 21~30 31~40 0.500000000023939 0.500015318846052 0.5000000005302936 0.4999853512481308 0.499999995873979 0.4989418602397619 0.5000000226109451 0.5090229246013306 41~50

4.结论

P93.1,2,3

1）利用矩阵乘的方法可对已知三维图形进行坐标变换，实现该图形在空间的旋转。 2)矩阵乘一般不满足交换律，由于利用矩阵乘的方法对三维图形进行坐标变换，例如V=R1*U;W=R2*V，此时W?R2*R1*U;故在进行多次旋转时，需要分步进行，不可一步完成旋转。 P108.1,7

通过对系数矩阵的增广矩阵进行高斯消元和回带容易得到线性方程组的解，同时，利用这种方法可以求得矩阵的逆。 P120.3,4

如果系数矩阵能分解为LU的形式，其中L为下三角矩阵，U为上三角矩阵，通过对系数矩阵的分解再求解可应用简单的迭代进行求解x。

《数值计算方法》实验报告

P129.4.

求解带状线性方程组的解可使用高斯-赛德尔迭代法。

附件（代码）： P93.1:

X=zeros(8,3);

X([5:8,11,12,15,16,18,20,22,24])=1;

d=[1 2 4 3 1 5 6 8 7 5 6 2 4 8 7 3];

for i=0:1:100

r1=[cos(i*pi/600) -sin(i*pi/600) 0;0 1 0;-sin(i*pi/600) 0 cos(i*pi/600)];

U=X*r1';

plot3(U(d,1),U(d,2),U(d,3)); drawnow end

for i=0:1:100

r2=[cos(i*pi/400) -sin(i*pi/400) 0;sin(i*pi/400) cos(i*pi/400) 1];

W=U*r2';

plot3(W(d,1),W(d,2),W(d,3)); drawnow end

subplot(2,2,1);

plot3(X(d,1),X(d,2),X(d,3)); subplot(2,2,2);

plot3(U(d,1),U(d,2),U(d,3)); subplot(2,2,3);

plot3(W(d,1),W(d,2),W(d,3)); xlabel('x') ylabel('y') zlabel('z')

view(3); rotate3d;

23

0 0;0 《数值计算方法》实验报告 24

p93.2

X=zeros(8,3);

X([5:8,11,12,15,16,18,20,22,24])=1;

d=[1 2 4 3 1 5 6 8 7 5 6 2 4 8 7 3];

for i=0:1:100

r1=[1 0 0;0 cos(i*pi/1200) -sin(i*pi/1200) ;0 sin(i*pi/1200) cos(i*pi/1200)];

U=X*r1';

plot3(U(d,1),U(d,2),U(d,3)); drawnow end

for i=0:1:100

r2=[cos(i*pi/600) -sin(i*pi/600) 0;sin(i*pi/600) cos(i*pi/600) 0;0 0 1];

W=U*r2';

plot3(W(d,1),W(d,2),W(d,3)); drawnow

xlabel('x') ylabel('y') zlabel('z') end

subplot(2,2,1);

plot3(X(d,1),X(d,2),X(d,3)); subplot(2,2,2);

plot3(U(d,1),U(d,2),U(d,3)); subplot(2,2,3);

plot3(W(d,1),W(d,2),W(d,3)); view(3); rotate3d; xlabel('x') ylabel('y') zlabel('z') p93.3

X=zeros(4,3);