数值计算方法试验指导Matlab版

1、数值计算方法实验指导版）(Matlab数值计算方法实验指导（Matlab版）肇庆学院数学与统计学学院计算方法课程组数值计算方法实验1报告班20xx级XXXXx班学20xx2409xxxx级：号：成绩:姓XXX名：1 .实验名称实验1算法设计原则验证(之相近数相减、大数吃小数和简化计算步骤)2 .实验题目(1)取z1016,计算7T75和1/(，T7石)，验证123）与两个相近的数相减会造成有效数字的损失.(2)按不同顺序求一个较大的数(1000个较小的数(数的现象.15310)的和，验证大数吃小(3)分别用直接法和秦九韶算法计算多项式n_n1P(x)a0xa1xan1xan在x=1.00037

2、处的值.验证简化计算步骤能减少运算时间.对于第(3)题中的多项式P(x),直接逐项计算需要n(n1)21:1次乘法和n次加法，使用秦九韶算法P(x)(a0Xai)xa2)xani)xan则只需要n次乘法和n次加法.3 .实验目的验证数值算法需遵循的若干规则.4 .基础理论设计数值算法时，应避免两个相近的数相减、防止大数吃小数、简化计算步骤减少运算次数以减少运算时间并降低舍入误差的积累.两相近的数相减会损失有效数字的个数，用一个大数依次加小数，小数会被大数吃掉，乘法运算次数太多会增加运算时间.5 .实验环境操作系统：Windowsxp;程序设计语言：Matlab6 .实验过程(1)直接计算并比较

3、；(2)法1:大数逐个加1000个小数，法2:先把1000个小数相加再与大数加；(3)将由高次项到低次项的系数保存到数组An中，其中n为多项式次数.7 .结果与分析(1)计算的vT7<z=,1/(0_1Jz).分析：(2) 123逐次加1000个3106的和是,先将1000个3106相加，再用这个和与123相加得.分析：(3)计算次的多项式：直接计算的结果是,用时;用秦九韶算法计算的结果是,用时.分析：8 .附录：程序清单(1)两个相近的数相减.%*%*程序名：ex1_1.m*%*程序功能：验证两个相近的数相减会损失有效数字个数*%*z=1e16;%x=d”VZ;%y=1/kVZ)；x,

4、y(2)大数吃小数%*%*ex12.m*%*程序功能：验证大数吃小数的现象.*%*clc;%清屏clearall;%释放所有内存变量formatlong;%按双精度显示浮点数z=123;%大数t=3e-15;%小数x=z;%大数依次加小数%重复1000次给x中加上ty=0;%先累加小数%重复1000次给y中加上ty=z+y;%再加到大数x,y=(3)秦九韶算法%*%*程序名：ex1_3.m*%*程序功能：验证秦九韶算法可节省运行时间.*%*clc;*%清屏clearall;formatlong;%释放所有内存变量%按双精度显示浮点数A=8,4,-1,-3,6,5,3,2,1,3,2,-1,4,

