19次方系数对方程根的影响课程设计报告书.docx

1、19次方系数对方程根的影响课程设计报告书西安郵電學院数值分析课程设计报告书系部名称：应 用 数 理 系学生姓名：专业名称：信息与计算科学班 级： 时间：2 实验X 19次方系数对方程根的影响一、问题的提出考虑代数多项式,显然它的全部根为1,2,20;如果19次方的系数有扰动,那么对方程的解有何影响.二、实验内容直接利用Matlab中的roots()和poly()函数,选择不同的扰动做数值计算.主要步骤:ve=zeros(1,21);ve(2)=ess;roots(poly(1:20)+ve).Matlab程序:见附件bingtaiwenti.m三、实验结果1.取扰动分别为所得结果如下:27.0

2、817+5.03812i21.3025+1.56717i19.952+0i27.0817-5.03812i21.3025-1.56717i19.2293+0i19.5337+9.1664i18.5028+3.6004i17.6573+0.692896i19.5337-9.1664i18.5028-3.6004i17.6573-0.692896i13.8235+7.77167i15.1651+3.76125i15.4524+0.875524i13.8235-7.77167i15.1651-3.76125i15.4524-0.875524i10.7211+5.4609i12.4866+2.8827

3、8i13.3527+0.486992i10.7211-5.4609i12.4866-2.88278i13.3527-0.486992i8.91282+3.47317i10.5225+1.71958i11.8578+0i8.91282-3.47317i10.5225-1.71958i11.0427+0i7.69268+1.89884i9.04487+0.59455i9.9916+0i7.69268-1.89884i9.04487-0.59455i9.00201+0i6.75761+0.654714i7.94891+0i7.99952+0i6.75761-0.654714i7.00247+0i7.

4、00009+0i5.95208+0i5.99995+0i5.99999+0i5.00061+0i5+0i5+0i4+0i4+0i4+0i3+0i3+0i3+0i2+0i2+0i2+0i1+0i1+0i1+0i2.绘制相应的图形,即将以上结果可视化.四、实验结果分析观测现象:由结果可以观察到误差是不可避免的,并且扰动的减小不能使得所有根除的误差减小,还可以观察到处的误差总是最大,而靠近处的误差较小.误差分析:原问题数学上描述为:的求解问题,也可等价为:那么要考虑的问题是对方程的某个解的影响,我们不妨将看成的函数,从而可以求得此问题的条件数为:可以求得20个根处的条件数为:19次方系数变化对根的影

5、响问题条件数X=1,5X=6,10X=11,15X=16,20.172633e-142038.69.886669e9.315986e11.859841e-876273.1.347381e10.235250e11.114368e-3.156701e7.978210e10.116152e11.114955.195843e8.199945e11.341541e1025.5252.159475e9.296669e11.452548e9可见而靠近处的条件数非常小,而处的条件数最大.这与我们的观测结果是一致的.由于条件数太大所以此问题是病态问题.实验Y 用Jacobi法求对称矩阵的特征值及特征向量一、实验

6、内容已下列矩阵为例,求对称矩阵的全部特征值及特征向量:1. 2. 3. 二、方法步骤1.在的非主对角线元素中,找出按模最大的元素;2.计算平面旋转矩阵,其中的及由计算;3.计算, (的初始值取单位阵);4.如果(其中为的元素),则停止计算,所求特征值为:,特征向量:即得第列;否则令,重复以上各步.三、实验结果讨论的矩阵为:A = 5 2 4 2 5 3 4 3 1要求误差为:err=1e-005迭代次数:7第1个特征值为:9.8057相应的特征向量为:0.65053 0.57817 0.49248第2个特征值为:3.0604相应的特征向量为:-0.62538 0.77572 -0.084615

7、第3个特征值为:-1.8661相应的特征向量为:-0.43095 -0.25294 0.8662讨论的矩阵为:A = 3 5 -1 5 2 9 -1 9 3要求误差为:err=1e-006迭代次数:7第1个特征值为:3.8441相应的特征向量为:0.87749 0.052794 -0.47668第2个特征值为:-8.2319相应的特征向量为:-0.3691 0.70896 -0.60095第3个特征值为:12.3879相应的特征向量为:0.30622 0.70327 0.64159讨论的矩阵为:A = 2 1 3 1 2 6 3 6 1要求误差为:err=1e-009迭代次数:8第1个特征值为

