数值计算方法实验指导(Matlab版)

1、数值计算方法实验指导（Matlab版）肇庆学院 数学与统计学学院计算方法课程组数值计算方法实验1报告班级：20xx级XXXXx班学号：20xx2409xxxx姓名：XXX成绩：1. 实验名称 实验1 算法设计原则验证（之相近数相减、大数吃小数和简化计算步骤）2. 实验题目 (1) 取，计算和，验证两个相近的数相减会造成有效数字的损失 (2) 按不同顺序求一个较大的数（123）与1000个较小的数（）的和，验证大数吃小数的现象 (3) 分别用直接法和秦九韶算法计算多项式在x=1.00037处的值验证简化计算步骤能减少运算时间对于第(3)题中的多项式P(x)，直接逐项计算需要次乘法和n次加法，使用

2、秦九韶算法则只需要n次乘法和n次加法3. 实验目的 验证数值算法需遵循的若干规则4. 基础理论 设计数值算法时，应避免两个相近的数相减、防止大数吃小数、简化计算步骤减少运算次数以减少运算时间并降低舍入误差的积累两相近的数相减会损失有效数字的个数，用一个大数依次加小数，小数会被大数吃掉，乘法运算次数太多会增加运算时间5. 实验环境 操作系统：Windows xp； 程序设计语言：Matlab 6. 实验过程 (1) 直接计算并比较； (2) 法1：大数逐个加1000个小数，法2：先把1000个小数相加再与大数加；(3) 将由高次项到低次项的系数保存到数组An中，其中n为多项式次数7. 结果与分析

3、 (1) 计算的= ， 分析： (2) 123逐次加1000个的和是 ，先将1000个相加，再用这个和与123相加得 分析： (3) 计算 次的多项式：直接计算的结果是 ，用时 ；用秦九韶算法计算的结果是 ，用时 分析：8. 附录：程序清单 (1) 两个相近的数相减%*%* 程 序 名：ex1_1.m *%* 程序功能：验证两个相近的数相减会损失有效数字个数 *%*% x = ；% y = ；z=1e16; x,y=(2) 大数吃小数%*%* 程 序 名：ex1_2.m *%* 程序功能：验证大数吃小数的现象. *%*clc; % 清屏clear all; % 释放所有内存变量format l

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

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

