《矩阵与数值分析实验.docx》由会员分享,可在线阅读,更多相关《矩阵与数值分析实验.docx(12页珍藏版)》请在第一文库网上搜索。
1、数值实验题目一、线性方程组求解对于线性方程组-4x + x + 刍 + 4 = 1xl - 4x2 + x3 + x4 = 1x1 + x2 - 4x3 +x4 = 1%1 + x2 + x3 - 4x4 = 1(1)用直接法求解;(2)用Jaco瓦迭代法求解;(3)分别取力= 0.75,1.0,1.25,1.5,用SOR方法求解.比较迭代结果(与精确解比较).用gauss消去程序如下(mallab):function x=Gauss, (, b)A=-4, 1, 1, 1; 1,-4, 1, 1;1, 1,-4, 1;1, 1, 1,-4;n=length (b);for k=l: n-1
2、%消去nT次index = k+l:n;% 用向量m = -A (index, k)A(k, k);A(index, index) = A(index, index) + m*A(k, index) ;%矩阵A消元b (index) = b(index) + m*b(k);%矩阵b消元endx = zeros (n, 1);x(n) = b(n)A(n, n) ;%回代for i - nl:-1:1x(i) = ( b(i) - A(i, i+l:n)*x(i+l:n) )A(i,i)lendx结果:X 二-1.00000000000000-1.00000000000000-1.0000000
3、0000000-1.00000000000000(2) Jaco抗迭代法A=-4, 1, 1,1 ; 1, -4, 1, 1; 1, 1, -4,1; 1, 1,1, 4;x0=0 0 0 0, ; %初始n=length(B);M=151%迭代次数for k=l:Mfor j=l:nx(j) = (B(j)-A(j,j+l:n)*xO(l:j-l, j+l:n)/A(j, j);endx=,;endJacobi_X=x结果:Jacobi_X 二-0.98663653898984-0.98663653898984-0. 98663653898984-0. 98663653898984SOR:f
4、unction tx= sor ( , b, imax, x, tol, w)A=-4, 1, 1, 1; 1,-4, 1, 1;1, 1,-4, 1;1, 1, 1,-4;imax=33;x0=0 0 0 0;tol=0.001;w=l. 1;n=length (x);del=le-10 ; %主封角的元素不能太小必须大於deltx=x;for i=l:ndg=A(i, i);if abs(dg) delreturnendend%for k = 1 :imaxx=x ;为保留前一个近似值for i = l:nsm=b(i);for j = l:nif J =1sm = sm -A(i, j)
5、*x(j);endend %for jx(i)=smA(i, i) ; %下一个近似解x(i)=w*x(i) + (l-w)*x(i);end %for itx=tx ;x ; %吧下一个近似值放到表里面if norm(-O)二、三次样条插值已知函数值A012345678910yk00.791.532.192.713.033.272.893.063.193.29和边界条件:s(0) = 0.8, (10) = 0.2. 求三次样条插值函数y = s(x)并画出其图形.解:程序:(chazhi.m)format long g;xk=f() 123456789 1();% 给定节点y=0 0.79
6、 1.53 2.19 2.71 3.03 3.27 2.89 3.06 3.19 3.29;% 节点处的值h=zeros(l,10);%分别设置初始状态为0向量c=zeros(l,9); d=zeros(l,9); g=zeros(l,9); b=zeros(l,9); % 设置初始状态均为 0 向量m2=zeros(9,l); ml =zeros(9,l);m0=0.8; ml0=0.2;%已知的两个边界条件for i=l:l:10h(i)=xk(i+l)-xk(i); % 算出 hendfor i=l:l:9c(i)=h(i+l)/(h(i+l)+h(i);d =h (h(i+l)+h )
7、;endfor i=l:l:9g(i)=3*(d(i)*(y(i+2)-y(i+ l)h(i)+c(i)*(y(i+1 )y(i)h );endv=2222 22222;% 对角元均为 2A=diag(v);%生成矩阵for i=l:l:8A(i,i+l)=d ;A(i+l,i)=c(i+l);b(i)=g(i);endb(l)=g(l)-c(l)*mO;b(9)=g(9)-d(9)*ml0;L,U=lu(A);%将A进行LU分解,解方程,得出m向量ml=Lb;% LU 法求解m2=Um 1;m(l)=m0; m( 11 )=m 10;for i=2:l:10m(i)=m2(i-l);endf
8、or i=l:l:10;% 计算 S(x)的格式x=i-l:0.1:i; sx(i,:)=(l-2.*(x-xk(i)./(-h(i).*(x-xk(i+l)./(-h(i).A(2).*y(i)+.(L2.*(x-xk(i+ l)./h(i).*(x-xk(i)./h(i).A(2).*y(i+1)+ .(x-xk(i).*(x-xk(i+1 )-h(i) (2).*m(i)+ .(x-xk(i+1 ).*(x-xk(i).h(i).(2).*m(i+1);endsxl=zeros(l,101);xl=0:0.1:10;forj=l:10;for i=l:ll;sxl(l,il O*(j-l
9、)=sx(j,i);endendplot(x1,sxl);% 画出图形sx;m=m%输出m向量运行结果: chazhim =0.80.7714859631499640.7040561474001430.6122894472494650.3867860636019980.360566298342541-0.149051256972164-0.1843612704538850.2564963387877020.05837591530307440.2图形如下:三、求的近似值,要求绝对误差小于i(r4.解:程序:(fuhua.m)%输入为积分上下限,输出为积分值function y=fuhua(a,b
10、)n=l;tl =(b-a)/2*(exp(-aA2)/a+exp(-bA2)/b);while 1h=0;for i=0:n- lh=exp(-(a+(2*i+1 )*(b-a)2n)(2)(a+(2*i+1 )*(b-a)/2/n)+h;% 复化中矩形公式endh2=(b-a)/n*h;y=l/(b-a)*(tl+h2);ifabs(y-tl)0.5e-4;% 判断绝对误差break;endtl=y;n=2*n;end运行结果: y=fuhua(l, 3)y =0.10969135555231x2四、构造迭代法求方程sinx- = 0的所有根,要求误差计算不超过105).2解:Sinx的值域为L +1, x八2/2的值域为0:无穷,所以方程的根应在0,2八0.5之间,且2根。将有根区间0,pi二分搜索单根区间:经搜索两单根区间为0,pi4J, pi4,pi2J(1)在(),pi4上时:-5.916456789157589e-030J Figures - Figure 1ommand Window-0.51Shortcuts 2 How to Add 回 Wat,s New.,.一_0.150.05023456u = 2f-3u,“(0) = 1,0rl MATLABFile Edit View Insert Tools Debug Desktop VMndow Help