8、:1.1984相应的特征向量为:0.88975 -0.4478 -0.088467第2个特征值为:8.6956相应的特征向量为:0.39037 0.64607 0.6559第3个特征值为:-4.894相应的特征向量为:-0.23656 -0.61812 0.74965以上为程序:tezhjcobi.m运行结果,Matlab程序见附件数值分析试验报告矩阵的LU分解1. 题目：求4阶矩阵的LU分解2. 方法：杜里特尔分解法3. 程序：function f=LU_decom(A)m,n=size(A)L=eye(n);U=zeros(n);flag=ok;for i=1:n U(1,i)=A(1,i

9、); endfor r=2:n L(r,1)=A(r,1)/U(1,1); end for i=2:n for j=i:n z=0; for r=1:i-1 z=z+L(i,r)*U(r,j); end U(i,j)=A(i,j)-z; end if abs(U(i,i) LU_decom(A)m = 4n = 4L = 1 0 0 0 2 1 0 0 1 2 1 0 3 3 2 1U = 2 4 2 6 0 1 2 3 0 0 3 6 0 0 0 15.拓展在编写程序过程中由于角标较多因此在运行过程中出现了不少角标不对的错误题目；给出函数f(x)=1/(1+25x2)，求f（x在-1,1上取

10、5个和9个等距节点，做最小二乘拟合，得出均方误差方误差。五个节点时，matlab编码为：首先建立M文件，并保存function y=f(x)y=1/(1+25*x2);endx=-1 -0.5 0 0.5 1;for i=1:5y(i)=f(x(i);enda=polyfit(x,y,3)syms xf1=a(1)*x3+a(2)*x2+a(3)*x+a(4) x=-1 -0.5 0 0.5 1;for i=1:5E(i)=(f(x(i)-(a(1)*x(i)3+a(2)*x(i)2+a(3)*x(i)+a(4)2;endsum(E)输出结果为a = -0.0000 -0.6063 -0.00

11、00 0.5737 f1 = -4878849915647781/1298074214633706907132624082305024*x3-1600/2639*x2-5348847520430703/649037107316853453566312041152512*x+1514/2639 （拟合的多项式） ans = 0.3534（均方误差）九个点的时候，matlab编码为：x=-1 -0.75 -0.5 -0.25 0 0.25 0.5 0.75 1;for i=1:9y(i)=f(x(i);enda=polyfit(x,y,3)syms xf2=a(1)*x3+a(2)*x2+a(3)

12、*x+a(4) x=-1 -0.75 -0.5 -0.25 0 0.25 0.5 0.75 1;for i=1:5E1(i)=(f(x(i)-(a(1)*x(i)3+a(2)*x(i)2+a(3)*x(i)+a(4)2;endsum(E1)输出结果为：a = -0.0000 -0.5609 0.0000 0.4855 f2 = -728732707776039/2535301200456458802993406410752*x3-1263030908712921/2251799813685248*x2+4915246442354361/20282409603651670423947251286

13、016*x+1093229300764671/2251799813685248 （最小二乘拟合多项式）ans =0.3350（均方误差） 用复合梯形公式求积分的值。function i=combinetraprl(f,a,b,eps)%复化梯形公式求函数f在区间a，b上的定积分%函数名：f%积分下限：a%积分上限：b%积分精度：eps%积分值：i%积分划分的子区间个数：stepn=1;h=(b-a)/2;i1=0;i2=(subs(sym(f),findsym(sym(f),a)+subs(sym(f),findsym(sym(f),b)/h;while abs(i2-i1)eps n=n+1

14、; h=(b-a)/n; i1=i2; i2=0; for i=0:n-1 x=a+h*i; x1=x+h; i2=i2+(h/2)*(subs(sym(f),findsym(sym(f),x)+. subs(sym(f),findsym(sym(f),x1); endend四阶龙格-库塔法分别求解下列初值问题；functionx,y=runge_kutta(fun,x0,xt,y0,pointnum,varargin)if nargin3 y0=0;endy(1,:)=y0(:);h=(xt-x0)/(pointnum-1);x=x0+0:pointnum*h;for k=1:pointnu

15、m f1=h*feval(fun,x(k),y(k,:),varargin:); f1=f1(:); f2=h*feval(fun,x(k)+h/2,y(k,:)+f1/2,varargin:); f2=f2(:); f3=h*feval(fun,x(k)+h/2,y(k,:)+f2/2,varargin:); f3=f3(:); f4=h*feval(fun,x(k)+h,y(k,:)+f3,varargin:); f4=f4(:); y(k+1,:)=y(k,:)+(f1+2*(f2+f3)+f4)/6;endfunctionx,y=euler(fun,x0,xt,y0,pointnum)

16、if narginx0=0;xt=1;fun=inline(9*y/(1+2*x),x,y);y0=1;pointnum=10;x1,yr=runge_kutta(fun,x0,xt,y0,pointnum)一 运行结果：1） Runge_kutta方法：x1 =0 0.1111 0.2222 0.3333 0.4444 0.5556 0.6667 0.7778 0.8889 1.0000 1.1111yr =1.0000 2.4611 5.2136 9.9212 17.4198 28.7302 45.0700 67.8650 98.7600 139.6286 192.5834小行星轨道问题：

17、一天文学家要确定一颗小行星绕太阳运行的轨道，他在轨道平面内建立以太阳为原点的直角坐标系，在五个不同的对小行星作了五次观察，测得轨道上五个点的坐标数据（单位：万公里）如下表所示：P1P2P3P4P5X坐标5360558460628596666268894Y坐标602611179169542349268894由开普勒第一定律知，小行星轨道为一椭圆，椭圆的一般方程可表示为：a1x2+2a2xy+a3y2+2a4x+2a5y+1=0现需要建立椭圆的方程以供研究。（1）分别将五个点的数据代入椭圆一般方程中，写出五个待定系数满足的等式，整理后写出线性方程组AX=R以及方程组的系数矩阵和右端项b；（2）用M

18、ARLAB求低阶方程的指令Ab求出待定系数a1,a2,a3,a4,a5.答：1)2)一源程序：A=536052 2*53605*6026 60262 2*53605 2*6026; 584602 2*58460*11179 111792 2*58460 2*11179; 628592 2*62859*16954 169542 2*62859 2*16954; 666622 2*66662*23492 234922 2*66662 2*23492; 688942 2*68894*68894 688942 2*68894 2*68894;b=-1,-1,-1,-1,-1;x=Ab二运行结果： 数值

19、分析课程设计报告 摘要（中）：本文建立在数值分析的理论基础上，能够在Matlab环境中运行，给出了理论分析、程序清单以及计算结果。更重要的是，还有详细的对算法的框图说明。首先运用Romberg积分方法对给出定积分进行积分，然后对得到的结果用插值方法，分别求出Lagrange插值多项式和Newton插值多项式，再运用最小二乘法的思想求出拟合多项式，最后对这些不同类型多项式进行比较，找出它们各自的优劣。Summary(Enlish):This essay is based on numerical analysis ,and could be operated in Matlab environm

20、ent, including theroy analysis, program and results. Whats more, there are detailed diagram which shows how the algorithm works. It first uses the integrating method of Romberg , which is an improved trapezoidal integration, to solve the given definite integral,then we create Lagranges interpolation

21、 polynomial and Newtons interpolation polynomial. And according to least square method, curve fitting polynomial is created. At the last part of the essay, we compare these different patterns of polynomial, founding their distinctive advantages and disadvantages.主题词：Romberg积分，插值方法，Langrange插值多项式，New

22、ton插值多项式，拟合多项式。一. 问题提出： 已知椭圆的周长可以表示成s=a (01),取a=1, 针对从0.1到0.9（步长h=0.1）分别求出周长s;(用Romberg积分方法)对于以上数据，求出的插值多项式；对于中数据，试用最小二乘法的思想求作拟合多项式（要求是偶次），并对这些多 项式的优劣进行比较。二问题解决： 依据问题出现先后顺序，对三个小问进行讨论。 Romberg算法是在复化梯形求积公式的基础上，应用理查逊外推构造的一种数值积分方法。 由复合梯形公式的展开定理，得到如下关系式：T1(h)-I=aah2+a2h4+a3h6+amh2m+其中，I=,T1(h)=Tn.利用Richardson外推定理对T1(h)进行加速，注意这里取m=1,q=,有利用T0()和T0()可以得到实际上T1(h)就是复化抛物线求积公式，一般的计算公式(m=1,2, k=0,1,2,) 即.由于Romberg求积过程是每次把区间缩小一半，所以Romberg积分方法也叫做逐次分半加速收敛法。Robermg求积算法的计算过程如下：取k=0,h=b-a,求= 令1k(k记区间a,b的二分次数).求梯形值T0，即按递推公式计算.求加速值，按公式逐个求出如表的第k行其余各元素Tj(k-j)(j=1,2,k).若0.0001&ik)|i4 i=i+1; h=h/2;