数值分析代码.doc

《数值分析代码.doc》由会员上传分享，免费在线阅读，更多相关内容在应用文档-天天文库。

1、第二章求解线性方程组的数值方法2.1Gauss消去法functionx=gauss(a,b)n=length(b);fork=1:nfori=k+1:nforj=k+1:na(i,j)=a(i,j)-a(i,k)*a(k,j)/a(k,k);endb(i)=b(i)-a(i,k)*b(k)/a(k,k);endfori=n:-1:1%以下为回代过程s=b(i);ifi==n;x(n)=b(n)/a(n,n);elseendforj=i+1:n;s=s-a(i,j)*x(j);endx(i)=s/a(i,i);endendend例子：输入矩阵a=[2,-1,3;4

2、,2,5;1,2,0];b=[1;4;7];调用函数x=gauss(a,b)答案为9-1-6列主元Gauss消去法functionx=gausslie(a,b)n=length(b);fori=1:nmax=abs(a(i,i));m=i;forj=i+1:nifmax

3、素法）fork=1:nfori=k+1:nforj=k+1:na(i,j)=a(i,j)-a(i,k)*a(k,j)/a(k,k);endb(i)=b(i)-a(i,k)*b(k)/a(k,k);endfori=n:-1:1s=b(i);ifi==n;x(n)=b(n)/a(n,n);elseendforj=i+1:n;s=s-a(i,j)*x(j);endx(i)=s/a(i,i);endendend%后面这段就是Gauss消去法的代码例子可以用上面的，结果是一样的2.2迭代法需注意迭代格式是怎样推出来的（逐次逼近法），注意L-U分解是怎样实现的Jacobi迭

4、代法functionx=jacobi(A,b)n=length(A);D=zeros(n,n);L=zeros(n,n);U=zeros(n,n);fori=1:n;D(i,i)=A(i,i);endfori=2:nforj=1:i-1L(i,j)=-A(i,j);U(j,i)=-A(j,i);endendB=inv(D)*(L+U);f=inv(D)*b;x0=ones(n,1);x1=B*x0+f;epsin=1e-5;whilemax(abs(x1-x0))>epsinx0=x1;x1=B*x0+f;endx=x1;例子A=[8,-3,2;4,11,-1;

5、6,3,12]b=[20;33;36]x=jacobi(A,b)答案为321Gauss-Seidel迭代法functionx=gaussseidel(A,b,p,eps)%p为x的零矩阵,eps为精度L=-tril(A,-1);D=diag(diag(A));U=-triu(A,1);G=(inv(D-L))*U;f=(inv(D-L))*b;x=G*p+f;n=1;whileabs(x-p)>=epsp=x;x=G*p+f;n=n+1;endn%n为迭代次数的记录例子A=[8,-3,2;4,11,-1;6,3,12]b=[20;33;36]p=[0;0;0]x

6、=gaussseidel(A,b,p,0.00001)答案与上一样第三章非线性方程组的数值解法3.1求非线性方程实根的对分法二分法求根clearall;f=@(x)x^3-exp(-x);%这是匿名函数，@是函数句柄，意思是这个函数已经被写进程序中a=0;b=1;%这个函数是f=x^3-e^(-x)epsin=1e-8;iff(a)*f(b)<0x1=(a+b)/2;endwhileabs(x1-b)>epsiniff(a)*f(x1)<0;b=x1;x1=(a+b)/2;elsea=x1;x1=(a+b)/2;endendxstar=x1;disp(x1);这

7、个不需要调用，直接运行就有结果了3.2单个非线性方程的迭代法注：迭代格式必须是收敛的epsin=1e-5;x0=3.5;x1=1/2*(log(x0)+7);k=0;whileabs(x1-x0)>=epsinx0=x1;x1=1/2*(log(x0)+7);disp(x1);x=max(x1);k=k+1;ifk==100breakdisp('nosolution');endenddisp(x);3.3单个非线性方程的Newton法Newton法书上例子P101，x-e^(x)-2=0epsin=1e-5x0=0;x1=x0-(x0+exp(x0)-2)/(1

8、+exp(x0));k=