第一章 误差与范数
第一章 误差与范数
1.1 误差的来源例1.1.1 用差商ha f h a f a f )()()(-+≈'求x x f ln )(=在3=x 处导数的近似值.取1.0=h ，1000.0=h ，h =0.000 000 000 000 001和h =0.000 000 000 000 000 1分别用MATLAB 软件计算，取十五位数字计算.解 在MATLAB 工作窗口输入下面程序>>a=3;h=0.1;y=log(a+h)-log(a);yx=y/h运行后得yx = 0.32789822822991 将此程序中h 改为0.000 1，运行后得yx = 0.33332777790385后者比前者好.再取h = 0.000 000 000 000 001，运行后得yx = 0.44408920985006不如前者好.取h = 0.000 000 000 000 000 1，运行后得yx = 0算出的结果反而毫无价值.例1.1.2 分别求方程组b AX =在下列情况时的解，其中A ⎪⎪⎭⎫⎝⎛=011111.. （1）⎪⎪⎭⎫ ⎝⎛=22b ;（2）⎪⎪⎭⎫⎝⎛=0122.b . 解 （1） 首先将方程组b AX =化为同解方程b A X 1-=，然后在MATLAB 工作窗口输入程序>> b=[2,2]';A=[1,1;1,1.01]; X=A\b运行后输出当⎪⎪⎭⎫ ⎝⎛=22b 时，b AX =的解为⎪⎪⎭⎫⎝⎛=02X ；（2）同理可得，当⎪⎪⎭⎫ ⎝⎛=0122.b 时，b AX =的解为⎪⎪⎭⎫⎝⎛=11X .例1.1.3 计算e 的近似值. 解 泰勒级数e +++++++=!!4!3!21432 n x x x x x n x)(∞<<-∞x , 取1=x ，得e +++++++=!1!41!31!2111 n . （1.2）这是一个无限过程，计算机无法求到精确值.只能在（1.2）取有限项时计算，再估计误差.如果取有限项!!!!)( n s n 1413121111++++++=作为e 的值必然会有误差，根据泰勒余项定理可知其截断误差为e !)1()1( +=-n e s n θ)10(<<θ.如果取（1.2）的前九项，输入程序>> n =8;s=1;S =1; for k=1:ns=s*k;S=S+1/s, ends,S,R=3/(s*(n+1)) 或>>S1=1+1+1/2+1/(1*2*3)+1/(1*2*3*4)+1/(1*2*3*4*5)+1/(1*2*3*4*5*6)+1/(1*2*3*4*5*6*7)+1/(1*2*3*4*5*6*7*8),R1=3/(1*2*3*4*5*6*7*8*9)运行后结果S = R =2.71827876984127 8.267195767195768e-006 因为截断误差为e ),10(101968.267!93!)18()1(6-8<<⨯≈<+=-θθ e s 所以e 的近似值e ≈≈++++++++=!81!71!61!51!41!31!2111)1(8 s 2.718 28.1.2 误差和有效数字例1.2.1 取282.718作为e 的四舍五入近似值时，求其绝对误差和相对误差. 解 在MATLAB 工作窗口输入程序>>juewu=exp(1)-2.71828运行后输出结果为juewu = 1.828 459 045 505 326e-006例1.2.2 计算⎰π20sin x xd x 的近似值，并确定其绝对误差和相对误差.解 因为被积函数xxsin 的原函数不是初等函数，故用泰勒级数求之.++-+-=!!!!sin ９ ７５ ３８６x x x x x x 421 )(∞<<-∞x , （1.5） 这是一个无限过程，计算机无法求到精确值.可用（1.5）的前四项!!!７５ ３６x x x -+-421代替被积函数xxsin ，得 ⎰π=20sin x x y d ⎰π≈20(x !!!142７５ ３６x x x -+-)d x =!7)2(!5)2(!3)2(2753７５ ３⋅π-⋅π+⋅π-π=y ˆ. 根据泰勒余项定理和交错级数收敛性的判别定理，得到绝对误差!99)2(ˆ9⋅<-=πyy R = WU ， 在MA TLAB 命令窗口输入计算程序如下：>>syms xf=1-x^2/(1*2*3)+x^4/(1*2*3*4*5)-x^6/(1*2*3*4*5*6*7)y=int(f,x,0,pi/2),y1=double(y)y11=pi/2-(pi/2)^3/(3*3*2)+(pi/2)^5/(5*5*4*3*2)-(pi/2)^7/(7*7*6*5*4*3*2)inf=int(sin(x)/x,x,0,pi/2) ,infd=double(inf) WU =(pi/2)^9/(9*9*8*7*6*5*4*3*2), R =infd-y11因为运行后输出结果为： =y 1.370 762 168 154 49，yˆ=1.370 744 664 189 38，=R 1.750 396 510 491 47e-005， WU = 1.782 679 830 970 664e-005410-<.所以，yˆ的绝对误差为=ε410-，故⎰π=20sin x xy d 7 1.370≈x .yˆ的相对误差为 =r ε71.37010ˆ4-=y ε<0.007 3%.1.3 误差估计的基本方法例1.3.4 设计三种算法求方程01522=-+x x 在)3,2(的一个正根*x 的近似值，并研究每种算法的误差传播情况.解 为解已知方程，我们可以设计如下三种算法，然后将计算结果与此方程的精确解5.2*=x 比较，观察误差的传播.算法1 将已知方程化为同解方程=x 2215x -.取初值20=x ，按迭代公式21215k k x x -=+依次计算 ,,,,21n x x x ，结果列入表1–3中.算法2 将已知方程化为同解方程1215+=x x .取初值20=x ，按迭代公式 12151+=+k k x x依次计算 ,,,,21n x x x ，结果列入表1–3中.算法3 将已知方程化为同解方程141522+-+-=x x x x x .取初值20=x ，按迭代公式为1415221+-+-=+k k kk k x x x x x 依次计算 ,,,,21n x x x ，结果列入表1–3中.我们为这三种算法的计算编写两套MATLAB 程序如下： （1）MATLAB 主程序function [k,juecha,xiangcha,xk]= liti112(x0,x1,limax) % 输入的量--x0是初值, limax 是迭代次数和精确值x; % 输出的量--每次迭代次数k 和迭代值xk,% --每次迭代的绝对误差juecha 和相对误差xiangcha ， x(1)=x0;for i=1:limaxx(i+1)=fl(x(i));%程序中调用的fl.m juecha = abs(x(i)-x1);xiangcha = juecha /( abs(x(i))+eps);xk=x(i);k=i-1;[(i-1),juecha,xiangcha,xk] end（2）MATLAB 调用函数程序及其计算结果①算法2的MATLAB 调用函数程序function y1=fl(x)y1=15/(2*x+1);② 将MATLAB 主程序和调用函数程序分别命名liti112.m 和fl.m ，分别保存为M 文件，然后在MATLAB 工作窗口输入命令>> [k,juecha,xiangcha,xk]= liti112(2,2.5,100) ③运行后输出计算结果列入表1–3和表 1-4中.④将算法2的MATLAB 调用函数程序的函数分别用y1=15-2*x^2和y1=x-(2*x^2+x-15)/(4*x+1)代替，得到算法1和算法3的调用函数程序，将其保存，运行后将三种算法的前8个迭代值821,,,x x x 列在一起（见表 1-3），进行比较.将三种算法的821,,,x x x 对应的绝对误差和相对误差的值列在一起（见表 1-4），进行比较.1.4 数值计算中应注意的问题例1.4.1 求数）181(71915-+⨯=-x 的近似值. 解 （1）直接用MATLAB 命令>> x=(7^15)*(sqrt(1+8^(-19))-1)运行后输出结果x = 0问题出现在两个相近的数1981-+与1相减时，计算机运行程序>>sqrt(1+8^(-19))-1运行后输出结果ans = 0由于计算机硬件只支持有限位机器数的运算，因此在计算中可能引入和传播舍入误差.因为有效数字的严重损失，导致输出18119-+-的结果为0，计算机不能再与数157继续进行真实的计算，所以，最后输出的结果与x 的精确值不符.（2）如果化为18187)181(71919151915++⨯=-+⨯=---x ，再用MATLAB 命令>> x=(7^15)*( (8^(-19))/(sqrt(1+8^(-19))+1))运行后输出结果x = 1.6471e-005 这是因为18119-+-化为18181919++--后，计算机运行程序>> x= (8^(-19))/(sqrt(1+8^(-19))+1)运行后的结果为x =3.4694e-018 由于有效数字的损失甚少，所以运算的结果-18103.4694⨯再与157继续计算，最后输出的结果与x 的精确值相差无几.例1.4.2 求数)13030ln(2--=y 的近似值. 解 （1）直接用MATLAB 程序>> x=30;x1= sqrt(x^2-1)运行后输出结果x1 = 29.9833 输入MATLAB 程序>> x=30; x1=29.9833;y=log(x-x1)运行后输出结果y = -4.0923（2）因为)13030ln(2--中的30=x 很大，如果采用倒数变换法111221-+=--=x x x x z ，即130301ln)13030ln(22-+=--)190030ln(-+-=.输入MATLAB 程序>> x=30;y=-log(x+sqrt(x^2-1))运行后输出结果y = -4.0941（3）输入MA TLAB 程序>> x=30; y=log(x-sqrt(x^2-1))运行后输出结果y = -4.0941 可见，（2）计算的近似值比（1）的误差小.参加计算的数，有时数量级相差很大.如果不注意采取相应的措施，在它们的加减法运算中，绝对值很小的那个数经常会被绝对值较大的那个数“吃掉”，不能发挥其作用，造成计算结果失真.例1.4.4 请在16位十进制数值精度计算机上利用软件MATLAB 计算下面的两个数0.30.1111111111111111*++=x 和0.30.11111111111111111*++=y将计算结果与准确值比较，解释计算结果.解 在MATLAB 工作窗口输入下面程序>> x=111111*********+0.1+0.3, y=1111111111111111+0.1+0.3运行后输出结果x = 1.111111*********e+014，y =1.111111*********e+015 从输出的结果可以看出，x *x =，而y *y ≠.为什么*y 仅仅比*x 多一位1，而y *y ≠呢？这是因为计算机进行运算时，首先要把参加运算的数写成绝对值小于1而“阶码”相同的数，这一过程称为数的“对阶”.在16位十进制数值精度计算机上利用软件MATLAB 计算这两个数，把运算的数*x 写成浮点规格化形式为，151515*103000**********.0001010000000000000.000100111111111111111.0⨯+⨯+⨯=x在16位十进制数值精度计算机上，三项的数都表示为小数点后面16位数字的数与1510之积，所以，计算机没有对数进行截断，而是按原来的三个数进行计算.因此，计算的结果x *x =.而161616*10030000000000000.00010010000000000000.000101111111111111111.0⨯+⨯+⨯=y三项的数都表示写成绝对值小于1而“阶码”都为1610的数以后，第一项的纯小数的小数点后面有16位数字.但是，后两项数的纯小数的小数点后面有17位数字，超过了16位十进制数值精度计算机的存储量，计算机对后两项的数都进行截断最后一位，即后两项的数都是16位机上的零，再进行计算，所以计算结果与实际不符.五、向量和矩阵的范数例1．5.1: 用matlab 求下列向量的2,1,5，,-∞∞范数。
第一章 误差与范数
习 题1.11. 分别求方程组b AX=在下列情况时的解，其中⎪⎪⎭⎫⎝⎛=01.2222A .（1）⎪⎪⎭⎫ ⎝⎛=22b ；（2）⎪⎪⎭⎫ ⎝⎛=01.22b . 2. 用差商h a f h a f a f )()()('-+≈，求x x f ln )(=在5=x 处导数)5('f 的近似值，取十五位数字计算，h 分别取0.000 1和0.000 000 01.习 题1.21. 为了使 4.420=的近似值的相对误差小于%1.0，试问取几位有效数字？2. 写出光速10*10)001000.0925997.2(⨯±=c s /cm 的绝对误差和相对误差. 3. 取655923.141作为π的四舍五入近似值时，求其有效数字、绝对误差和相对误差.4. 问722,141.3,142.3分别作为π的近似值各具有几位有效数字？习 题1.31.测得电压)2110(±=V V ，电流)5.020(±=I A ，则由欧姆定律得Ω==5.5I VR ，求R 的绝对误差和相对误差.2. 设nx y =，求y 的相对误差与x 的相对误差的关系式.3. 设计三种算法求方程0322=-+x x 在)2,1(的一个正根*x 的近似值，并研究每种算法的误差传播情况.4. 已测得某场地长*x 的值110=x m ，宽*y 的值80=y m ，已知x 的绝对误差是2.0m ，y 的绝对误差是1.0m ，求场地面积xy S =的绝对误差和相对误差.5. 正方形的边长约为200cm ，怎样测量才能使其面积误差不超过2cm 1？1. 求数）274(41915-+⨯=-x 的近似值.2. 利用四位数学用表求2cos 1-的近似值，要求至少有2位有效数字.3. 请在16位十进制数值精度计算机上利用软件MATLAB 计算下面的两个数0.30.1333333333333333*++=x 和0.30.133333333333333333*++=y将计算结果与准确值比较，解释计算结果.4. 利用恒等变换，使下列表达式的计算结果比较精确.（1）1,11>>--+x x x x x ； （2）e12-x ； （3）x xsin cos 1-. 5. 设x y ln =.当)0(>≈a a x 时，如果已知对数a ln 的绝对误差限为n-⨯1021，试估计真数a 的相对误差限及有效数字位数.6. 指出下列各题的合理计算途径（对给出具体数据的，请算出结果）：（1）o1cos 1-（三角函数值取四位有效数字）；（2）)13030ln(2--；（3）x xsin cos 1- （其中x充分小）； （4）127x .习 题 1.51.设⎪⎪⎪⎪⎪⎭⎫⎝⎛=1095791068565778710A ,当b T)31,33,23,32(=有微小误差)1.0,1.0,1.0,1.0(--=δb T 时，估计方程组b Ax =解的变化.2. 用MATLAB 求矩阵A 的2范数条件数、1范数条件数、∞范数条件数、Frobenius 范数条件数、条件数倒数的估计值,其中3. 用MATLAB 求第2题中矩阵A 的2范数、1范数、∞范数、Frobenius 范数、谱半径.4. 用MATLAB 求下列向量的2范数、1范数、∞范数、5范数、-∞范数:(1))5,2,1,1(-=X ;(2) )12,45,06.0,89,4.0(--=Y T. 5. 线性方程组求解与性态讨论.求b Ax =的解向量x ，其中然后把b 扰动为b ˆT )99.30,01.33,99.22,01.32(=，再求解b xA ˆˆ=. 计算xx x -ˆ（使用1范数或∞范数），讨论方程组性态..98755.03.02.064321543210⎪⎪⎪⎪⎪⎭⎫⎝⎛---=A .31332332,9105791068565778710⎪⎪⎪⎪⎪⎭⎫ ⎝⎛=⎪⎪⎪⎪⎪⎭⎫⎝⎛=b A。
误差与范数
第一章误差与范数1.1 误差的来源1.2 误差和有效数字1.3 误差估计的基本方法1.4 数值计算中应注意的问题第二章非线性方程（组）的数值解法2.1 方程（组）的根及其MATLAB命令2.1.2 求解方程（组）的solve命令2.1.4 求解方程（组）的fsolve命令2.2 搜索根的方法及其MATLAB程序2.2.1 作图法及其MATLAB程序2.2.2 逐步搜索法及其MATLAB程序2.3 二分法及其MATLAB程序2.3.1 二分法求方程根的迭代次数的MATLAB命令2.3.2 二分法的MATLAB程序2.4 迭代法及其MATLAB程序2.4.2 迭代法的MATLAB程序12.4.5 迭代法的MATLAB程序22.5 迭代过程的加速方法及其MATLAB程序2.5.2 加权迭代的MATLAB程序2.5.4 艾特肯(Aitken)加速方法的MATLAB程序2.6 牛顿(Newton)切线法及其MATLAB程序2.6.2 牛顿切线法的收敛性及其MATLAB程序2.6.3 牛顿切线法的MATLAB程序2.6.6 牛顿切线法的加速及其两种MATLAB程序2.7 割线法及其MATLAB程序2.7.2 割线法的MATLAB程序2.8 抛物线法及其MATLAB程序2.8.2 抛物线法的MATLAB程序2.9 求解非线性方程组的牛顿法及其MATLAB程序2.9.1 求解二阶非线性方程组的牛顿法及其MATLAB程序2.9.2 求解n阶非线性方程组的牛顿法及其MATLAB程序2.9.3 求解n阶非线性方程组的拟牛顿法及其MATLAB程序第三章解线性方程组的直接方法3.1 方程组的逆矩阵解法及其MATLAB程序3.1.3 线性方程组有解的判定条件及其MATLAB程序3.2 三角形方程组的解法及其MATLAB程序3.2.2 解三角形方程组的MATLAB程序3.3 高斯消元法和列主元素消元法及其MATLAB程序3.3.1 高斯消元法及其MATLAB程序3.3.2 列主元素消元法及其MATLAB程序3.4 LU分解法及其MATLAB程序3.4.1 判断矩阵LU分解的充要条件及其MATLAB程序3.4.2 直接LU分解法及其MATLAB程序3.4.3 含交换阵的LU分解法及其MATLAB程序3.4.4 正定对称矩阵的Cholesky分解及其MATLAB程序3.4.5 正定对称矩阵的Cholesky分解及其MATLAB程序3.5 求解线性方程组的LU方法及其MATLAB实现3.5.1 求解线性方程组的Cholesky分解及其MATLAB实现3.5.2 求解线性方程组的直接LU分解法及其MATLAB实现3.5.3 解线性方程组的选主元的LU方法及其MATLAB实现3.6 误差分析及其两种MATLAB程序3.6.1 用MATLAB软件作误差分析AX=解的性态的MATLAB程序3.6.2 求P条件数和讨论bAX=解和A的性态的MATLAB程序3.6.3 用P范数讨论b第四章解线性方程组的迭代法4.1 迭代法和敛散性及其MATLAB程序4.1.2 用谱半径判别迭代法产生的迭代序列的敛散性的MATLAB程序4.2 雅可比（Jacobi）迭代及其MATLAB程序4.2.2 判别雅可比迭代收敛性的MATLAB程序4.2.3 雅可比迭代的两种MATLAB程序4.3 高斯—塞德尔（Gauss-Seidel）迭代及其MATLAB程序4.3.3 高斯—塞德尔迭代两种MATLAB程序4.4 解方程组的超松弛迭代法及其MATLAB程序4.4.2 超松弛迭代法收敛性及其MATLAB程序4.4.3 超松弛迭代法的MATLAB程序第五章矩阵的特征值与特征向量的计算5.2 幂法及其MATLAB程序5.2.2 幂法的MATLAB程序5.3 反幂法和位移反幂法及其MATLAB程序5.3.3 原点位移反幂法的两种MATLAB程序5.4 雅可比(Jacobi)方法及其MATLAB程序5.4.3 雅可比方法的MATLAB程序5.5 豪斯荷尔德方法及其MATLAB程序5.5.1 豪斯荷尔德方法及其MATLAB程序5.5.2 矩阵约化为上Householder阵及其MATLAB程序5.5.3 实对称矩阵的三对角化及其MATLAB 程序5.6 QR 方法及其MATLAB 程序5.6.2 QR 方法的MATLANB 程序5.6.5 最末元素位移QR 法计算实对称矩阵特征值及其MATLAB 程序5.6.6 求根位移QR 方法计算实对称矩阵A 的特征值及其MATLAB 程序5.7 广义特征值问题及其MATLAB 程序5.7.2 用MATLAB 计算BX AX λ=型的广义特征值和特征向量第六章 函数的插值方法6.1 插值问题及其误差6.1.2 与插值有关的MATLAB 函数6.2 拉格朗日（Lagrange ）插值及其MATLAB 程序6.2.1 线性插值及其MATLAB 程序6.2.2 抛物线插值及其MATLAB 程序 6.2.3 n 次拉格朗日（Lagrange ）插值及其MATLAB 程序6.2.5 拉格朗日多项式和基函数的MATLAB 程序6.2.6 拉格朗日插值及其误差估计的MATLAB 程序6.3 牛顿（Newton ）插值及其MATLAB 程序6.3.3 牛顿插值多项式、差商和误差公式的MATLAB 程序6.3.4 牛顿插值及其误差估计的MATLAB 程序6.3.5 牛顿插值法的MATLAB 综合程序6.4 埃尔米特（Hermite ）插值及其MATLAB 程序6.4.3 埃尔米特插值多项式的MATLAB 程序6.5 分段插值及其MATLAB 程序6.5.1 高次插值的振荡和MATLAB 程序6.5.3 分段线性插值的MATLAB 程序6.5.4 作有关分段线性插值图形的MATLAB 程序6.5.5 用MATLAB 计算有关分段线性插值的误差6.6 分段埃尔米特（Hermite ）插值及其MATLAB 程序6.6.2 Hermite 插值的MATLAB 程序6.6.3 作有关分段Hermite 插值图形的MATLAB 程序6.6.4 用MATLAB 计算有关分段Hermite 插值的误差6.7 三次样条插值及其MATLAB 程序6.7.4 用一阶导数计算的几种样条函数和MATLAB 程序6.7.6 用MATLAB 计算三次样条6.7.7 几种作三次样条有关图像的MATLAB 程序6.7.8 用MATLAB 计算有关分段三次样条的误差6.8 高元插值及其MATLAB 程序6.8.1 MESHGRID 命令的功能和调用格式6.8.2 单调数据点上的二元插值及其MATLAB 程序6.8.3 三元插值及其MATLAB 程序第七章 函数逼近与曲线（面）拟合7.1 曲线拟合、误差及其MATLAB 实现7.2 曲线拟合的线性最小二乘法及其MATLAB 程序7.3 函数)(x r k 的选取及其MATLAB 程序7.4 多项式拟合及其MATLAB 程序7.5 拟合曲线的线性变换及其MATLAB 程序7.6 函数逼近及其MATLAB 程序7.7 三角多项式逼近及其MATLAB 程序7.8 随机数据点上的二元拟合及其MATLAB 程序7.9 随机数据点上的n 元拟合及其MATLAB 程第八章 数值微分8.2 一阶导数的数值计算及其MATLAB 程序 (117)8.2.1 差商求导及其MATLAB 程序8.2.2 中心差商公式求导四种MATLAB 程序8.2.3 Richardson 外推法求导及其MATLAB 程序8.2.4 牛顿多项式求导及其MATLAB 程序8.2.5 diff 和gradient 函数在数值求导中的应用8.3 高阶导数的数值计算及其MATLAB 程序8.3.1 插值或拟合高阶数值导数及其MATLAB 程序8.3.2 高阶泰勒数值导数及其MATLAB 程序8.4 数值梯度和数值偏导数的计算及其MATLAB 程序8.4.1 梯度和偏导数的数值计算及其MATLAB 程序8.4.2 计算雅克比矩阵及其行列式的MATLAB 方法第九章 数值积分9.1 积分的符号计算9.1.1 定积分的MATLAB 符号计算9.1.2 变限积分的MATLAB 符号计算9.2 数值积分的思想及其MATLAB 程序9.2.3 矩形公式的MATLAB 运算9.3 插值型数值积分及其MATLAB 程序9.3.2 梯形公式的四种MATLAB 程序9.3.4 辛普生（Simpson ）数值积分的MATLAB 程序9.3.6 牛顿-科茨（Newton-Cotes ）数值积分和误差分析的三种MATLAB 程序9.3.7 利用三次样条求表格型数值积分的MATLAB 方法9.3.8 利用拉格朗日插值等方法求表格型数值积分的MATLAB 方法9.4 龙贝格公式及其MATLAB 程序9.4.2 龙贝格积分的MATLAB 程序9.5 自适应积分及其MATLAB 程序9.6 高斯（Gauss ）型积分公式及其MATLAB 程序9.6.2 在]1,1[ 上的高斯-- 勒让德积分公式及其MATLAB 程序9.6.3 在],[b a 上的高斯-- 勒让德积分公式及其MATLAB 程序9.6.4 Radau 和Lobatto 积分公式及其MATLAB 程序9.7 反常积分的计算及其MATLAB 程序9.7.1 无穷积分的符号计算及其MATLAB 程序9.7.2 无穷积分的近似计算及其MATLAB程序9.7.3 无界函数反常积分的符号计算的五种MATLAB程序9.7.4 无界函数反常积分的近似计算的两种MATLAB程序9.8 多重积分的计算及其MATLAB程序9.8.1 二重积分的符号计算及其MATLAB程序9.8.2 二重积分的梯形公式及其MATLAB程序9.8.3 矩形域上的辛普生公式及其MATLAB程序9.8.4 一般域上二重积分的数值计算及其MATLAB程序9.8.5 三重积分的计算及其MATLAB程序第十章常微分方程（组）求解10.1 常微分方程（组）的MATLAB符号求解10.1.1 MATLAB求常微分方程（组）的通解10.1.2 MATLAB求常微分方程（组）的特解10.1.3 线性常微分方程组的MATLAB解法10.3 欧拉方法及其MATLAB程序10.3.2 向前欧拉方法的三种MATLAB程序10.3.4 向后欧拉方法的MATLAB程序10.4 改进的欧拉方法及其MATLAB程序10.4.2 梯形公式的两个MATLAB程序10.4.4 改进的欧拉公式的MATLAB程序10.5 龙格—库塔方法10.5.2 二阶龙格—库塔方法及其MATLAB程序10.5.3 三阶龙格—库塔方法及其MATLAB程序10.5.4 四阶龙格—库塔方法及其MATLAB程序10.5.5 自适应龙格—库塔方法及其MATLAB库函数10.6 线性多步法及其MATLAB程序10.6.2 亚当斯显式公式及其MATLAB程序10.6.3 亚当斯隐式公式及其MATLAB程序10.6.4 米尔恩公式及其MATLAB程序10.6.5 汉明公式及其MATLAB程序10.6.6 预测-校正系统及其MATLAB程序10.7 一阶（高）阶微分方程（组）的数值解及其MATLAB程序10.7.1 一阶微分方程组的数值解及其MATLAB程序10.7.2 高阶微分方程（组）的数值解及其MATLAB实现10.8 边值问题的数值解及其MATLAB程序10.8.1 打靶法及其MATLAB程序10.8.2 有限差分方法及其MATLAB程序10.8.3 求解常微分方程（组）边值问题数值解的MATLAB库函数。
误差理论第一章绪论
①在近似加减运算中，各运算数据以小数位数最少的数据为 准，其余各数据可多取一位小数，但最后结果应与小数位数 最少的数据位数相同。
②在乘除运算中，各运算数据以有效位数最少的数据位数为 准，其余各数据要比有效位数最少的数据位数多取一位数字， 而最后结果应与有效位数最少的数据位数相同。
③平方或开方运算，可按乘除运算处理。 12
三、误差分类 按误差的特点和性质，误差可分为系统误差、差：在同一条件下，多次测量同一量值时， 误差的绝对值符号保持不变，或在条件改变时，按一定规 律变化的误差。（如常用的杆秤）
①按对误差掌握的程度分：已定系统误差（误差的绝对值 和符号已确定）；未定系统误差（误差的绝对值和符号未 能确定，但范围可估计出）。
15000 15080.3 80.3 0.4%
20000
20000
二、误差来源
（一）测量装置误差：①标准量具误差（如标准电阻、标准
砝码）； ②仪器误差（如天平、压力表、温度计）；③附件
误差（如测长仪的标准环规）
5
（二）环境误差：各种环境因素与规定的标准状态不一致 而引起的误差（如温度、湿度、振动等） （三）方法误差：由测量方法不完善而引起的。（如间接 测量圆直径） （四）人员误差：由测量人员的习惯或疲劳原因等引起的 误差。
§1-1 研究误差的意义
误差存在的必然性和普遍性：由于实验方法、实验设备的不 完善、周围环境的影响、人们认识能力的限制，使得测量和 实验所得的数据和被测量的真值之间，不可避免存在差异。 尽管科技发展和人们认识水平的提高可使误差控制的很小， 但终究不能完全消除，这种必然性和普遍性已为大量实践所 证明。
（三）粗大误差：超出在规定条件下预期的误差，又称 “寄生误差”，此误差值较大，明显歪曲测量结果（如 人员因素、有缺陷的仪器等）
第1章引言、误差(Y)
二、误差(error)和误差限
1.绝对误差 设 x为准确值, x* 为 x的一个近似值, 称
e( x* ) = x* − x
为近似值 x*的 绝对误差, 简称误差,可简记为 e.
16
Wuhan University
二、误差和误差限
2.绝对误差限 因为准确值 x 往往是未知甚至是无法 知道的 * * 因此 e(x ) = x − x 往往也无法求出
x* 具有 n 位有效数字 则便说 近似值
23
Wuhan University
三、有效数字
例3.
求下列四舍五入近似值的有效数字个数.
x * = 3587 . 64 = 0 . 358764 × 10 4
若
1 1 −2 x − x ≤ ×10 = ×10 4−6 2 2
*
x* = 3587.64具有6位有效数字. 则
π = 3.141 59 0.000 002 65 L≤ 0.5 × 10
* *
−5
6位有效数字
0.000 000 04 L≤ 0.5 × 10 −7 8位有效数字 π = 3.141 592 7
π = 3.1415
*
≤ 0.5 × 10 −3 4位有效数字 0.000092 L
22
Wuhan University
12
Wuhan University
第一章 误差
一、误差的种类及来源
1.模型误差 在建立数学模型过程中,要将复杂的现 象抽象归结为数学模型,往往要忽略一 些次要因素的影响,而对问题作一些简 化,因此和实际问题有一定的区别. 2.观测误差 在建模和具体运算过程中所用的数据往 往是通过观察和测量得到的,由于精度的 限制,这些数据一般是近似的,即有误差
张宏伟《矩阵与数值分析》第一章 绪论
第一章 误差分析与向量与矩阵的范数一、内容提要本章要求掌握绝对误差、相对误差、有效数字、误差限的定义及其相互关系；掌握数值稳定性的概念、设计函数计算时的一些基本原则和误差分析；熟练掌握向量和矩阵范数的定义及其性质。
1．误差的基本概念和有效数字 1）．绝对误差和相对误差的基本概念设实数x 为某个精确值，a 为它的一个近似值，则称a x -为近似值a 的绝对误差，简称为误差． 当0≠x 时，xax -称为a 的相对误差．在实际运算中，精确值x 往往是未知的，所以常把aax -作为a 的相对误差． 2）．绝对误差界和相对误差界的基本概念设实数x 为某个精确值，a 为它的一个近似值，如果有常数a e ，使得 a e a x ≤- 称a e 为a 的绝对误差界，或简称为误差界．称ae a 是a 的相对误差界．此例计算中不难发现，绝对误差界和相对误差界并不是唯一的，但是它们越小，说明a 近似x 的程度越好，即a 的精度越好．3）．有效数字设实数x 为某个精确值，a 为它的一个近似值，写成n k a a a a 21.010⨯±=它可以是有限或无限小数的形式，其中),2,1( =i a i 是9,,1,0 中的一个数字，k a ,01≠为整数．如果n k a x -⨯≤-1021则称a 为x 的具有n 位有效数字的近似值．如果a 有n 位有效数字，则a 的相对误差界满足：n a a a x -⨯≤-111021。
4）．函数计算的误差估计如果),,,(21n x x x f y =为n 元函数，自变量n x x x ,,,21 的近似值分别为n a a a ,,,21 ，则)(),,,(),,,(12121k k n k akn n a x x fa a a f x x x f -⎪⎪⎭⎫⎝⎛∂∂≈-∑= 其中),,,(21n kak a a a f x x f ∂∂=⎪⎪⎭⎫ ⎝⎛∂∂，所以可以估计到函数值的误差界，近似地有 k a n k aka n n e x fe a a af x x x f ∑=⎪⎪⎭⎫⎝⎛∂∂≈≤-12121),,,(),,,( 如果令2=n ，设21,x x 的近似值分别为21,a a ，其误差界为111a e a x ≤-和≤-22a x 2a e ，取),(21x x f y =为21,x x 之间的四则运算，则它们的误差估计为，1121a a a a e e e +≈±；112121a a a a e a e a e +≈⋅；22211121a e a e a e a a a a +≈，02≠a 。
误差与范数概念附习题与答案
误差的来源例1.1.1 用差商ha f h a f a f )()()(-+≈'求x x f ln )(=在3=x 处导数的近似值.取1.0=h ，1000.0=h ，h =0.000 000 000 000 001和h =0.000 000 000 000 000 1分别用MATLAB 软件计算，取十五位数字计算.解 在MATLAB 工作窗口输入下面程序>>a=3;h=0.1;y=log(a+h)-log(a);yx=y/h运行后得yx = 0.32789822822991 将此程序中h 改为0.000 1，运行后得yx = 0.33332777790385后者比前者好.再取h = 0.000 000 000 000 001，运行后得yx = 0.44408920985006不如前者好.取h = 0.000 000 000 000 000 1，运行后得yx = 0算出的结果反而毫无价值.例1.1.2 分别求方程组b AX =在下列情况时的解，其中A ⎪⎪⎭⎫⎝⎛=011111.. （1）⎪⎪⎭⎫ ⎝⎛=22b ;（2）⎪⎪⎭⎫⎝⎛=0122.b . 解 （1） 首先将方程组b AX =化为同解方程b A X 1-=，然后在MATLAB 工作窗口输入程序>> b=[2,2]';A=[1,1;1,1.01]; X=A\b运行后输出当⎪⎪⎭⎫ ⎝⎛=22b 时，b AX =的解为⎪⎪⎭⎫⎝⎛=02X ；（2）同理可得，当⎪⎪⎭⎫ ⎝⎛=0122.b 时，b AX =的解为⎪⎪⎭⎫⎝⎛=11X .例1.1.3 计算e 的近似值. 解 泰勒级数e +++++++=!!4!3!21432 n x x x x x n x)(∞<<-∞x , 取1=x ，得e +++++++=!1!41!31!2111 n . （1.2）这是一个无限过程，计算机无法求到精确值.只能在（1.2）取有限项时计算，再估计误差.如果取有限项!!!!)( n s n 1413121111++++++=作为e 的值必然会有误差，根据泰勒余项定理可知其截断误差为e !)1()1( +=-n e s n θ)10(<<θ.如果取（1.2）的前九项，输入程序>> n=8; s=1;S =1; fork=1:n s=s*k; S=S+1/s, end s, S,R=3/(s*(n+1)) 或>>S1=1+1+1/2+1/(1*2*3)+1/(1*2*3*4)+1/(1*2*3*4*5)+1/(1*2*3*4*5*6)+1/(1*2*3*4*5*6*7)+1/(1*2*3*4*5*6*7*8),R1=3/(1*2*3*4*5*6*7*8*9)运行后结果S = R =2.71827876984127 8.267195767195768e-006 因为截断误差为e ),10(101968.267!93!)18()1(6-8<<⨯≈<+=-θθ e s 所以e 的近似值e ≈≈++++++++=!81!71!61!51!41!31!2111)1(8 s 2.718 28.1.2 误差和有效数字例1.2.1 取282.718作为e 的四舍五入近似值时，求其绝对误差和相对误差. 解 在MATLAB 工作窗口输入程序>>juewu=exp(1)-2.71828运行后输出结果为juewu = 1.828 459 045 505 326e-006例1.2.2 计算⎰π20sin x xd x 的近似值，并确定其绝对误差和相对误差.解 因为被积函数xxsin 的原函数不是初等函数，故用泰勒级数求之.++-+-=!!!!sin ９ ７５ ３８６x x x x x x 421 )(∞<<-∞x , （1.5） 这是一个无限过程，计算机无法求到精确值.可用（1.5）的前四项!!!７５ ３６x x x -+-421代替被积函数xxsin ，得 ⎰π=20sin x x y d ⎰π≈20(x !!!142７５ ３６x x x -+-)d x =!7)2(!5)2(!3)2(2753７５ ３⋅π-⋅π+⋅π-π=y ˆ. 根据泰勒余项定理和交错级数收敛性的判别定理，得到绝对误差!99)2(ˆ9⋅<-=πyy R = WU ， 在MA TLAB 命令窗口输入计算程序如下：syms xf=1-x^2/(1*2*3)+x^4/(1*2*3*4*5)-x^6/(1*2*3*4*5*6*7) y=int(f,x,0,pi/2),y1=double(y)y11=pi/2-(pi/2)^3/(3*3*2)+(pi/2)^5/(5*5*4*3*2)-(pi/2)^7/(7*7*6*5*4*3*2)inf=int(sin(x)/x,x,0,pi/2) ,infd=double(inf) WU =(pi/2)^9/(9*9*8*7*6*5*4*3*2), R =infd-y11因为运行后输出结果为： =y 1.370 762 168 154 49，yˆ=1.370 744 664 189 38，=R 1.750 396 510 491 47e-005， WU = 1.782 679 830 970 664e-005410-<.所以，yˆ的绝对误差为=ε410-，故⎰π=20sin x xy d 7 1.370≈x .yˆ的相对误差为 =r ε71.37010ˆ4-=y ε<0.007 3%.1.3 误差估计的基本方法例1.3.4 设计三种算法求方程01522=-+x x 在)3,2(的一个正根*x 的近似值，并研究每种算法的误差传播情况.解 为解已知方程，我们可以设计如下三种算法，然后将计算结果与此方程的精确解5.2*=x 比较，观察误差的传播.算法1 将已知方程化为同解方程=x 2215x -.取初值20=x ，按迭代公式21215k k x x -=+依次计算 ,,,,21n x x x ，结果列入表1–3中.算法2 将已知方程化为同解方程1215+=x x .取初值20=x ，按迭代公式 12151+=+k k x x依次计算 ,,,,21n x x x ，结果列入表1–3中.算法3 将已知方程化为同解方程141522+-+-=x x x x x .取初值20=x ，按迭代公式为1415221+-+-=+k k kk k x x x x x 依次计算 ,,,,21n x x x ，结果列入表1–3中.我们为这三种算法的计算编写两套MATLAB 程序如下： （1）MATLAB 主程序function [k,juecha,xiangcha,xk]= liti112(x0,x1,limax) % 输入的量--x0是初值, limax 是迭代次数和精确值x; % 输出的量--每次迭代次数k 和迭代值xk,% --每次迭代的绝对误差juecha 和相对误差xiangcha ， x(1)=x0;for i=1:limaxx(i+1)=fl(x(i));%程序中调用的fl.m juecha = abs(x(i)-x1);xiangcha = juecha /( abs(x(i))+eps);xk=x(i);k=i-1;[(i-1),juecha,xiangcha,xk] end（2）MATLAB 调用函数程序及其计算结果 ①算法2的MATLAB 调用函数程序function y1=fl(x)y1=15/(2*x+1);② 将MATLAB 主程序和调用函数程序分别命名liti112.m 和fl.m ，分别保存为M 文件，然后在MATLAB 工作窗口输入命令>> [k,juecha,xiangcha,xk]= liti112(2,2.5,100) ③运行后输出计算结果列入表1–3和表 1-4中.④将算法2的MATLAB 调用函数程序的函数分别用y1=15-2*x^2和y1=x-(2*x^2+x-15)/(4*x+1)代替，得到算法1和算法3的调用函数程序，将其保存，运行后将三种算法的前8个迭代值821,,,x x x 列在一起（见表 1-3），进行比较.将三种算法的821,,,x x x 对应的绝对误差和相对误差的值列在一起（见表 1-4），进行比较.1.4 数值计算中应注意的问题例1.4.1 求数）181(71915-+⨯=-x 的近似值. 解 （1）直接用MATLAB 命令>> x=(7^15)*(sqrt(1+8^(-19))-1)运行后输出结果x = 0问题出现在两个相近的数1981-+与1相减时，计算机运行程序>>sqrt(1+8^(-19))-1运行后输出结果ans = 0由于计算机硬件只支持有限位机器数的运算，因此在计算中可能引入和传播舍入误差.因为有效数字的严重损失，导致输出18119-+-的结果为0，计算机不能再与数157继续进行真实的计算，所以，最后输出的结果与x 的精确值不符.（2）如果化为18187)181(71919151915++⨯=-+⨯=---x ，再用MATLAB 命令>> x=(7^15)*( (8^(-19))/(sqrt(1+8^(-19))+1))运行后输出结果x = 1.6471e-005 这是因为18119-+-化为18181919++--后，计算机运行程序>> x= (8^(-19))/(sqrt(1+8^(-19))+1)运行后的结果为x =3.4694e-018 由于有效数字的损失甚少，所以运算的结果-18103.4694⨯再与157继续计算，最后输出的结果与x 的精确值相差无几.例1.4.2 求数)13030ln(2--=y 的近似值. 解 （1）直接用MATLAB 程序>> x=30;x1= sqrt(x^2-1)运行后输出结果x1 = 29.9833 输入MATLAB 程序>> x=30; x1=29.9833;y=log(x-x1)运行后输出结果y = -4.0923（2）因为)13030ln(2--中的30=x 很大，如果采用倒数变换法111221-+=--=x x x x z ，即130301ln)13030ln(22-+=--)190030ln(-+-=.输入MATLAB 程序>> x=30;y=-log(x+sqrt(x^2-1))运行后输出结果y = -4.0941（3）输入MA TLAB 程序>> x=30; y=log(x-sqrt(x^2-1))运行后输出结果y = -4.0941 可见，（2）计算的近似值比（1）的误差小.参加计算的数，有时数量级相差很大.如果不注意采取相应的措施，在它们的加减法运算中，绝对值很小的那个数经常会被绝对值较大的那个数“吃掉”，不能发挥其作用，造成计算结果失真.例1.4.4 请在16位十进制数值精度计算机上利用软件MATLAB 计算下面的两个数0.30.1111111111111111*++=x 和0.30.11111111111111111*++=y将计算结果与准确值比较，解释计算结果.解 在MATLAB 工作窗口输入下面程序>> x=111111*********+0.1+0.3, y=1111111111111111+0.1+0.3运行后输出结果x = 1.111111*********e+014，y =1.111111*********e+015 从输出的结果可以看出，x *x =，而y *y ≠.为什么*y 仅仅比*x 多一位1，而y *y ≠呢？这是因为计算机进行运算时，首先要把参加运算的数写成绝对值小于1而“阶码”相同的数，这一过程称为数的“对阶”.在16位十进制数值精度计算机上利用软件MATLAB 计算这两个数，把运算的数*x 写成浮点规格化形式为，151515*103000**********.0001010000000000000.000100111111111111111.0⨯+⨯+⨯=x在16位十进制数值精度计算机上，三项的数都表示为小数点后面16位数字的数与1510之积，所以，计算机没有对数进行截断，而是按原来的三个数进行计算.因此，计算的结果x *x =.而161616*10030000000000000.00010010000000000000.000101111111111111111.0⨯+⨯+⨯=y三项的数都表示写成绝对值小于1而“阶码”都为1610的数以后，第一项的纯小数的小数点后面有16位数字.但是，后两项数的纯小数的小数点后面有17位数字，超过了16位十进制数值精度计算机的存储量，计算机对后两项的数都进行截断最后一位，即后两项的数都是16位机上的零，再进行计算，所以计算结果与实际不符.。
数值分析(计算方法)总结
第一章绪论误差来源：模型误差、观测误差、截断误差（方法误差）、舍入误差是的绝对误差,是的误差，为的绝对误差限(或误差限）为的相对误差，当较小时，令相对误差绝对值得上限称为相对误差限记为：即:绝对误差有量纲，而相对误差无量纲若近似值的绝对误差限为某一位上的半个单位，且该位直到的第一位非零数字共有n位，则称近似值有n位有效数字，或说精确到该位。
例:设x==3。
1415926…那么,则有效数字为1位，即个位上的3，或说精确到个位.科学计数法：记有n位有效数字，精确到。
由有效数字求相对误差限：设近似值有n位有效数字，则其相对误差限为由相对误差限求有效数字：设近似值的相对误差限为为则它有n位有效数字令1.x+y近似值为和的误差（限）等于误差(限）的和2.x-y近似值为3.xy近似值为4.1．避免两相近数相减2．避免用绝对值很小的数作除数3．避免大数吃小数4．尽量减少计算工作量第二章非线性方程求根1。
逐步搜索法设f (a) <0, f (b）〉 0，有根区间为 (a, b）,从x0=a出发，按某个预定步长（例如h=（b-a）/N）一步一步向右跨，每跨一步进行一次根的搜索，即判别f(x k）=f（a+kh）的符号，若f（x k)〉0(而f（x k-1）<0)，则有根区间缩小为［x k-1，x k] （若f（x k）=0,x k即为所求根)，然后从x k—1出发，把搜索步长再缩小，重复上面步骤，直到满足精度：｜x k—x k-1｜< 为止，此时取x*≈(x k+x k-1）/2作为近似根.2。
二分法设f（x）的有根区间为［a,b]= [a0，b0], f（a）<0，f（b)〉0。
将［a0,b0]对分,中点x0= （(a0+b0)/2），计算f（x0）。
3.比例法一般地,设 [a k，b k］为有根区间，过（a k，f（a k）)、 (b k, f（b k）)作直线，与x轴交于一点x k，则:1.试位法每次迭代比二分法多算一次乘法，而且不保证收敛.2。
- 1、下载文档前请自行甄别文档内容的完整性，平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益，请联系客服反馈,我们会尽快为您处理(人工客服工作时间：9:00-18:30)。
1.1 误差的来源例1.1.1 用差商ha f h a f a f )()()(-+≈'求x x f ln )(=在3=x 处导数的近似值.取1.0=h ，1000.0=h ，h =0.000 000 000 000 001和h =0.000 000 000 000 000 1分别用MATLAB 软件计算，取十五位数字计算.解 在MATLAB 工作窗口输入下面程序>>a=3;h=0.1;y=log(a+h)-log(a);yx=y/h运行后得yx = 0.32789822822991 将此程序中h 改为0.000 1，运行后得yx = 0.33332777790385后者比前者好.再取h = 0.000 000 000 000 001，运行后得yx = 0.44408920985006不如前者好.取h = 0.000 000 000 000 000 1，运行后得yx = 0算出的结果反而毫无价值.例1.1.2 分别求方程组b AX =在下列情况时的解，其中A ⎪⎪⎭⎫⎝⎛=011111.. （1）⎪⎪⎭⎫ ⎝⎛=22b ;（2）⎪⎪⎭⎫⎝⎛=0122.b . 解 （1） 首先将方程组b AX =化为同解方程b A X 1-=，然后在MATLAB 工作窗口输入程序>> b=[2,2]';A=[1,1;1,1.01]; X=A\b运行后输出当⎪⎪⎭⎫ ⎝⎛=22b 时，b AX =的解为⎪⎪⎭⎫⎝⎛=02X ；（2）同理可得，当⎪⎪⎭⎫ ⎝⎛=0122.b 时，b AX =的解为⎪⎪⎭⎫⎝⎛=11X .例1.1.3 计算e 的近似值. 解 泰勒级数e +++++++=!!4!3!21432 n x x x x x n x)(∞<<-∞x , 取1=x ，得e +++++++=!1!41!31!2111 n . （1.2）这是一个无限过程，计算机无法求到精确值.只能在（1.2）取有限项时计算，再估计误差.如果取有限项!!!!)( n s n 1413121111++++++=作为e 的值必然会有误差，根据泰勒余项定理可知其截断误差为e !)1()1( +=-n e s n θ)10(<<θ.如果取（1.2）的前九项，输入程序>> n=8; s=1;S =1; fork=1:n s=s*k; S=S+1/s, end s, S,R=3/(s*(n+1)) 或>>S1=1+1+1/2+1/(1*2*3)+1/(1*2*3*4)+1/(1*2*3*4*5)+1/(1*2*3*4*5*6)+1/(1*2*3*4*5*6*7)+1/(1*2*3*4*5*6*7*8),R1=3/(1*2*3*4*5*6*7*8*9)运行后结果S = R =2.71827876984127 8.267195767195768e-006 因为截断误差为e ),10(101968.267!93!)18()1(6-8<<⨯≈<+=-θθ e s 所以e 的近似值e ≈≈++++++++=!81!71!61!51!41!31!2111)1(8 s 2.718 28.1.2 误差和有效数字例1.2.1 取282.718作为e 的四舍五入近似值时，求其绝对误差和相对误差. 解 在MATLAB 工作窗口输入程序>>juewu=exp(1)-2.71828运行后输出结果为juewu = 1.828 459 045 505 326e-006例1.2.2 计算⎰π20sin x xd x 的近似值，并确定其绝对误差和相对误差.解 因为被积函数xxsin 的原函数不是初等函数，故用泰勒级数求之.++-+-=!!!!sin ９ ７５ ３８６x x x x x x 421 )(∞<<-∞x , （1.5） 这是一个无限过程，计算机无法求到精确值.可用（1.5）的前四项!!!７５ ３６x x x -+-421代替被积函数xxsin ，得 ⎰π=20sin x x y d ⎰π≈20(x !!!142７５ ３６x x x -+-)d x =!7)2(!5)2(!3)2(2753７５ ３⋅π-⋅π+⋅π-π=y ˆ. 根据泰勒余项定理和交错级数收敛性的判别定理，得到绝对误差!99)2(ˆ9⋅<-=πyy R = WU ， 在MA TLAB 命令窗口输入计算程序如下：syms xf=1-x^2/(1*2*3)+x^4/(1*2*3*4*5)-x^6/(1*2*3*4*5*6*7) y=int(f,x,0,pi/2),y1=double(y)y11=pi/2-(pi/2)^3/(3*3*2)+(pi/2)^5/(5*5*4*3*2)-(pi/2)^7/(7*7*6*5*4*3*2)inf=int(sin(x)/x,x,0,pi/2) ,infd=double(inf) WU =(pi/2)^9/(9*9*8*7*6*5*4*3*2), R =infd-y11因为运行后输出结果为： =y 1.370 762 168 154 49，yˆ=1.370 744 664 189 38，=R 1.750 396 510 491 47e-005， WU = 1.782 679 830 970 664e-005410-<.所以，yˆ的绝对误差为=ε410-，故⎰π=20sin x xy d 7 1.370≈x .yˆ的相对误差为 =r ε71.37010ˆ4-=y ε<0.007 3%.1.3 误差估计的基本方法例1.3.4 设计三种算法求方程01522=-+x x 在)3,2(的一个正根*x 的近似值，并研究每种算法的误差传播情况.解 为解已知方程，我们可以设计如下三种算法，然后将计算结果与此方程的精确解5.2*=x 比较，观察误差的传播.算法1 将已知方程化为同解方程=x 2215x -.取初值20=x ，按迭代公式21215k k x x -=+依次计算 ,,,,21n x x x ，结果列入表1–3中.算法2 将已知方程化为同解方程1215+=x x .取初值20=x ，按迭代公式 12151+=+k k x x依次计算 ,,,,21n x x x ，结果列入表1–3中.算法3 将已知方程化为同解方程141522+-+-=x x x x x .取初值20=x ，按迭代公式为1415221+-+-=+k k kk k x x x x x 依次计算 ,,,,21n x x x ，结果列入表1–3中.我们为这三种算法的计算编写两套MATLAB 程序如下： （1）MATLAB 主程序function [k,juecha,xiangcha,xk]= liti112(x0,x1,limax) % 输入的量--x0是初值, limax 是迭代次数和精确值x; % 输出的量--每次迭代次数k 和迭代值xk,% --每次迭代的绝对误差juecha 和相对误差xiangcha ， x(1)=x0;for i=1:limaxx(i+1)=fl(x(i));%程序中调用的fl.m juecha = abs(x(i)-x1);xiangcha = juecha /( abs(x(i))+eps);xk=x(i);k=i-1;[(i-1),juecha,xiangcha,xk] end（2）MATLAB 调用函数程序及其计算结果 ①算法2的MATLAB 调用函数程序function y1=fl(x)y1=15/(2*x+1);② 将MATLAB 主程序和调用函数程序分别命名liti112.m 和fl.m ，分别保存为M 文件，然后在MATLAB 工作窗口输入命令>> [k,juecha,xiangcha,xk]= liti112(2,2.5,100) ③运行后输出计算结果列入表1–3和表 1-4中.④将算法2的MATLAB 调用函数程序的函数分别用y1=15-2*x^2和y1=x-(2*x^2+x-15)/(4*x+1)代替，得到算法1和算法3的调用函数程序，将其保存，运行后将三种算法的前8个迭代值821,,,x x x 列在一起（见表 1-3），进行比较.将三种算法的821,,,x x x 对应的绝对误差和相对误差的值列在一起（见表 1-4），进行比较.1.4 数值计算中应注意的问题例1.4.1 求数）181(71915-+⨯=-x 的近似值. 解 （1）直接用MATLAB 命令>> x=(7^15)*(sqrt(1+8^(-19))-1)运行后输出结果x = 0问题出现在两个相近的数1981-+与1相减时，计算机运行程序>>sqrt(1+8^(-19))-1运行后输出结果ans = 0由于计算机硬件只支持有限位机器数的运算，因此在计算中可能引入和传播舍入误差.因为有效数字的严重损失，导致输出18119-+-的结果为0，计算机不能再与数157继续进行真实的计算，所以，最后输出的结果与x 的精确值不符.（2）如果化为18187)181(71919151915++⨯=-+⨯=---x ，再用MATLAB 命令>> x=(7^15)*( (8^(-19))/(sqrt(1+8^(-19))+1))运行后输出结果x = 1.6471e-005 这是因为18119-+-化为18181919++--后，计算机运行程序>> x= (8^(-19))/(sqrt(1+8^(-19))+1)运行后的结果为x =3.4694e-018 由于有效数字的损失甚少，所以运算的结果-18103.4694⨯再与157继续计算，最后输出的结果与x 的精确值相差无几.例1.4.2 求数)13030ln(2--=y 的近似值. 解 （1）直接用MATLAB 程序>> x=30;x1= sqrt(x^2-1)运行后输出结果x1 = 29.9833 输入MATLAB 程序>> x=30; x1=29.9833;y=log(x-x1)运行后输出结果y = -4.0923（2）因为)13030ln(2--中的30=x 很大，如果采用倒数变换法111221-+=--=x x x x z ，即130301ln)13030ln(22-+=--)190030ln(-+-=.输入MATLAB 程序>> x=30;y=-log(x+sqrt(x^2-1))运行后输出结果y = -4.0941（3）输入MA TLAB 程序>> x=30; y=log(x-sqrt(x^2-1))运行后输出结果y = -4.0941 可见，（2）计算的近似值比（1）的误差小.参加计算的数，有时数量级相差很大.如果不注意采取相应的措施，在它们的加减法运算中，绝对值很小的那个数经常会被绝对值较大的那个数“吃掉”，不能发挥其作用，造成计算结果失真.例1.4.4 请在16位十进制数值精度计算机上利用软件MATLAB 计算下面的两个数0.30.1111111111111111*++=x 和0.30.11111111111111111*++=y将计算结果与准确值比较，解释计算结果.解 在MATLAB 工作窗口输入下面程序>> x=111111*********+0.1+0.3, y=1111111111111111+0.1+0.3运行后输出结果x = 1.111111*********e+014，y =1.111111*********e+015 从输出的结果可以看出，x *x =，而y *y ≠.为什么*y 仅仅比*x 多一位1，而y *y ≠呢？这是因为计算机进行运算时，首先要把参加运算的数写成绝对值小于1而“阶码”相同的数，这一过程称为数的“对阶”.在16位十进制数值精度计算机上利用软件MATLAB 计算这两个数，把运算的数*x 写成浮点规格化形式为，151515*103000**********.0001010000000000000.000100111111111111111.0⨯+⨯+⨯=x在16位十进制数值精度计算机上，三项的数都表示为小数点后面16位数字的数与1510之积，所以，计算机没有对数进行截断，而是按原来的三个数进行计算.因此，计算的结果x *x =.而161616*10030000000000000.00010010000000000000.000101111111111111111.0⨯+⨯+⨯=y三项的数都表示写成绝对值小于1而“阶码”都为1610的数以后，第一项的纯小数的小数点后面有16位数字.但是，后两项数的纯小数的小数点后面有17位数字，超过了16位十进制数值精度计算机的存储量，计算机对后两项的数都进行截断最后一位，即后两项的数都是16位机上的零，再进行计算，所以计算结果与实际不符.。