5、3,1,-2,4,6,8,9,50,-80,12,35,7,-6,42,5,6,23,74,65,55,80,78,77,98,56;A(10001)=0;的都是分量0%扩展到10001项，后面%A为多项式系数，从高次项到低次项%n为多x=1.00037;n=9000;项式次数%直接计算begintime=clock;%开始执行的时间的i次品%累加多项式的i次项endtime=clock;%结束执行的时间time1=etime(endtime,begintime);%运行时间disp('直接计算');disp('p(',num2str(x),')=

6、9;,num2str(p);disp('运行时间：,num2str(time1),'秒');%秦九韶算法计算begintime=clock;%开始执行的时间p=Ci；fori0:ap=x*p-KA(i+1);end%累加秦九韶算法中的一项endtime=clock;%结束执行的时间time2=etime(endtime,begintime);%运行时间disp(',);disp('秦九韶算法计算,);disp('p(',num2str(x),')=',num2str(p);disp('运行时间：',num2

7、str(time2),秒1)；10数值计算方法实验1报告班20xx级XXXXx班学20xx2409xxxx姓XXX名：号:成绩:1.实验名称2.实验题目计算定积分；xnex1dx,n0,1,10）分别用教材例实验1算法设计原则验证（之数值稳定性）1-7推导出的算法A和B,其中：算法A：100M1算法B：1In1(1In)nI100验证算法不稳定时误差会扩大.3 .实验目的验证数值算法需遵循的若干规则.4 .基础理论设计数值算法时，应采用数值稳定性好的算法.数值稳定的算法，误差不会放大，甚至会缩小；而数值不稳定的算法会放大误差.115 .实验环境操作系统：Windowsxp;程序设计语言:Mat

8、lab6 .实验过程分别用数组IA和旧口保存两种算法计算的结果.分析:7 .结果与分析n算法A算法B精确值012345678910运行结果：（或拷屏）128 .附录：程序清单%*%*程序名：ex1_4.m*%*程序功能：验证数值稳定性算法可控制误差.*%*clc;%清屏clearall;%释放所有内存变量formatlong;%按双精度显示浮点数I=0.63212055882856,0.36787944117144,0.26424111765712,0.20727664702865,.0.17089341188538,0.14553294057308,0.12680235656154,0.11

9、238350406938,.0.10093196744492,0.09161229300662,130.08387707010843;%保留14位小数的精确值，是Matlab中的续行符%算法AIA(1)=0.6321;%Matlab下标从1开始，所以要用IA(n+1)表示原问题中的I(n)forn-1:10end%算法BIBflij'O:forn=10:-I:1=(I-3(口-ruenddisp('n算法A算法B精确值')；forn=1:11fprintf('%2d%14.6f%14.6f%14.6fn',n-1,IA(n),IB(n),I(n);14e

10、nd%n显示为2位整数,其它显示为14位其中小数点后显示6位的小数15数值计算方法实验1报告班20xx级XXXXx班学20xx2409xxxx级：号：成绩:姓XXX名：1 .实验名称实验1算法设计原则（除数绝对值不能太小）2 .实验题目将线性方程组增广矩阵利用初等行变换可化为aiiai2a21a22由此可解得bib2Xir2皆1aiia12b1ri会口aii00a22'b2'b'0a22'b2bi'/aii，x2b2'/a22.分别解增广矩阵为1021612和i02612的方程组，验证除数绝对值远小于被除数绝对值的除法会导致结果失真.3 .实验目

11、的验证数值算法需遵循的若干规则.4 .基础理论设计数值算法时，应避免除数绝对值远小于16被除数绝对值的除法，否则绝对误差会被放大，使结果失真.5 .实验环境操作系统：Windowsxp;程序设计语言：Matlab6 .实验过程用二维数组A和B存放方程组的增广矩阵，利用题目所给初等行变换求解方程组.7 .结果与分析第1种顺序的方程组的解为x=,y=;第2种顺序的方程组的解为x=,y=.分析：8 .附录：程序清单17*%*程序名：ex15.m*%*程序功能：验证除数的绝对值太小可能会放大误差.*%*clc;A=1e-16,1,1;2,1,2;B=2,1,2;1e-16,1,1;%增广矩阵%方程组A

12、AQ3+：);m=A(t2)/Aa3;AQjHAVm*A(Z工81.为=川3)/A(L1)；AQ3)=用23)AR2工%m=-a_21/a_11是第2行加第1行的倍数%消去a_21%m=-a_12/a_22是第1行加第2行的倍数%消去a_12,系数矩阵成对角线18%未知数x1的值%未知数x2的值disp('方程组A的解：x1=',num2str(A(1,3),x2=',num2str(A(2,3);%万程组B1);BC2.=:)+m*B(k:);m2)B(2,2);5(1.=+B(la1);Bf2T2).disp(',);%m=-b_21/b_11是第2行加第1

13、行的倍数%消去b_21%m=-b_12/b_22是第1行加第2行的倍数%消去b_12,系数矩阵成对角线%未知数x1的值%未知数x2的值disp('方程组B的解：x1=',num2str(B(1,3),x2=',num2str(B(2,3);19数值计算方法实验2报告班20xx级XXXXx班学20xx2409xxxx级：号：姓XXX成名：绩：1 .实验名称实验2非线性方程的迭代解法(之简单迭代法)2 .实验题目用简单迭代法求方程x34x2100在区间1,2内的一个实根，取绝对误差限为104.3 .实验目的掌握非线性方程的简单迭代法.4 .基础理论简单迭代法：将方程f(x)

14、0改写成等价形式x(x),从初值x0开始，使用迭代公式xk1(xlj可以得到一个数列，若该数列收敛，则其极限即为原方程的解.取数列中适当的项可作为近似解.5.实验环境操作系统：Windowsxp;程序设计语言：Matlab206.实验过程7.结果与分析8.附录：程序清单21数值计算方法实验2报告班20xx级XXXXx班学20xx2409xxxx级：号：成绩:姓XXX名：1 .实验名称(之Newton实验2非线性方程的迭代解法迭代法)2 .实验题目100在区间104.用Newton迭代法求方程x34x21,2内的一个实根，取绝对误差限为3.实验目的迭代法.4.基础理论Newton代公式为xk1掌

15、握求解非线性方程的Newtonf'(Xk)迭代法：解方程f(x)0的Newton迭xk3.5 .实验环境操作系统：Windowsxp;程序设计语言:Matlab6 .实验过程227.结果与分析8.附录：程序清单23数值计算方法实验2报告班20xx级XXXXx班学20xx2409xxxx号:姓XXX成名：绩：1 .实验名称实验2非线性方程的迭代解法(之对分区间法)2 .实验题目用对分区间法求方程x3x10在区间1,1.5内的一个实根，取绝对误差限为104.3 .实验目的掌握求解非线性方程的对分区间法.4 .基础理论对分区间法：取a,b的中点p,若f(p)个或b-a<&则p为

16、方程f(x)0的近似解；若f(a)f(p)<0,则说明根在区间取a,p中；否则，根在区间取p,b中.将新的有根区间记为a，b1,对该区间不断重复上述步骤，即可得到方程的近似根.5 .实验环境24操作系统：Windowsxp;程序设计语言：Matlab6 .实验过程用宏定义函数f(x);为了循环方便，得到的新的有根区间始终用a,b表示;由于新的有根区间可能仍以a为左端点，这样会反复使用函数值f(a),为减少运算次数，将这个函数值保存在一个变量fa中；同样在判断新的有根区间时用到函数值f(p),若新的有根区间以p为左端点，则下一次用到的f(a)实际上就是现在的f(p),为减少运算次数，将这个

17、函数值保存在一个变量fp中.算法的伪代码描述：Input:区间端点a,b;精度要求(即误差限)邑函数f(x);最大对分次数NOutput:近似解或失败信息行号伪代码注释1n-1;对分次数计数器2faf(a);左端点的函数25值3whilen<Ndo4p(a+b)区间中点/2;5fpf(P)中点的函数值6iffp=0or(b-a)/2<then函数值为0或半区间长不超£7returnp；输出近似解并退出程序8endif9nn+1;计数器加一10iffafp>0then若中点与左端点函数值同号11a-p新区间取右半区间12fa-fp13else否则14b-p新区间取左右

18、半区间15endif16enddo17return错误信息输出错误信息26并结束程序7 .结果与分析8 .附录：程序清单说明：源程序中带有数字的空行，对应着算法描述中的行号%*%*程序名：Bisection.m*%*程序功能：使用二分法求解非线性方程*%*f=inline('xA3-x-1');%定义函数f(x)a=input（'有根区间左端点:a='）;27b=input('右端点:b=');epsilon=input('误差限:epsilona=');N=input('最大对分次数：N=');1%对分次数计数器

19、n置12%左端点的函数值给变量fafprintf('nkpf(p)a(k)f(a(k)1);fprintf('b(k)b-an');%显示表头fprintf('%2d%36.6f%12.6f%12.6f%12.6fn',0,a,fa,b,b-a);%占2位其中0位小数显示步数0,共12位其中小数6位显示各值%whilen<N28区间中点p5%求p点函数值给变量fpfprintf('%2d%12.6f%12.6f',n,p,fp);%输出迭代过程中的中点信息p和f(p)6%如果f(p)=0或b-a的一半小于误差限efprintf(&#

20、39;nn近似解为:%fn',p);%则输出近似根p(7)return;%并结束程序(7)89%计数器加110%若f(a)与f(p)同号1129%则取右半区间为新的求根区间,即a取作P12%保存新区间左端点的函数值13%否则14%左半区间为新的求根区间，即b取作p15fprintf（'%12.6f%12.6f%12.6f%12.6fn',a,fa,b,b-a）；%显示新区间端点及左端函数值、区间长度16fprintf（'nn经过%d次迭代后未达到精度要求.n',N）;%输出错误信息（行17）30数值计算方法实验2报告班20xx级XXXXx班学20xx24

21、09xxxx级：号：姓XXX成名：绩：1.实验名称实验2非线性方程的迭代解法(之Aitken-Steffensen力口速法)2.实验题目用Aitken-Steffensen加速法求方程x34x2100在区间1,2内的一个实根，取绝对误差限为104.3 .实验目的熟悉求解非线性方程的Aitken-Steffensen加速法.4 .基础理论将方程f(x)0改写成等价形式x(x),得到从初值x0开始的迭代公式xk1(xk)后，基于迭代公式xk1闻)的Aitken-Steffensen加速法是通过(ykzjZk2ykxk迭代-再迭代-加速”完成迭代的，具体过程为yk(xk),Zk(yk),xk1xk3

22、15 .实验环境操作系统：Windowsxp;程序设计语言：Matlab6 .实验过程为了验证Aitken-Steffensen加速法可以把一些不收敛的迭代加速成迭代收敛，我们使用将方程组变形为X3而工3,取迭代函数（x）3万耳,并利用宏定义出迭代函数.由于不用保存迭代过程，所以用X0表示初值同时也存放前一步迭代的值，y和z是迭代过程中产生的yk和zk,x存放新迭代的结果.算法的伪代码描述：Input:初值x。；精度要求（即误差限）e;迭代函数0（x）;最大迭代次数NOutput:近似解或失败信息行号伪代码注释1n-1;迭代次数计数器2whilen<Ndo3y&x。）；迭代324

23、z-&y);再迭代5xx(y-xo)2/(z-2y+xo)加速6if|xxo|<£then如果达到精居要求7returnx;则输出近似4并退出程序18endif9nn+1;计数器加一10xox;新近似值给xo做下次的初值11enddo12return错误信息输出错误信卡并结束程序d7 .结果与分析8 .附录：程序清单%*%*程序名：Aitken_Steffensen.m33%*程序功能：用Aitken-Steffensen加速法求方程.*%*clc;clearall;phi=inline('0.5*sqrt(10-xA3)');%迭代函数x0=input

24、('初值：x0=');epsilon=input('误差限:epsilon=');N=input('最大迭代次数：N=');disp('n迭代中间值y(n-1)再迭代结构z(n-1)加速后的近似值x(n)');fprintf('%2d%54.6fn',0,x0);%占2位整数显示步数0,为了对齐占54位小数6位显示x01%n是计数器34%whilen<=Ny=3;%迭代z=3;%再迭代x=3;%力口速%x0初值及前一步的近似值，y和z是中间变量，x是下一步的近似值fprintf('%2d%18.6f%

25、18.6f%18.6fn',n,y,z,x);%显示中间值和迭代近似值6%如果与上一步近似解差的绝对值不超过误差限fprintf('nn近似解x=x(%d)=%f,n,x);%则输出近似根(7),可简35略为:fprintf('nn近似解x=%f',x);return;%并结束程序(7)8%相当于endif9%计数器加110%新近似值x作为下一次迭代的初值11fprintf('n迭代%d次还不满足误差要求.nn',N);%输出错误信息(12)36数值计算方法实验2报告班20xx级XXXXx班学20xx2409xxxx级：号：姓XXX成名：绩：1

26、.实验名称实验2非线性方程的迭代解法(之Newton下山法)2 .实验题目用Newton下山法求方程x34x2100在区间2 1,2内的一个实根，取绝对误差限为104.3 .实验目的熟悉非线性方程的Newton下山法.4 .基础理论Newton下山法：Newton下山法公式为xk1xkkf(xk),使If(xk1)IIf(xk)I,其中0k1.f'(xk)5 .实验环境操作系统：Windowsxp;程序设计语言：Matlab6 .实验过程37定义函数f(x)和df(x),其中df(x)是f(x)的导函数.每步迭代时先取下山因子为1,尝试迭代，判断尝试结果是否满足下山因子，若满足则作为这

27、步的迭代结果；否则将下山因子减半，然后再尝试.为防止当前的xk是极小值点，附近不会有满足下述条件的其它点，使尝试陷入死循环，同时计算机中能表示出的浮点数也有下界，因此我们设置了最大尝试次数.当超过最大尝试次数时，不再进行下山尝试.由于反复尝试迭代且要判断下山条件，所以f(x0)和f乂0)会反复使用，为避免重复计算浪费运行时间，将这两个值分别保存在变量fx0和dfx0.而尝试产生的节点，判断下山条件时要用到它的函数值，若尝试成功，这个点会作为下一步的初值再使用，所以把该点的函数值也保存在变量fx中.算法的伪代码描述：mput:初值xo；精度要求(即误差限)£；函数及其导函数f(x)和f

28、'(x)最大迭代次数N;K下山尝试最大次数38Output:近似解或失败信息行号伪代码注释1n-1;迭代次数计娄器攵2Fo<f(xo);3whilen<Ndo4Fo-f0);(x5ifF。'=0then若该点导数为06returnFalse;则无法进行关代，结束程序d7endif8lambda_1;下山因子入从试起19k-1尝试次数计数器10whilek<Kdo(11xX0-lambda*F0/F0'Newton下山式口12Fxf(x);13if|Fx|<|F0|then判断下山条件14break;退出尝试循环3915endif16lambda一

29、lambda/2;下山因子减半17kk+1;尝试次数计数器加118endwhile19ifk>Kthen如果因超过尝试次数退出循环20returnFalse;则提示错误信息并结束程序21endif否则时尝试成功退出上边循环22if|xX0|<£then如果达到精度要求23returnx;则输出近似值并退出程序24endif25X0x;新近似值给x。做下次的初值26F。Fx;所求函数值卜次也用到4027nn+1;计数器加一28enddo29return错误信息输出错误信息并结束程序7 .结果与分析8 .附录：程序清单%*%*程序名：NewtonDownhill.m*%*程序

30、功能：用Newton下山法求解非线性方程.*%*clc;clearall;f=inline('xA3-x-1');41%函数f(x)df=inline('3*xA2-1');%函数f(x)的导函数x0=input('初值：x0=');epsilon=input('误差限:epsilon=');N=input('最大迭代次数：N=');K=input('最大下山尝试次数：K=');1 %迭代次数计数器2 %存x0点函数值fprintf('nnnx(n)f(x(n)n');%显示表头fp

31、rintf('%2d%14.6f%14.6fn',0,x0,fx0);%2位整数显示0,共14位小数6位显示x0和fx042whilen<Ndisp(',);程的表头disp('%换行显示下山尝试过下山因子尝试满足下山条件);%存x0点x(n)对应f(x(n)disp('');导数值，每次下山尝试不用重新计算ifdfx0=0%导数为0不能迭代disp(无法进行Newton迭代');return;endlambda=1.0;%下山因子从1开始尝试k=1;k下山尝试次数计数器whilek<=K43下山最多尝试K次x=xO-lambda*&0dfxO;%下山公式fx=f(x);%函数值fprintf('%22.6f%14.6f%14.6f',lambda,x,fx);%显示尝试结果if(abs(fx)<abs(fx0)%判断是否满足下山条件fprintf('满足n');break;%是，则退出下山尝试的循环elsefprintf('不满足n');endlambda=lambda/2;%不是，则下山因子减半k=k+1;%计数器加1endifk>Kfpr