1、1第二讲 MATLAB 的数值分析2-1 矩阵运算与数组运算 矩阵运算和数组运算是 MATLAB 数值运算的两大类型,矩阵运算是按矩阵的运算规则进行的,而数组运算则是按数组元素逐一进行的。因此,在进行某些运算(如乘、除)时,矩阵运算和数组运算有着较大的差别。在 MATLAB 中,可以对矩阵进行数组运算,这时是把矩阵视为数组,运算按数组的运算规则。也可以对数组进行矩阵运算,这时是把数组视为矩阵,运算按矩阵的运算规则进行。1、矩阵加减与数组加减矩阵加减与数组加减运算效果一致,运算符也相同,可分为两种情况:(1)若参与运算的两矩阵(数组)的维数相同,则加减运算的结果是将两矩阵的对应元素进行加减,如A
2、=1 1 1;2 2 2;3 3 3;B=A;A+Bans=2 2 24 4 46 6 6(2)若参与运算的两矩阵之一为标量(1*1 的矩阵) ,则加减运算的结果是将矩阵(数组)的每一元素与该标量逐一相加减,如A=1 1 1;2 2 2;3 3 3;A+2ans=3 3 34 4 45 5 52、矩阵乘与数组乘(1)矩阵乘矩阵乘与数组乘有着较大差别,运算结果也完全不同。矩阵乘的运算符为“*” ,运算是按矩阵的乘法规则进行,即参与乘运算的两矩阵的内维必须相同。设 A、B 为参与乘运算的两矩阵,C 为 A 和 B 的矩阵乘的结果,则它们必须满足关系 Cmn=AmkBkn。因此,参与运算的两矩阵的顺
3、序不能任意调换,因为 A*B 和 B*A 计算结果很可能是完全不一样的。如:A=1 1 1;2 2 2;3 3 3;2B=A;A*Bans=6 6 612 12 1218 18 18F=ones(1,3);G=ones(3,1);F*Gans3G*Fans=1 1 11 1 11 1 1(2)数组乘数组乘的运算符为“.*” ,运算符中的点号不能遗漏,也不能随意加空格符。参加数组乘运算的两数组的大小必须相等(即同维数组) 。数组乘的结果是将两同维数组(矩阵)的对应元素逐一相乘,因此,A.*B 和 B.*A 的计算结果是完全相同的,如:A=1 1 1 1 1;2 2 2 2 2;3 3 3 3 3
4、;B=A;A.*Bans=1 1 1 1 14 4 4 4 49 9 9 9 9B.*Aans=1 1 1 1 14 4 4 4 49 9 9 9 9由于矩阵运算和数组运算的差异,能进行数组乘运算的两矩阵,不一定能进行矩阵乘运算。如A=ones(1,3);B=A;A.*Bans=1 1 1 A*A3?Error using= =Inner matrix dimensions must agree.3、矩阵除与数组除矩阵除分为矩阵右除和矩阵左除两种情况。矩阵右除的运算符为“/” ,设 A、B 为两矩阵,则“A/B”是指方程 X*B=A 的解矩阵 X。显然,矩阵右除运算对参与运算两矩阵的维数是有一
5、定要求的,即矩阵 A 和 B 的列数必须相等。如A=1 1 1 1;2 2 2 2;3 3 3 3;B=1 1 1 1;X=A/BX=12X*Bans=1 1 1 12 2 2 2矩阵右除允许参与右除运算的矩阵 B 为标量,这时矩阵右除运算的结果是将矩阵 A 的每一元素逐一与该标量进行除法运算。如:A=2 4 6 8;8 6 4 2;B=2;A/Bans=1 2 3 44 3 2 1矩阵左除运算符为“” ,设 A、B 为两矩阵,则“AB”是指方程 B*X=A 的解矩阵 X。显然,矩阵左除运算对参与运算两矩阵的维数也有一定要求的,即矩阵 A 和 B 的行数必须相等。如数组右除的运算符为“./”
6、,左除的运算符为“.” 。数组右除和左除的运算结果是完全等效的。设 A、B 为两同维矩阵,则“A./B”的运算结果是将矩阵 A 的每一个元素与矩阵 B的对应元素相除。应注意的是,参与数组运算的两矩阵(数组)的大小必须相等。A=2 2 3 3 4 4;1 1 2 2 3 3;4 4 5 5 6 6;B=1 2 3 3 2 2;1 1 1 1 1 1;2 2 5 5 3 3;A./Bans=2 2 1 1 2 21 1 2 2 3 32 2 1 1 2 2B./Aans=40.5000 0.5000 1.0000 1.0000 0.5000 0.50001.0000 1.0000 0.5000 0
7、.5000 0.3333 0.33330.5000 0.5000 1.0000 1.0000 0.5000 0.50004、常用的矩阵运算函数(1)用 size()函数计算矩阵 A 的维数,调用格式:d=size(A) %将矩阵 A 的行数和列数赋给变量 dm,n=size(A) %将矩阵 A 的行数赋给变量 m、列数赋给变量 n(2)用 rand()函数产生随机矩阵, ,调用格式:rand(n) %产生值在 01 之间随机分布的 n*n 的随机方阵rand(m,n) %生值在 01 之间随机分布的 n*m 的随机矩阵(3)计算矩阵长度(列数)的函数 length(),调用格式:a=lengt
8、h(B) %将矩阵 B 的列数赋值给变量 a(4)矩阵元素的求积运算函数 prod(),调用格式:prod(A) %若 A 为向量,将计算矩阵 A 所有元素之积;若 A 为矩阵,将产生一行向量,其元素分别为矩阵 A 的各列元素之积。prod(A,k) %将对矩阵 A 按 k 定义的方向进行示积运算,若 k=1 则按列的方向求积,若 k=2 则按行的方向求积。(5)矩阵元素的求和运算函数 sum(),调用格式同 prod()函数。2-2 多项式及其运算 MATLAB 提供了标准多项式的常用函数,包括求根、相乘、相除等。1多项式的表达与创建MATLAB 采用将多项式系数按幂次序排列形式的行向量来表
9、征一多项式。设多项式为: 011assasAnn)(则表征该多项式的行向量表示为: 。P因此,在 MATLAB 中,创建多项式即可创用建行向量的方法,直接输入按顺序排序的多项式系数即可。如,输入语句:A=1 2 2 1;即表示创建多项式 ,并赋值给变量 A。123ss2多项式求根函数 roots()用于对多项式求根,调用格式:p=roots(A)5其中 A 为表征多项式的行向量,p 返回该多项式的根(用列向量表示) 。如B=1 3 2; %创建多项式 232sroots(B) %求多项式的根ans=-2 -1A=1 2 2 1;roots(A)ans=-1.0000 -0.5000+0.866
10、0i-0.5000-0.8660i3由指定根求多项式函数 poly 用于由给定根求多项式系数向量,调用格式为:A=poly(p)其中 p 为多项式根(行或列向量表示) ,A 为返回的多项式系数(行向量表示) 。如:p=2 1;poly (p)ans=1 -3 2可见 roots()与 poly()互为逆运算。4多项式相乘(卷积)函数 conv()用于求两个多项式的乘积多项式,调用格式:R=conv(A,B)其中 A、B 分别为表征两个多项式的行向量,R 为返回的每间积多项式的系数向量(按降幂次序排列) 。如A=1 3 2;B=1 2 1;R=conv(A,B)R=1 5 9 7 25多项式相除
11、(解卷)函数 deconv()用于进行两个多项式的相除运算,它是相乘运算(conv )的逆运算,其调用格式为:B,t=deconv(R,A)其中 R 为被除多项式,A 为除数多项式,B 为商多项式,t 为余多项式。即多项式 R 除以多项式 A 后得商多项式 B 和余多项式 t。如:R=1 5 9 7 2;6A=1 3 2;B,t=deconv(R,A)B=1 2 1t=0 0 0 0 0若多项式系数向量为零向量,则表示 R 能被 A 除尽。2-3 线性代数方程(组)求解我们习惯将上组方程式以矩阵方式表示如下:AX=B 其中 A 为等式左边各方程式的系数项,X 为欲求解的未知项,B 代表等式右边
12、之已知项,要解上述的联立方程式,我们可以利用矩阵左除 做运算,即是 X=AB。 如果将原方程式改写成 XA=B 其中 A 为等式左边各方程式的系数项,X 为欲求解的未知项,B 代表等式右边之已知项。注意上式的 X, B 已改写成列向量,A 其实是前一个方程式中 A 的转置矩阵。上式的 X 可以矩阵右除 / 求解,即是 X=B/A。 若以反矩阵运算求解 AX=B, X=B,即是 X=inv(A)*B,或是改写成 XA=B, X=B,即是X=B*inv(A)。 我们直接以下面的例子来说明这三个运算的用法: A=3 2 -1; -1 3 2; 1 -1 -1; % 将等式的左边系数键入 B=10 5
13、 -1; % 将等式右边之已知项键入,B 要做转置 X=AB % 先以左除运算求解 X = % 注意 X 为行向量 -2 5 6 C=A*X % 验算解是否正确 C = % C=B 10 5 -1 A=A; % 将 A 先做转置 B=10 5 -1; X=B/A % 以右除运算求解的结果亦同 X = % 注意 X 为列向量 10 5 -1 7X=B*inv(A); % 也可以反矩阵运算求解 2-4 信号的 MATLAB 表示举例1、指数信号 atAeMATLAB 中表示:y=A*exp(a*t)%example2-1 decaying exponential signalA=1;a=-0.4;
14、t=0:0.001:10;yt= A*exp(a*t);plot(t,yt)2、正弦信号 )sin(tA0MATLAB 中表示:y=A*sin(w0*t+phi)%example2-2 sinusoidal signalA=1;w0=2*pi;phi=pi/6;t=0:0.001:8;yt= A*sin(w0*t+phi);plot(t,yt)3、抽样函数 )(tSaMATLAB 中抽样函数用 sinc 函数表示:y=sinc(t)%example2-3 Sample functiont=-3*pi:pi/100:3*pi;yt=sinc(t/pi);plot(t,yt)4、矩形脉冲信号MAT
15、LAB 中矩形脉冲信号用 rectpuls 函数表示:y=rectpuls (t,width) %width 缺省值为 1%example2-4 产生一幅度为 1,宽度为 width 以零点为对称的矩形波t=0:0.001:4;T=1;yt=rectpuls (t-2*T,T);8plot(t,yt)axis(0,4,-0.5,1.5)5、三角波脉冲信号MATLAB 中三角波脉冲信号用 tripuls 函数表示:y=tripuls (t,width,skew) 用以产生一幅度最大为 1,宽度为 width 的三角波,非零范围( -width/2, width/2) 。skew定义 skew=2
16、tmax/ width,t max 为三角波顶点坐标,skew 取值:-1+1 ,决定三角波的形状,缺省时为 0。%example2-5 t=-4:0.001:4;yt=tripuls (t,4,0.5);plot(t,yt)6、周期矩形脉冲信号y=square(w0*t,duty_cycle);用以产生一个幅度为+1 和-1,基波频率为 w0 的矩形脉冲信号。duty_cycle 指一个周期内正脉冲的宽度和负脉冲宽度的百分比,缺省值为 1。%example2-6 square wavet=0:0.0001:5;A=1;T=1;w0=2*pi/T;ft=A*square(w0*t,20);pl
17、ot(t,ft)axis(0,5,-1.5,1.5)7、周期三角波信号y=sawtooth(w0*t,width);产生一个基波频率为 w0 的三角波信号,其幅度为 +1 和-1,-1 出现在 nT 点。设在0,T区间内+1 出现的位置为 tmax,则 width =tmax/T。width 缺省值为 0.5.%example2-7 triangular wavet=0:0.0001:5;A=1;T=1;w0=2*pi/T;yt=sawtooth(w0*t,1);plot(t,yt)9axis(0,5,-1.5,1.5)8、信号的尺度变换、翻转、平移将信号 x(t)写成函数,函数名为 x2_1,程序如下:function yt=x2_1(t)yt=0.5*(t+2).*(t=-2调用格式:f=Heaviside(t)2、已知 x(t)的波形如下图,(1)编写实现该信号波形的函数 x(t)的 MATLAB 函数;(2)画出 x(t),x(0.5t) 和 x(2-0.5t)。01235t)(tx