Matlab与数值分析：摘要

矩阵元素可以是任何matlab表达式，可以是实数，也可以是复数，复数可用特殊函数I，j 输入 a=1 2 3;4 5 6 x=2 pi/2;sqrt(3)3+5i 矩阵元素:矩阵的标识和引用w 矩阵的子矩阵通过向量、标量的标识来引用和赋值。1.向量标识 A（u,v)A中u行，v列元素 2.“01”向量标识 A（L）A的第L行首元规则：规则：w所有字符串都用单引号括起来；所有字符串都用单引号括起来；w字符串中的每个字符都是字符串变量中字符串中的每个字符都是字符串变量中的一个元素；的一个元素；w字符串中的字符以字符串中的字符以ASCIIASCII码形式存储并码形式存储并区分大小区分大小字符串数组字符串数组w 字符串数组的建立。字符串数组的建立。例：例：a=This is a test.A=中国中国 成都成都 B=China 中國中國 w把matlab工作空间中一些有用的数据长久保存下来的方法是生成mat数据文件。save 将工作空间中所有的变量存到matlab.mat文件中。数据的保存与获取默认文件名save filename将工作空间中所有的变量存到filename.mat文件中。save filename a b 将工作空间中a和b变量存到filename.mat文件中。load load filename load filename a b mat文件是标准的二进制文件，还可以ASCII码形式保存。即可恢复保存过的所有变量 矩阵除的运算在线性代数中没有，有矩阵逆的运算，在matlab中有两种矩阵除运算。矩阵左除：AB 定义为：X=AB AX=B 矩阵右除：B/A 定义为：X=B/A XA=B (B/A)=(AB)a p a 自乘p次幂 方阵方阵整数整数3.矩阵乘方当p为正整数时，表示方阵自乘p次；当p为负整数时，表示方阵自乘p次后的逆；当p为零时，为方阵同维的单位阵。winv 矩阵求逆wdet 行列式的值weig 矩阵的特征值wdiag 对角矩阵 w 矩阵转置w sqrt 矩阵开方4.矩阵的其它运算 a./b=b.aa.b=b./aa./b=b.a 都是都是a的元素被的元素被b的对应元的对应元 素除素除a.b=b./a 都是都是b的元素被的元素被a的对应元的对应元 素除素除例例:a=1 2 3;b=4 5 6;c1=a.b;c2=b./ac1=4.0000 2.5000 2.0000c2=4.0000 2.5000 2.0000 给出a,b对应元素间的商.例例:a=1 2 3;b=4 5 6;z=a.2z=1.00 4.00 9.00z=a.bz=1.00 32.00 729.003.数组乘方数组乘方(.)元素对元素的幂元素对元素的幂关系操作和逻辑操作关系操作和逻辑操作w 关系操作符：关系操作符：w、=、=(等于等于)、=(不等于不等于)w 关系运算规则：关系运算规则：w 两个变量都是标量，则结果为真两个变量都是标量，则结果为真(1)或假或假(0)。w 两个变量都是数组，则必须大小相同，结果也两个变量都是数组，则必须大小相同，结果也是同样大小的数组，数组的元素为是同样大小的数组，数组的元素为0或或1。关系操作关系操作w 一个数组和一个标量，则把数组的每个元一个数组和一个标量，则把数组的每个元素分别与标量比较，结果为与数组大小相素分别与标量比较，结果为与数组大小相同的数组，数组的元素为同的数组，数组的元素为0或或1。w、=，仅对参加比较变量的实，仅对参加比较变量的实部进行比较，部进行比较，w=和和=，则同时对实部和虚部进行，则同时对实部和虚部进行比较。比较。逻辑操作逻辑操作w 逻辑操作符：w&(与)、|(或)、(非)和xor(异或)。逻辑操作逻辑操作w 逻辑操作符：w&(与)、|(或)、(非)和xor(异或)。多项式的表示和创建 在Matlab中，任意的多项式都是用一个行向量表示，将多项式的系数按降幂排列存放在行向量中。多项式p(x)=a0 xn+a1xn-1+an-1x+an对应的行向量为:P=a0 a1 an。通过根创建多项式polyw其调用格式如下：lP=poly(A)l如A为向量，则P是多项式的系数行向量，该多项式的根A。l注意：l如A为nn的矩阵，则P是一个长度为n+1的行向量，它是A的特征多项式的系数；例：A=6-8 6;1 0 0;0 2 0;p=poly(A),roots(p),D=eig(A)w a=1-6 15-20 15-6 1w a=w 1 -6 15 -20 15 -6 1w roots(a)w ans=w 1.0042+0.0025iw 1.0042-0.0025iw 1.0000+0.0049iw 1.0000-0.0049iw 0.9958+0.0024iw 0.9958-0.0024ilC=conv(a,b)。例：求a=x3+3x2+2x1+1和b=4x3+3x2+9x1+10的乘积。a=1 3 2 1;b=4 3 9 10;c=conv(a,b)2.多项式的除法dconvw多项式的除法的运算实质就是多项式系数的解卷积运算。向量a对向量c进行解卷积得到的商向量q和余向量r，满足：knc(k)-r(k)=a(i)q(k+1-i)i=1n在Matlab，完成此功能的命令是dconv：lq,r=dconv(c,a)其中：c=conv(a,q)+r4.多项式的求导polyder 不仅可以计算单个多项式的导数，还可以计算两个多项式相乘和相除的导数。其调用方式如下：lpolyder(p):返回多项式系数向量p的导数；lpolyder(a,b):返回多项式a*b的导数；lq,d=polyder(b,a):返回b/a的导数，用q/d表示。其中q是结果的分子多项式，d是结果的分母多项式。5.多项式的求值polyvaln在Matlab中，多项式求值有两种方式：按数组规则和按矩阵规则计算多项式值。ly=polyval(p,x):按数组规则运算。用来计算多项式p在自变量x处的值，可以是数、向量和矩阵。如果x是向量或矩阵，则该命令将对x的每个元素计算p的值，结果大小和x的大小相同；ly=polyvalm(p,x):按矩阵规则运算。自变量x是一方阵。l他的计算方法是：lY=P(1)*XN+P(2)*X(N-1)+.+P(N)*X+P(N+1)*I6.多项式的部分分式展开没有多重根:w B(s)R(1)R(2)R(n)w -=-+-+.+-+K(s)w A(s)s-P(1)s-P(2)s-P(n)w 有多重根时，则包含有：w R(j)R(j+1)R(j+m-1)w -+-+.+-w s-P(j)(s-P(j)2 (s-P(j)mw w 在Matlab中完成此任务的指令是residue:lr,p,k=residue(b,a):求b(s)/a(s)的部分分式展开，返回值r是留数，p是极点，k是整式向量，a,b分别是分母和分子；lb,a=residue(r,p,k):将展开式还原为两个多项式b(s)和a(s)相除的形式。plot 最基本的二维图形指令最基本的二维图形指令lplot的功能：的功能：l plot命令自动打开一个图形窗口命令自动打开一个图形窗口Figurel 用直线连接相邻两数据点来绘制图形用直线连接相邻两数据点来绘制图形l根据图形坐标大小自动缩扩坐标轴，将数根据图形坐标大小自动缩扩坐标轴，将数据标尺及单位标注自动加到两个坐标轴上据标尺及单位标注自动加到两个坐标轴上，可自定坐标轴，可把，可自定坐标轴，可把x,y 轴用对数坐标轴用对数坐标表示表示plot(x,y):参数x,y为向量，以x为横坐标，以y为纵坐标。注意注意：x,y的长度必须相等例：t=0:0.001:2*pi;y=cos(t);plot(t,y)lplot(Y)参数Y为mn的矩阵。该命令把矩阵Y转化为长度为m的n个向量，然后对n个向量划出n条线。例：y=5 4 4 8 9 10;3 4 4 5 8 2;8 5 7 9 6 25;9 8 8 9 6 7 plot(y)l复向量式 plot(Z):当z为一个复向量时，该命令等价于plot(real(Y),imag(Y);其它情况忽略z的虚部。但会报错。例：t=0:0.1:3*pi;x=t+i*cos(t);plot(x);figure;plot(t,x)l综合调用方式plot(x1,y1,x2,y2,)该命令可以在同一窗口绘制多条曲线。例：t=0:0.1:3*pi;x=t+i*cos(t);plot(real(x),imag(x),t,x)三维绘图可分为三个过程：三维绘图可分为三个过程：l 生成平面网格点生成平面网格点l 生成和绘制三维网格点生成和绘制三维网格点l 处理三维表面。处理三维表面。三维曲面绘图三维曲面绘图生成平面网格点 X,Y=meshgrid(x,y)lx 是区间是区间x0,xm上划分点组成的向量；上划分点组成的向量；ly 是区间是区间y0,yn上划分点组成的向量；上划分点组成的向量；lX,Y是输出变量矩阵，矩阵是输出变量矩阵，矩阵X的行向量都的行向量都是是x，矩阵，矩阵y的列向量都是的列向量都是y。例：例：X,Y=meshgrid(-2:2:2,-2:1:2)l生成三维网格点 mesh(x,y,z,c)：x,y,z,c是同维矩阵是同维矩阵 x,y,z对应确定空间上的网格点对应确定空间上的网格点 c为颜色矩阵；为颜色矩阵；例：例：X,Y=meshgrid(-2:.2:2,-2:.2:2);Z=X.*exp(-X.2-Y.2);mesh(Z)l三维表面处理 对网格曲面的网格块区域着色就得对网格曲面的网格块区域着色就得到实曲面。到实曲面。mesh命令仅对网格线着色，命令仅对网格线着色，而对网格片着色使用而对网格片着色使用surf命令。命令。surf的的使用方式与使用方式与mesh相同。相同。例：例：X,Y=meshgrid(-2:.2:2,-2:.2:2);Z=X.*exp(-X.2-Y.2);surf(Z);lshading faceted 网格修饰，缺省方式lshading flat 去掉黑色线条，根据小方块的值确定颜色lshading interp 颜色整体改变，根据小方块四角的值差补过度点的值确定颜色X,Y,Z=peaks(30);surfl(X,Y,Z)shading interp;colormap(cool);axis offhidden onl命令文件格式命令文件格式注解运行命令注解运行命令将matlab命令构成的程序存储成以m为扩展名的文件，然后再执行该程序文件，这种工作模式称为命令文件模式。程序文件不能在指令窗口下建立，因为指令窗口只允许一次执行一行上的一个或几个语句。误差的来源误差的来源现现 实实 世世 界界 研 究对 象测量数据数学模型的建立计算方法数值运算测量误差模 型误 差方 法误 差舍 入误 差绝对误差绝对误差:l设x是精确值，p是近似值，则定义两者之差是绝对误差绝对误差:为近似数p的绝对误差,简称误差误差axp 绝对误差限绝对误差限a-a-a+a+a aAn由于精确值一般是未知的,因而不能求出来,但可以根据测量误差或计算情况估计它的上限|-|x p称为绝对误差限。相对误差相对误差l相对误差定义为绝对误差与精确值之比 l表示单位长度中含有多少误差arx n由于x一般不知道,实际计算时用近似值p代替。arp n绝对误差是有量纲的量，相对误差没有量纲,有时亦用百分比、千分比表示相对误差限相对误差限arx 称为相对误差限误差来源：误差来源：舍入误差舍入误差 A=a0 a1 am.am+1 am+n am+n+1 高位部分高位部分低位部分低位部分表示成表示成A=a0 a1 am.am+1 am+n +0 0.0 0 am+n+1 n将无限位字长的精确数处理成有限位字长近似数的处理方法称为舍入方法n如果数A的字长要求保留到小数后n位,可以将数切分成高位部分和低位部分l舍去低位部分，即取 a=a0 a1 am.am+1 am+n，那么截断近似法截断近似法|=0.0 0 am+n+1 0.0 0 999.0.0 0 1=1x10-nn n位位n n位位n n位位n如果近似数a的绝对误差限不超过最末位的1个单位，称a为准确到小数后第n位的可靠数，每位数字均为可靠数字四舍五入法四舍五入法|=0.0 0 am+n+1 0.0 0 499.0.0 0 5=0.5x10-n1,2,3,4n n位位n n位位n n位位n四舍情况，当am+n+1=1,2,3,4时，取 a=a0 a1 am.am+1 am+n，那么四舍五入法四舍五入法l五入情况l当am+n+1=5,6,7,8,9时，取 a=a0 a1 am.am+1 (am+n+1),那么|=0.0 0 1-0.0 0 am+n+1 0.0 0 1-0.0 0 50.=0.5x10-n5,6,7,8,9n n位位n n位位n n位位n n位位n四舍五入法的|0.5x10-n，在a的最末一位上有半个单位误差有效数字概念有效数字概念l对于 a=a0 a1 am.am+1 am+n(a00)的近似数，若|0.5x10-n，则称a为具有m+n+1位有效数字的有效数，其中每一位数字都叫做a的有效数字。有效数和可靠数的最末位数字称为可疑数字n有效数位的多少直接影响到近似值的绝对误差与相对误差的大小关于有效数字的推论关于有效数字的推论mnaaax10.021n推论1 对于给出的有效数，其绝对误差限不大于其最末数字的半个单位1102m nxx 关于有效数字的推论关于有效数字的推论15()10nrxan推论2 对于给出的一个有效数，其相对误差限可估计如下：mnaaax10.021:例：例：的有限位数的有限位数取取 x1 1=3=3，误差限不超过，误差限不超过0.5;0.5;取取 x2 2=3.14,=3.14,误差限不超过误差限不超过0.0050.005;取取 x3 3=3.1416,=3.1416,误差限不超过误差限不超过0.000050.00005;分析：分析：算式系解析方程导出算式系解析方程导出 计算过程中未见有相近大数相减计算过程中未见有相近大数相减 剩下的可能性是从初始值起的舍入误差剩下的可能性是从初始值起的舍入误差数值计算数值计算问题的适定性问题的适定性l“良态良态”问题和问题和“病态病态”问题问题l在适定的情况下，若对于原始数据很小的变化在适定的情况下，若对于原始数据很小的变化XX，对应的参数误差对应的参数误差yy也很小，则称该数学问题是也很小，则称该数学问题是良态良态问题；若问题；若yy很大，则称为病态问题很大，则称为病态问题l病态问题中解对于数据的变化率都很大，因此数据微病态问题中解对于数据的变化率都很大，因此数据微小变化必将导致参数模型精确解的很大变化小变化必将导致参数模型精确解的很大变化l数学问题的性态完全取决于该数学问题本身的属性，数学问题的性态完全取决于该数学问题本身的属性，在采用数值方法求解之前就存在，与数值方法无关在采用数值方法求解之前就存在，与数值方法无关数值计算数值计算问题的适定性问题的适定性l稳定算法和不稳定算法稳定算法和不稳定算法l如果用数值方法计算时如果用数值方法计算时,初始误差对结果影响小初始误差对结果影响小的算法称为稳定算法。否则称为不稳定算法的算法称为稳定算法。否则称为不稳定算法例例:计算积分值计算积分值 1011010)1exp()1exp()1exp(dxxxnxxdxxxnnn,3,2,11nAnAnn,2,1,)1exp(10ndxxxAnneexddxxxA1111exp1)1exp(11101001经经A1出发出发,递推计算得到：递推计算得到：123456780.3 6 7 8 7 90.2 6 4 2 4 20.2 0 7 2 7 40.1 7 0 9 0 40.1 4 5 4 8 00.1 2 7 1 2 00.1 1 0 1 6 00.1 1 8 7 2 0AAAAAAAA?068480.09A?068480.09A90,1e x p(1)0 xxx 分析：分析：算式系解析方程导出算式系解析方程导出 计算过程中未见有相近大数相减计算过程中未见有相近大数相减 剩下的可能性是从初始值起的舍入误差剩下的可能性是从初始值起的舍入误差4.注意事项注意事项几点注意事项几点注意事项1.避免相近二数相减避免相近二数相减例：例：a1=0.12345，a2=0.12346，各有，各有5位有效数字。位有效数字。而而 a2 a1=0.00001，只剩下，只剩下1位有效数字。位有效数字。几点注意事项几点注意事项 几种经验性避免方法：几种经验性避免方法：;xxxx ;1lnlnln xxx当当|x|1 时：时：;2sin2cos12xx .6121112xxxex4 Remarks2.避免小分母避免小分母:分母小会造成浮点溢出分母小会造成浮点溢出3.避免大数避免大数吃吃小数小数4 Remarks例：用单精度计算例：用单精度计算 的根。的根。010)110(992 xx精确解为精确解为110291 x,x 算法算法1 1：利用求根公式利用求根公式aacbbx242 229124410,022bbacbbacxxaa 4 Remarks在计算机内，在计算机内，109存为存为0.1 1010，1存为存为0.1 101。做加法时，做加法时，两加数的指数先向大指数对齐，再将浮点部分相加。即两加数的指数先向大指数对齐，再将浮点部分相加。即1 的指数部分须变为的指数部分须变为1010，则：，则：1=0.0000000001 1010，取，取单精度时就成为：单精度时就成为：109+1=0.10000000 1010+0.00000000 1010=0.10000000 1010大数大数吃吃小数小数算法算法2：先解出先解出 再利用再利用9211024)(aacbbsignbx11010991221 xacxacxx求和时求和时从小到大从小到大相加，可使和的误差减小。相加，可使和的误差减小。例：按从小到大、以及从大到小的顺序分别计算例：按从小到大、以及从大到小的顺序分别计算1+2+3+40+1094.先化简再计算，减少步骤，避免误差积累。先化简再计算，减少步骤，避免误差积累。一般来说，计算机处理下列运算的速度为一般来说，计算机处理下列运算的速度为 exp ,5.选用稳定的算法。选用稳定的算法。