ch3数值阵列及向量化运算

1第3章数值阵列及向量化运算3.1数值计算的特点和地位符号计算的局限性：1）有很多问题无法解，2）求解时间过长数值计算：适用范围广，能处理各种复杂的函数关系，计算速度快，容量大。【例3.1-1】已知tttfcos)(2，求dttfxsx0)()(。（1）符号计算解法symstxft=t^2*cos(t)sx=int(ft,t,0,x)ft=t^2*cos(t)sx=sin(x)*(x^2-2)+2*x*cos(x)（2）数值计算解法dt=0.05;t=0:dt:5;Ft=t.^2.*cos(t);Sx=dt*cumtrapz(Ft);%由一个个宽度为dt的小梯形面积的累加求Ft曲线下的面积t(end-4:end)Sx(end-4:end)plot(t,Sx,'.k','MarkerSize',12)xlabel('x'),ylabel('Sx'),gridonans=4.80004.85004.90004.95005.0000ans=-20.1144-19.9833-19.7907-19.5345-19.213100.511.522.533.544.55-25-20-15-10-505xSx图3.1-1在区间[0,5]采样点上算得的定积分值2【例3.1-2】已知)sin()(tetf，求40)()(dttfxs。（1）符号计算解法symstxft=exp(-sin(t))sx=int(ft,t,0,4)ft=exp(-sin(t))Warning:Explicitintegralcouldnotbefound.sx=int(exp(-sin(t)),t==0..4)（2）数值计算解法dt=0.05;t=0:dt:4;Ft=exp(-sin(t));Sx=dt*cumtrapz(Ft);Sx(end)plot(t,Ft,'*r','MarkerSize',4)holdonplot(t,Sx,'.k','MarkerSize',15)holdoffxlabel('x')legend('Ft','Sx')ans=3.063200.511.522.533.5400.511.522.533.5xFtSx3.2数值数组的创建和寻访3.2.1一维数组的创建x=[1,3,5,7,9]逐个元素输入法x=a:inc:b冒号生成法，inc缺省时步长为1x=linspace(a,b,n)线性定点法x=logspace(a,b,n)对数定点法运用标准数组生成函数。3【例3.2-1】一维数组的常用创建方法举例。a1=1:6a2=0:pi/4:pia3=1:-0.1:0a1=123456a2=00.78541.57082.35623.1416a3=Columns1through61.00000.90000.80000.70000.60000.5000Columns7through110.40000.30000.20000.10000b1=linspace(0,pi,4)b2=logspace(0,3,4)%创建数组[100101102103]b1=01.04722.09443.1416b2=1101001000c1=[2pi/2sqrt(3)3+5i]c1=Columns1through32.00001.57081.7321Column43.0000+5.0000irand('state',0)%利用标准数组生成函数产生均匀分布随机数组c2=rand(1,5)c2=0.95010.23110.60680.48600.89133.2.2二维数组的创建一小规模数组的直接输入法【例3.2-2】在MATLAB环境下，用下面三条指令创建二维数组C。a=2.7358;b=33/79;C=[1,2*a+i*b,b*sqrt(a);sin(pi/4),a+5*b,3.5+i]C=1.00005.4716+0.4177i0.69090.70714.82443.5000+1.0000i三个要素：整个输入数组必须以方括号“[]”为其首尾；行与行间“;”或“Enter”；同行中元素间“,”或“空格”。4二中规模数组的数组编辑器创建法【例3.2-3】根据现有数据创建一个)83(的数组。图3.2-1利用数组编辑器创建中规模数组三中规模数组的M文件创建法【例3.2-4】创建和保存数组AM的MyMatrix.m文件。（1）打开文件编辑调试器，并在空白填写框中输入所需数组（见图3.2-2）。（2）最好，在文件的首行，编写文件名和简短说明，以便查阅（见图3.2-2）。（3）保存此文件，并且文件起名为MyMatrix.m。（4）以后只要在MATLAB指令窗中，运行MyMatrix.m文件，数组AM就会自动生成于MATLAB内存中。图3.2-2利用M文件创建数组四利用MATLAB函数创建数组【例3.2-5】利用最常用标准数组生成函数产生标准数组的演示。ones(2,4)%产生（2×4）全1数组ans=11111111randn('state',0)%把正态随机数发生器置0randn(2,3)%产生正态随机阵ans=-0.43260.1253-1.1465-1.66560.28771.19095D=eye(3)%产生3×3的单位阵D=100010001diag(D)%取D阵的对角元ans=111diag(diag(D))%外diag利用一维数组生成对角阵ans=100010001randsrc(3,20,[-3,-1,1,3],1)%在[-3,-1,1,3]上产生3×20均布随机数组，随机发生器的状态设置为1ans=Columns1through13313-1-33-3-3-13-1-1-31313-1111111133-1-3-11-13-1-111-1-3Columns14through2011-33-113-1-1-3-1-11-33-1-13-1333.2.3二维数组元素的编址和寻访6【例3.2-6】本例演示：数组元素及子数组的各种标识和寻访格式；冒号的使用；end的作用。A=zeros(2,6)A(:)=1:12%单下标法：单下标全元素寻访A=000000000000A=135791124681012A(2,4)%全下标法：指定行、指定列A(8)%单下标法：单下标寻访ans=8ans=8A(:,[1,3])%全下标法：全部行、指定列A([1,2,5,6]')%单下标法:生成指定的一维行（或列）数组ans=1526ans=1256A(:,4:end)%全下标法：全部行、指定列，end表示最后一列。ans=791181012A(2,1:2:5)=[-1,-3,-5]%全下标法：指定行、指定列A=1357911-14-38-512B=A([1,2,2,2],[1,3,5])%全下标法：指定行、指定列B=159-1-3-5-1-3-5-1-3-5相互转换指令：ind2sub：把单序号编址转换成全下标编址sub2ind：把全下标编址转换成到单序号编址73.2.4数组构作技法综合【例3.2-7】数组操作函数reshape,diag,repmat的用法；空阵[]删除子数组的用法。a=1:8A=reshape(a,4,2)A=reshape(A,2,4)%改变行数和列数a=12345678A=15263748A=13572468b=diag(A)%提取对角元素，。B=diag(b)%生成对角阵b=14B=1004D1=repmat(B,2,4)%排列B模块D1=10101010040404041010101004040404D1([1,3],:)=[]%删除指定行D1=0404040404040404【例3.2-8】函数flipud,fliplr,rot90对数组的操作体现着“矩阵变换”。A=reshape(1:9,3,3)A=147258369B=flipud(A)%上下对称交换B=3692581478C=fliplr(A)%左右对称交换C=741852963D=rot90(A,2)%逆时针旋转90度，2次D=9638527413.3数组、矩阵及向量化编程MATLAB面向数组/矩阵编程和运算：用“数组或矩阵运算”模式去处理那些“借助循环而反复执行的标量运算”显著提高程序执行速度书写简洁、便于阅读3.3.1数组概念和运算规则（3-2）数组运算符和矩阵运算符加减乘左除右除求幂数组运算符.+.-.*.\./.^矩阵运算符+-*\/^9服从数组运算规则的MATLAB初等函数及关系逻辑算符举例服从数组运算规则的初等函数三角、反三角sin,cos,tan,cot,sec,csc,asin,acos,atan,acos,asec,acsc双曲、反双曲sinh,cosh,…,asinh,acosh,…指数、对数exp,sqrt,pow2,log,log10,log2圆整、求余ceil,floor,fix,round,mod,rem模、角、虚实部abs,angle,real,imag,conj符号函数sign关系运算符==,~=,,,=,=逻辑运算符&,|,~3.3.2矩阵概念及矩阵运算规则一矩阵概念NMijMNMMNNNMaaaaaaaaaa212222111211A（3-1）二矩阵运算规则【例】已知矩阵42A，34B，采用两种不同的编程求这两个矩阵的乘积344232BAC。（1）以“标量”为基本处理单元的传统编程法njmibaclkkjikij,,1;,,11clearrng('default')%随机发生器恢复默认设置A=rand(2,4);B=rand(4,3);C1=zeros(size(A,1),size(B,2));forii=1:size(A,1)%6forjj=1:size(B,2)fork=1:size(A,2)C1(ii,jj)=C1(ii,jj)+A(ii,k)*B(k,jj);endendend%12C1C1=1.27261.38701.22812.29481.46601.8204（2）MATLAB求矩阵乘积C2=A*B%13C2=1.27261.38701.22812.29481.46601.8204103.3.3数组、矩阵算术运算规则的功能比较相同：加，减，与标量加、减、乘不同：转置，乘，除，求幂【例】本例演示：矩阵和数组乘法的不同。（1）Am=magic(3)%3阶魔方阵Aa=reshape(1:12,3,4)B=2*ones(3,4)Am=816357492Aa=147102581136912B=222222222222（2）Cm1=Am*BCa1=Aa.*BCm1=303030303030303030303030Ca1=281420410162261218243.3.4向量化编程尽可能用“数组或矩阵运算”指令【例】欧姆定律：iur，其中iur,,分别是电阻（欧姆）、电压（伏特）、电流（安培）。验证实验：据电阻两端施加的电压，测量电阻中流过的电流，然后据测得的电压、电流计算平均电阻值。（测得的电压电流具体数据见下列程序）。（1）非向量化程序clearvr=[0.89,1.20,3.09,4.27,3.62,7.71,8.99,7.92,9.70,10.41];ir=[0.028,0.040,0.100,0.145,0.118,0.258,0.299,0.257,0.308,0.345];%--------------------L=length(vr);fork=1:Lr(k)=vr(k)/ir(k);11end%---------------------------sr=0;fork=1:Lsr=sr+r(k);endrm=sr/Lrm=30.5247（2）向量化程序clearvr=[0.89,1.20,3.09,4.27,3.62,7.71,8.99,7.92,9.70,10.41];ir=[0.028,0.040,0.100,0.145,0.