当前位置:文档之家› 常州大学数值分析作业(共六章)

常州大学数值分析作业(共六章)

第一章:9.设2cos 1)(x xx f -=,给出计算函数值)012.0(f 的一个合适算法,并在字长m 给定的,十进制计算机上给出数值计算结果。

解:由 )2421(242)2421(1)cos(1224242x x x x x x x -=-=+--≈- 得 )2421(cos 1)(22x x x x f -≈-=10. 字长为5的十进制计算机上计算)015.0(f 和)015.0(g ,并与)015.0(f 的精确值1.0075376410479比较,说明差异存在理由,其中x e x f x 1)(-=,24621)(32x x x x g +++=。

clearf=@(x)1/2-x^2/24; f(0.012)ans =0.5000解:字长为5时的误差很大,这是因为设置的字长有限,就不可避免的使舍入误差不断积累。

把字长改为9时,误差已经大幅度减小。

这说明,加大字长可以显著减小误差。

11. 举例介绍数组矩阵常见运算。

解:举例如下clearf=@(x)digit(digit(exp(x)-1,5)/x,5);g=@(x)digit(digit(1,5)+digit(x/2,5)+digit... (digit(x^2,5)/6,5)+digit(digit(x^3,5)/24,5),5); exc=1.0075376410479; f(0.015) g(0.015)err1=f(0.015)-exc err2=g(0.015)-excans =1.0075 ans =1.0075 err1 =-3.7641e-05 err2 =-3.7641e-05A*Bans =30 30 30 30 70 70 70 70 110 110 110 110 150 150 150 150A.*B ans =1 2 3 4 10 12 14 16 27 30 33 36 52 56 60 64 A^2 ans =A.^2 ans =f=@(x)digit(digit(exp(x)-1,9)/x,9);g=@(x)digit(digit(1,9)+digit(x/2,9)+digit... (digit(x^2,9)/6,9)+digit(digit(x^3,9)/24,9),9); err1=f(0.015)-exc err2=g(0.015)-excerr1 =-1.0479e-09 err2 =-1.0479e-09clearA=[1:4;5:8;9:12;13:16]B=[1,1,1,1;2,2,2,2;3,3,3,3;4,4,4,4] A ’ A =1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 B =1 1 1 12 2 2 23 3 3 34 4 4 4 ans =%%编写m 文件使用digit 函数设置字长%% function y=digit(x,m) k=max(size(x)); y=x;for i=1:k if x(i)<0 sign=-1; elsesign=1; endx(i)=abs(x(i)); p=0;if x(i)<0.1&x(i)>eps while x(i)<0.1 x(i)=x(i)*10; p=p-1; end endif x(i)>=1while x(i)>=1 x(i)=x(i)/10; p=p+1; end endy(i)=round(x(i)*10^m)/10^m; y(i)=sign*y(i)*10^p; end return12.对任意给定的实数a 、b 、c 、试编写Matlab 程序,求方程02=++c bx ax 的根。

解:利用教材例11的方法: 当b>0时,a ac b b x 2421---=,bac b c x +--=4222。

当b<0时,bac b cx --=4221,a acb b x 2422-+-=。

13.利用1,753arctan 753<+-+-=x x x x x x 及()3/3arctan 6=π,给出一个计算π的方法,根据此方法编写程序,给出π的至少有10位有效数字的近似值。

解:根据题中所给公式,容易得到:A+B ans = 2 3 4 5 7 8 9 10 12 13 14 15 17 18 19 20 A-B ans =0 1 2 3 3 4 5 6 6 7 8 9 9 10 11 12cleara=input('输入a='); b=input('输入b='); c=input('输入c='); d=b^2-4*a*c if b>=0x1=(-b-d^0.5)/2*a x2=-2*c/(d^0.5+b) else b<0x1=-2*c/(d^0.5-b) x2=(-b+d^0.5)/2*a endclearx=3^(-0.5);y=atan(x) ; %精确值% s=0 ; %计算值% for k=1:2:100;s=s+(-1)^((k+1)/2)*(x^k)/k;err=y-s; m=abs(err) if m<=1e-11 输入a=2输入b=2 输入c=1 d = -4 x1 =-2.0000 - 2.0000i x2 =-0.5000 + 0.5000i 输入a=1 输入b=2 输入c=1 d = 0 x1 = -1 x2 = -1 输入a=1 输入b=5 输入c=2 d = 17 x1 =-4.5616 x2 =-0.4384()12)3/3(16)3/3arctan(61211--≈=-=+∑i i ni i π14.分别利用下式给出计算ln2的近似方法,编写相应程序并比较算法运行情况。

11,32)1()1ln(3211≤<-+++-=-=+∑∞=+x nx x x x n x x nn n n 11),1253(212211ln 1253112<<-+-++++=-=-+-∞=-∑x n x x x x n x x x n n n解:由运行结果可知,方法二的绝对误差比方法一的误差要小得多。

这是因为方法一给出的计算公式含有相近数相减项,损失了有效数字。

而方法二给出的计算公式避免了相近数相减,具有较好的精度。

第二章:20.(1)用Jacobi 迭代法解方程组AX=b.function [x,iternum,flag]=jacobi(A,b,x0,delta,max1) %检验输入参数,初始化if nargin<2,error('more augments are needed');end if nargin<3,x0=zeros(size(b));end if nargin<4,delta=1e-13;endif nargin<5,max1=100;endclear x=1/3; s=0;for k=1:2:100; s=s+((x)^k)/k; endy=2*s %计算值% g=log(2) %真实值% err=y-g %绝对误差% y =0.6931 g =0.6931 err =-2.2204e-16clearx=1; s=0;for k=1:100;s=s+((-1)^(k+1))*(x^k)/k; endy=s %计算值% g=log(2) %真实值% err=y-g %绝对误差% y =0.6882 g =0.6931 err =-0.005020.(1)用Gauss-Seidel 迭代法解方程组AX=b.function [x,iternum,flag]=gseid(A,b,x0,delta,max1) %检验输入参数,初始化if nargin<2,error('more augments are needed');end if nargin<3,x0=zeros(size(b));end if nargin<4,delta=1e-13;endif nargin<5,max1=100;endA=[1 2 -2;1 1 1;2 2 1]; b=[1;1;1]; [x,iternum,flag]=jacobi(A,b) 输出结果为: The Jacobi method converges. ans = -3 3 1 iternum =4 flag = 13.(1)将矩阵A 进行Crout 分解 A=[1 2 -2;1 1 1;2 2 1]; b=[1;1;1]; [x,iternum,flag]=gseid(A,b) 输出结果为: The Gauss-seidel method does not converge with 100 iterations ans = 1.0e+31 * -9.1905 9.2222 -0.0634 iternum = 4 flag = 1function [L,U]=Crout(A) %检验输入参数,初始化 b=size(A);n=b(1); if b(1)~=b(2)error('Matrix A should be a Square matrix.'); end3.(2)将矩阵A 进行Doolittle 分解。

A=[1 0 2 0;0 1 1 1;2 0 -1 1;0 0 1 1] A = 1 0 2 0 0 1 1 1 2 0 -1 1 0 0 1 1 [L,U]=Crout(A) L = 1.0000 0 0 0 0 1.0000 0 0 2.0000 0 -5.0000 0 0 0 1.0000 1.2000 U = 1.0000 0 2.0000 00 1.0000 1.0000 1.00000 0 1.0000 -0.20000 0 0 1.0000function [L,U]=Doolittle(A) %检验输入参数,初始化 b=size(A);n=b(1); if b(1)~=b(2)error('Matrix A should be a Square matrix.'); end7.用改进平方根法解方程组AX=bA=[1 0 2 0;0 1 1 1;2 0 -1 1;0 0 1 1]A =1 02 00 1 1 12 0 -1 10 0 1 1[L,U]=Doolittle(A) L = 1.0000 0 0 0 0 1.0000 0 0 2.0000 0 1.0000 0 0 0 -0.2000 1.0000 U = 1.0000 0 2.0000 00 1.0000 1.0000 1.00000 0 -5.0000 1.00000 0 0 1.2000function [x]=improvecholesky(A,b,n) %检验输入参数,初始化L=zeros(n,n);D=diag(n,0);S=L*D; for i=1:nL(i,i)=1;end8.(2)用追赶法求解方程组.function [L,U,x]=pursue(a,b,c,f) n=length(b);u(1)=b(1); for i=2:nif (u(i-1)~=0)l(i-1)=a(i-1)/u(i-1);u(i)=b(i)-l(i-1)*c(i-1); A=[4 1 -1 0;1 3 -1 0;-1 -1 5 2;0 0 2 4] b=[1;0;0;0] A = 4 1 -1 0 1 3 -1 0 -1 -1 5 2 0 0 2 4 b = 1 0 0 0n=4[x]=improvecholesky(A,b,n) n = 4 x =0.2821 -0.0769 0.0513 -0.0256第三章:1.设节点x 0=0,x 1=π/8,x 2=π/4,x 3=3π/8,x 4=π/2,试适当选取上述节点,用Lagrange 插值法分别构造cos x 在区间[0,π/2]上的一次、二次、四次差值多项式P 1(x ),P 2(x )和P 4(x ),并分别计算P 1(π/3),P 2(π/3)和P 4(π/3)。

相关主题