6、(' 运行时间: ',num2str(time1),'秒');% 秦九韶算法计算begintime=clock; % 开始执行的时间 % 累加秦九韶算法中的一项endtime=clock; % 结束执行的时间time2=etime(endtime,begintime); % 运行时间disp(' ');disp('秦九韶算法计算');disp('p(',num2str(x),')=',num2str(p);disp(' 运行时间: ',num2str(time2),'秒

7、9;);数值计算方法实验1报告班级：20xx级XXXXx班学号：20xx2409xxxx姓名：XXX成绩：1. 实验名称 实验1 算法设计原则验证（之数值稳定性）2. 实验题目 计算定积分，分别用教材例1-7推导出的算法A和B，其中： 算法A： 算法B：验证算法不稳定时误差会扩大3. 实验目的 验证数值算法需遵循的若干规则4. 基础理论 设计数值算法时，应采用数值稳定性好的算法数值稳定的算法，误差不会放大，甚至会缩小；而数值不稳定的算法会放大误差5. 实验环境 操作系统：Windows xp； 程序设计语言：Matlab 6. 实验过程分别用数组IA 和IB 保存两种算法计算的结果7. 结果与

8、分析 运行结果：（或拷屏）n算法A算法B精确值012345678910分析：8. 附录：程序清单%*%* 程 序 名：ex1_4.m *%* 程序功能：验证数值稳定性算法可控制误差. *%*clc; % 清屏clear all; % 释放所有内存变量format long; % 按双精度显示浮点数I=0.63212055882856, 0.36787944117144, 0.26424111765712, 0.20727664702865, .0.17089341188538, 0.14553294057308, 0.12680235656154, 0.11238350406938, .0.1

9、0093196744492, 0.09161229300662, 0.08387707010843; % 保留14位小数的精确值, 是Matlab中的续行符% 算法AIA(1) = 0.6321; % Matlab下标从1开始，所以要用IA(n+1)表示原问题中的I(n）% 算法Bdisp('n 算法A 算法B 精确值');for n=1:11 fprintf('%2d %14.6f %14.6f %14.6fn',n-1,IA(n),IB(n),I(n);end % n显示为2位整数, 其它显示为14位其中小数点后显示6位的小数数值计算方法实验1报告班级：20

10、xx级XXXXx班学号：20xx2409xxxx姓名：XXX成绩：1. 实验名称 实验1 算法设计原则（除数绝对值不能太小）2. 实验题目将线性方程组增广矩阵利用初等行变换可化为由此可解得分别解增广矩阵为和的方程组，验证除数绝对值远小于被除数绝对值的除法会导致结果失真3. 实验目的 验证数值算法需遵循的若干规则4. 基础理论 设计数值算法时，应避免除数绝对值远小于被除数绝对值的除法，否则绝对误差会被放大，使结果失真5. 实验环境 操作系统：Windows xp； 程序设计语言：Matlab 6. 实验过程用二维数组A和B存放方程组的增广矩阵，利用题目所给初等行变换求解方程组7. 结果与分析 第

11、1种顺序的方程组的解为x= ，y= ；第2种顺序的方程组的解为x= ，y= 分析：8. 附录：程序清单%*%* 程 序 名：ex1_5.m *%* 程序功能：验证除数的绝对值太小可能会放大误差. *%*clc;A=1e-16, 1, 1; 2, 1, 2;B=2, 1, 2; 1e-16, 1, 1; % 增广矩阵% 方程组A % m = - a_21/a_11 是第2行加第1行的倍数 % 消去a_21 % m = - a_12/a_22 是第1行加第2行的倍数 % 消去a_12, 系数矩阵成对角线 % 未知数x1的值 % 未知数x2的值disp('方程组A的解: x1=',n

12、um2str(A(1,3),', x2=',num2str(A(2,3);disp(' ');% 方程组B % m = - b_21/b_11 是第2行加第1行的倍数 % 消去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);数值计算方法实验2报告班级：20xx级XXXXx班学号：20xx2409xxxx姓名：XXX成绩：1

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

14、法(之Newton迭代法)2. 实验题目 用Newton迭代法求方程在区间1,2内的一个实根，取绝对误差限为3. 实验目的 掌握求解非线性方程的Newton迭代法4. 基础理论Newton迭代法：解方程的Newton迭代公式为5. 实验环境 操作系统：Windows xp； 程序设计语言：Matlab6. 实验过程 7. 结果与分析 8. 附录：程序清单数值计算方法实验2报告班级：20xx级XXXXx班学号：20xx2409xxxx姓名：XXX成绩：1. 实验名称 实验2 非线性方程的迭代解法(之对分区间法)2. 实验题目 用对分区间法求方程在区间1, 1.5内的一个实根，取绝对误差限为3.

15、实验目的 掌握求解非线性方程的对分区间法4. 基础理论对分区间法：取a，b的中点p，若f(p) 0或b a < ，则p为方程的近似解；若f(a) f(p) < 0，则说明根在区间取a，p中；否则，根在区间取p，b中将新的有根区间记为 a1，b1，对该区间不断重复上述步骤，即可得到方程的近似根5. 实验环境 操作系统：Windows xp； 程序设计语言：Matlab6. 实验过程用宏定义函数f(x)；为了循环方便，得到的新的有根区间始终用a，b表示；由于新的有根区间可能仍以a为左端点，这样会反复使用函数值f(a)，为减少运算次数，将这个函数值保存在一个变量fa中；同样在判断新的有根

16、区间时用到函数值f(p)，若新的有根区间以p为左端点，则下一次用到的f(a)实际上就是现在的f(p)，为减少运算次数，将这个函数值保存在一个变量fp中算法的伪代码描述：Input：区间端点a，b；精度要求（即误差限）；函数f（x）；最大对分次数NOutput：近似解或失败信息行号伪代码注释1n 1； 对分次数计数器2fa f（a）； 左端点的函数值3while n N do 4 p （a+b）/ 2； 区间中点5 fp f(p)； 中点的函数值6 if fp = 0 or (b-a)/2 < then 函数值为0或半区间长不超7 return p； 输出近似解并退出程序 8 endif

17、9 n n + 1； 计数器加一10if fa ·fp > 0 then 若中点与左端点函数值同号11a p； 新区间取右半区间12fa fp；13else 否则14b p；新区间取左右半区间15endif16enddo17 return 错误信息 输出错误信息并结束程序7. 结果与分析 8. 附录：程序清单说明: 源程序中带有数字的空行，对应着算法描述中的行号%*%* 程 序 名：Bisection.m *%* 程序功能：使用二分法求解非线性方程. *%*f=inline('x3-x-1'); % 定义函数f(x)a=input('有根区间左端点: a

18、=');b=input('右端点:b=');epsilon=input('误差限:epsilona=');N=input('最大对分次数: N='); 1 % 对分次数计数器n置1 2 % 左端点的函数值给变量fafprintf('n k p f(p) a(k) f(a(k)');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位其中小

19、数6位显示各值 3 % while n N 4 % 取区间中点p 5 % 求p点函数值给变量fp fprintf('%2d%12.6f%12.6f',n,p,fp);% 输出迭代过程中的中点信息p和f(p) 6 % 如果f(p)=0或b-a的一半小于误差限 fprintf('nn近似解为:%fn',p); % 则输出近似根p (7)return; % 并结束程序 (7) 8 9 % 计数器加1 10 % 若f(a)与f(p)同号 11 % 则取右半区间为新的求根区间, 即a取作p 12 % 保存新区间左端点的函数值 13 % 否则 14 % 左半区间为新的求根区

20、间, 即b取作p 15 fprintf('%12.6f%12.6f%12.6f%12.6fn',a,fa,b,b-a);% 显示新区间端点及左端函数值、区间长度 16 fprintf('nn经过%d次迭代后未达到精度要求.n',N);% 输出错误信息(行17)数值计算方法实验2报告班级：20xx级XXXXx班学号：20xx2409xxxx姓名：XXX成绩：1. 实验名称 实验2 非线性方程的迭代解法(之Aitken-Steffensen加速法)2. 实验题目 用Aitken-Steffensen加速法求方程在区间1,2内的一个实根，取绝对误差限为3. 实验目的

21、熟悉求解非线性方程的Aitken-Steffensen加速法4. 基础理论 将方程改写成等价形式，得到从初值开始的迭代公式后，基于迭代公式的Aitken-Steffensen加速法是通过“迭代-再迭代-加速”完成迭代的，具体过程为5. 实验环境 操作系统：Windows xp； 程序设计语言：Matlab6. 实验过程 为了验证Aitken-Steffensen加速法可以把一些不收敛的迭代加速成迭代收敛，我们使用将方程组变形为，取迭代函数，并利用宏定义出迭代函数由于不用保存迭代过程，所以用x0表示初值同时也存放前一步迭代的值，y和z是迭代过程中产生的yk和zk，x存放新迭代的结果算法的伪代码描

22、述：Input：初值x0；精度要求（即误差限）；迭代函数（x）；最大迭代次数NOutput：近似解或失败信息行号伪代码注释1n 1； 迭代次数计数器2while n N do 3 y （x0）； 迭代4 z （y）； 再迭代5 x x0 (y-x0)2 / (z - 2 y + x0) 加速6if | x x0 | < then 如果达到精度要求 7return x； 则输出近似值并退出程序 8 endif 9n n + 1；计数器加一10x0 x； 新近似值给x0做下次的初值11enddo12 return 错误信息 输出错误信息并结束程序7. 结果与分析 8. 附录：程序清单%*%*

23、 程 序 名：Aitken_Steffensen.m *%* 程序功能：用Aitken-Steffensen加速法求方程. *%*clc;clear all;phi=inline('0.5 * sqrt( 10 - x3)'); % 迭代函数x0=input('初值: x0 = ');epsilon=input('误差限: epsilon=');N=input('最大迭代次数: N=');disp(' n 迭代中间值y(n-1) 再迭代结构z(n-1) 加速后的近似值x(n)');fprintf('%2d%

24、54.6fn',0,x0); % 占2位整数显示步数0, 为了对齐, 占54位小数6位显示x0 1 % n是计数器 2 % while n<=N y= 3 ;% 迭代 z= 3 ;% 再迭代 x= 3 ;% 加速 % x0初值及前一步的近似值, y和z是中间变量, x是下一步的近似值 fprintf('%2d%18.6f%18.6f%18.6fn',n,y,z,x); % 显示中间值和迭代近似值 6 % 如果与上一步近似解差的绝对值不超过误差限 fprintf('nn 近似解 xx(%d)%fn',n,x); % 则输出近似根 (7)， 可简略为:

25、 fprintf('nn近似解 x=%f',x); return; % 并结束程序(7) 8 % 相当于endif 9 % 计数器加1 10 % 新近似值x作为下一次迭代的初值 11 fprintf('n 迭代%d次还不满足误差要求.nn',N); %输出错误信息(12)数值计算方法实验2报告班级：20xx级XXXXx班学号：20xx2409xxxx姓名：XXX成绩：1. 实验名称 实验2 非线性方程的迭代解法（之Newton下山法）2. 实验题目 用Newton下山法求方程在区间1,2内的一个实根，取绝对误差限为3. 实验目的 熟悉非线性方程的Newton下山

26、法4. 基础理论Newton下山法：Newton下山法公式为，使，其中5. 实验环境 操作系统：Windows xp； 程序设计语言：Matlab6. 实验过程定义函数f(x)和df(x)，其中df(x)是f(x)的导函数每步迭代时先取下山因子为1，尝试迭代，判断尝试结果是否满足下山因子，若满足则作为这步的迭代结果；否则将下山因子减半，然后再尝试为防止当前的xk是极小值点，附近不会有满足下述条件的其它点，使尝试陷入死循环，同时计算机中能表示出的浮点数也有下界，因此我们设置了最大尝试次数当超过最大尝试次数时，不再进行下山尝试由于反复尝试迭代且要判断下山条件，所以f(x0)和f(x0)会反复使用，

27、为避免重复计算浪费运行时间，将这两个值分别保存在变量fx0和dfx0而尝试产生的节点，判断下山条件时要用到它的函数值，若尝试成功，这个点会作为下一步的初值再使用，所以把该点的函数值也保存在变量fx中算法的伪代码描述：Input：初值x0；精度要求（即误差限）；函数及其导函数f（x）和f(x)；最大迭代次数N；K下山尝试最大次数Output：近似解或失败信息行号伪代码注释1n 1； 迭代次数计数器2F0 f(x0)； 3while n N do 4F0 f(x0)；5if F0= 0 then 若该点导数为06 return False； 则无法进行迭代，结束程序 7endif 8lambda

28、1；下山因子从1试起9k 1；尝试次数计数器10while k K do11x x0 lambda*F0 / F0Newton下山公式12 Fx f(x)；13if |Fx| < |F0| then判断下山条件14 break；退出尝试循环15 endif16lambda lambda/2；下山因子减半17k k+1；尝试次数计数器加118 endwhile19 if k > K then如果因超过尝试次数退出循环20 return False；则提示错误信息并结束程序21endif否则时尝试成功退出上边循环22if | x x0 | < then如果达到精度要求23retu

29、rn x；则输出近似值并退出程序24 endif 25 x0 x；新近似值给x0做下次的初值26F0 Fx；所求函数值下次也用到27n n + 1；计数器加一28enddo29return 错误信息 输出错误信息并结束程序7. 结果与分析 8. 附录：程序清单%*%* 程 序 名：NewtonDownhill.m *%* 程序功能：用Newton下山法求解非线性方程. * %*clc;clear all;f=inline('x3-x-1'); % 函数f(x)df=inline('3*x2-1'); % 函数f(x)的导函数x0=input('初值: x

30、0 = ');epsilon=input('误差限: epsilon=');N=input('最大迭代次数: N=');K=input('最大下山尝试次数: K='); 1 % 迭代次数计数器 2 % 存x0点函数值fprintf('nn n x(n) f(x(n)n');% 显示表头fprintf('%2d%14.6f%14.6fn',0,x0,fx0);% 2位整数显示0, 共14位小数6位显示x0和fx0 3 % while n N disp('');% 换行显示下山尝试过程的表头 d

31、isp(' 下山因子 尝试x(n) 对应f(x(n) 满足下山条件');disp(''); 4 % 存x0点导数值, 每次下山尝试不用重新计算if dfx0=0% 导数为0不能迭代disp(无法进行Newton迭代);return; end lambda=1.0; % 下山因子从1开始尝试k=1; % k下山尝试次数计数器while k<=K % 下山最多尝试K次 % 下山公式 fx=f(x);% 函数值 fprintf('%22.6f%14.6f%14.6f',lambda,x,fx);% 显示尝试结果if (abs(fx)<abs

32、(fx0)% 判断是否满足下山条件 fprintf(' 满足n'); break; % 是, 则退出下山尝试的循环 else fprintf(' 不满足n'); end lambda=lambda/2;% 不是, 则下山因子减半 k=k+1;% 计数器加1 end if k>K fprintf('n 下山条件无法满足, 迭代失败.nn'); return; endfprintf('%2d%14.6f%14.6fn',n,x,fx);% 2位整数显示步数n, 共14位小数6位显示下步迭代结果 22 % 达到精度要求否 fpri

33、ntf('nn 方程的近似解为: x%fnn',x);% (23) return;% 达到, 则显示结果并结束程序(23) end% (24) % 用x0,fx0存放前一步的近似值和它的函数值, 进行循环迭代 25 26 27 28 fprintf('n 迭代%d次还不满足误差要求.nn',N);数值计算方法实验2报告班级：20xx级XXXXx班学号：20xx2409xxxx姓名：XXX成绩：1. 实验名称 实验2 非线性方程的迭代解法（之弦截法）2. 实验题目 用弦截法求方程在区间1,2内的一个实根，取绝对误差限为3. 实验目的 熟悉非线性方程的弦截法4. 基

34、础理论将Newton迭代法中的导数用差商代替，得到弦截法(或叫正割法)公式5. 实验环境 操作系统：Windows xp； 程序设计语言：Matlab6. 实验过程不保存迭代过程，所以始终以x0和x1分别存放xk-1和xk，而x存放新产生的迭代值xk+1，这样，下一次迭代时需要把上一步的x1（即xk）赋值于x0（做新的xk-1）这些点的函数值会重复用到，在迭代公式中也要用到，上一步的x1作为下一步的x0也会再一次用它的函数值，为减少重新计算该点函数值的运行时间，将x1点的函数值保存在变量fx1中算法的伪代码描述：Input：初值x0，x1；精度要求（即误差限）；函数f（x）；最大迭代次数NOu

35、tput：近似解或失败信息行号伪代码注释1fx0 f(x0)； 初值点的函数值2fx1 f(x1)；3n 2； 迭代次数计数器(从2开始)4while n N do 5 x x1(x1-x0)*f(x1) / (f(x1 ) f(x0)； 迭代6if | x x0 | < then 如果达到精度要求 7return x； 则输出近似值并退出程序 8 endif 9x0 x1； （按顺序）最新的2个近似值，分别做下次的初值x0，x110x1 x；11 fx0 fx1； 下次初值点的函数值12fx1 f(x1)；13 n n+1；14enddo15return 错误信息 输出错误信息并结束程

36、序7. 结果与分析 8. 附录：程序清单%*%* 程 序 名：SecantMethod.m *%* 程序功能：用弦截法求解非线性方程. *%*clc;clear all;f=inline('2*x3-5*x-1'); % 函数f(x)x0=input('第一初值: x0 = ');x1=input('第二初值: x1 = ');epsilon=input('误差限: epsilon=');N=input('最大迭代次数: N=');fprintf('n n x(n)n'); % 显示表头fprint

37、f('%2d%14.6fn', 0, x0); % 占2位显示步数0, 共14位其中小数6位显示x0fprintf('%2d%14.6fn', 1, x1); % 占2位显示步数1, 共14位其中小数6位显示x1 1 % 存x0点函数值 2 % 存x1点函数值 3 % 迭代计数器 4 % while n N % 弦截法公式 fprintf('%2d%14.6fn', n, x); % 显示迭代过程 6 % 达到精度要求否 fprintf('nn 方程的近似解为: x%fnn', x); return; % 达到, 则显示结果并结束程序 8 9 % 原x1做x0为前两步的近似值 10 % 现x做x1为一两步的近似值 11 % x0点函数值 12 % 计算x1点函数值, 为下一次循环 13 % 计数器加1 14 fprintf('n 迭代%d次还不满足误差要求.nn',N); 数值计算方法实验3报告班级：20xx级XXXXx班学号：20xx2409xxxx姓名：XXX成绩：1. 实验名称 实验3 解线性方程组的直接法(之Gauss消去法)2. 实验题目 用Gauss消去法求解线性方程组3. 实验目的 掌握解线性方程组的Gauss消去法4. 基