矩阵及线型方程组.ppt

Matlab 与计算方法,华侨大学计算机学院张国亮 2010.9,第3章 矩阵及线型方程组,3.1.1 数值矩阵的生成，遵循下列几个基本步骤：用空格或者逗号来区分一行里不同的元素。用分号；来区分不同的行。用方括号来括住全体元素。,3.1 MATLAB矩阵表示和运算,MATLAB，即“矩阵实验室”，它是以矩阵为基本运算单元。,第3章 矩阵及线型方程组,3.1.1 数值矩阵的生成：例：Time=11 12 1 2 3 4 5 6 7 8 9 10 Time=11 12 1 2 3 4 5 6 7 8 9 10 X_Data=2.32 3.43;4.37 5.98X_Data=2.43 3.43 4.37 5.98 vect_a=1 2 3 4 5 vect_a=1 2 3 4 5 Null_M=%生成一个空矩阵,3.1 MATLAB矩阵表示和运算,第3章 矩阵及线型方程组,3.1.2 复数矩阵输入 两种方式：1 直接在数组元素中加入复数单位，a+i*b形式；2 输入数据矩阵，按照矩阵相乘的方式来表示矩阵第一种方式 a=2.7;b=13/25;C=1,2*a+i*b,b*sqrt(a);sin(pi/4),a+5*b,3.5+1C=1 27/5+13/25i 317/371 985/1393 53/10 9/2,3.1 MATLAB矩阵表示和运算,第3章 矩阵及线型方程组,3.1.2 复数矩阵输入 第2种方式 R=1 2 3;4 5 6,M=11 12 13;14 15 16R=1 2 3 4 5 6M=11 12 13 14 15 16 CN=R+i*MCN=1.0000+11.0000i 2.0000+12.0000i 3.0000+13.0000i 4.0000+14.0000i 5.0000+15.0000i 6.0000+16.0000i,3.1 MATLAB矩阵表示和运算,第3章 矩阵及线型方程组,3.1.3 符号矩阵的生成 在MATLAB中输入符号矩阵的方法和输入数值类型的矩阵在形式上很相像，只不过要用到符号矩阵定义函数sym，或者是用到符号定义函数syms，先定义一些必要的符号变量，再像定义普通矩阵一样输入符号矩阵。1用命令sym定义矩阵2.用命令syms定义矩阵,3.1 MATLAB矩阵表示和运算,第3章 矩阵及线型方程组,3.1.3 符号矩阵的生成 1用命令sym定义矩阵：这时的函数sym实际是在定义一个符号表达式，这时的符号矩阵中的元素可以是任何的符号或者是表达式，而且长度没有限制，只是将方括号置于用于创建符号表达式的单引号中。如下例：sym(a b c；Jack，Help Me!，NO WAY!，)sym_matrix=a b cJack Help Me!NO WAY!sym_digits=sym（1 2 3；a b c；sin（x）cos（y）tan（z）sym_digits=1 2 3a b csin（x）cos（y）tan（z）,3.1 MATLAB矩阵表示和运算,第3章 矩阵及线型方程组,3.1.3 符号矩阵的生成 2 用命令syms定义矩阵 先定义矩阵中的每一个元素为一个符号变量，而后像普通矩阵一样输入符号矩阵。例：syms a b c M1=sym(Classical);M2=sym(Jazz);M3=sym(Blues)syms_matrix=a b c;M1,M2,M3syms_matrix=a b cClassical Jazz Blues,3.1 MATLAB矩阵表示和运算,第3章 矩阵及线型方程组,包含全零阵、单位阵、全一阵、随机阵命令 全零阵,函数 zeros格式 B=zeros(n)%生成nn全零阵B=zeros(m,n)%生成mn全零阵B=zeros(m n)%生成mn全零阵B=zeros(d1,d2,d3)%生成d1d2d3全零阵或数组B=zeros(d1 d2 d3)%生成d1d2d3全零阵或数组B=zeros(size(A)%生成与矩阵A相同大小的全零阵命令 单位阵,函数 eye格式 Y=eye(n)%生成nn单位阵Y=eye(m,n)%生成mn单位阵Y=eye(size(A)%生成与矩阵A相同大小的单位阵,3.2 特殊矩阵的生成,第3章 矩阵及线型方程组,命令 全1阵,函数 ones 格式 Y=ones(n)%生成nn全1阵Y=ones(m,n)%生成mn全1阵Y=ones(m n)%生成mn全1阵Y=ones(d1,d2,d3)%生成d1d2d3全1阵或数组Y=ones(d1 d2 d3)%生成d1d2d3全1阵或数组Y=ones(size(A)%生成与矩阵A相同大小的全1阵,3.2 特殊矩阵的生成,第3章 矩阵及线型方程组,命令 均匀分布随机矩阵,函数 rand 格式 Y=rand(n)%生成nn随机矩阵，其元素在（0，1）内Y=rand(m,n)%生成mn随机矩阵Y=rand(m n)%生成mn随机矩阵Y=rand(m,n,p,)%生成mnp随机矩阵或数组Y=rand(m n p)%生成mnp随机矩阵或数组Y=rand(size(A)%生成与矩阵A相同大小的随机矩阵例：产生一个34随机矩阵 R=rand(3,4)R=0.9501 0.4860 0.4565 0.4447 0.2311 0.8913 0.0185 0.6154 0.6068 0.7621 0.8214 0.7919随机矩阵每次产生的数据是不同的。,3.2 特殊矩阵的生成,第3章 矩阵及线型方程组,3.3.1 加、减运算 运算符：“”和“”。运算规则：对应元素相加、减，即按线性代数中矩阵的“十”，“一”运算进行。A=1,1,1;1,2,3;1,3,6B=8,1,6;3,5,7;4,9,2A+BA-B 结果显示：A+B=9 2 74 7 105 12 8A-B=-7 0-5-2-3-4-3-6 4,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.2 矩阵相乘 运算符：*运算规则：按线性代数中矩阵乘法运算进行，即放在前面的矩阵的各行元素，分别与放在后面的矩阵的各列元素对应相乘并相加。1两个矩阵相乘X=2 3 4 5;1 2 2 1Y=0 1 1;1 1 0;0 0 1;1 0 0Z=X*Y结果显示为：Z=8 5 6 3 3 3,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.2 矩阵相乘 运算符：.*2 矩阵的数乘 例：a=2*X，则显示：a=4 6 8 10 2 4 4 2向量的点乘（内积）：维数相同的两个向量的点乘。数组乘法：A.*B表示A与B对应元素相乘。,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.2 矩阵相乘 3 矩阵(向量)的点乘（内积）函数 dot 格式 C=dot(A,B)%若A、B为向量，则返回向量A与B的点积，A与%B长度相同；若为矩阵，则A与B有相同的维数。C=dot(A,B,dim)%在dim维数中给出A与B的点积例:X=-1 0 2;Y=-2-1 1;Z=dot(X,Y)则显示：Z=4,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.2 矩阵相乘4 向量叉乘函数 cross格式 C=cross(A,B)%若A、B为向量，则返回A与B的叉 乘，即C=AB；矩阵同理。C=cross(A,B,dim)%在dim维数中给出向量A与B的叉积。A和B必须具有相同的维数，size(A,dim)和 size(B,dim)必须是3。例 计算垂直于向量(1,2,3)和(4,5,6)的向量。a=1 2 3;b=4 5 6;c=cross(a,b)结果显示：c=-3 6-3,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.2 矩阵相乘5 混合积混合积由cross、dot函数实现：例：计算向量a=(1,2,3)、b=(4,5,6)和c=(-3,6,-3)的混合积解：a=1 2 3;b=4 5 6;c=-3 6-3;x=dot(a,cross(b,c)结果显示：x=54注意：先叉乘后点乘，顺序不可颠倒。,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.2 矩阵相乘6 矩阵的卷积和多项式乘法函数 conv格式 w=conv(u,v)%u、v为向量，其长度可不相同。说明:长度为m的向量序列u和长度为n的向量序列v的卷积定义为：w(1)=u(1)*v(1)w(3)=u(1)*v(3)+u(2)*v(2)+u(3)*v(1)w(n)=u(1)*v(n)+u(2)*v(n-1)+u(n)*v(1),3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.2 矩阵相乘6 矩阵的卷积和多项式乘法例:展开多项式w=conv(1,2,2,conv(1,4,1,1)w=1 7 16 18 8P=poly2str(w,s)%将w表示成多项式P=s4+7 s3+16 s2+18 s+8注意：没有的项要用“0”来表示。,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.3 除法运算Matlab提供了两种除法运算：左除（）和右除（/）。一般情况下，x=ab是方程a*x=b的解，而x=b/a是方程x*a=b的解。例：a=1 2 3;4 2 6;7 4 9b=4;1;2;x=ab（区分x=b/a)则显示：x=-1.5000 2.0000 0.5000,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.3 除法运算 如果矩阵a为非奇异矩阵，则ab和b/a可通过a的逆矩阵与b阵得到：ab=inv(a)*b b/a=b*inv(a)注意:解方程组的时候inv()函数应用很少，因为很多情况下无法满足非奇异的条件。数组除法：A./B表示A中元素与B中元素对应相除。,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.4 矩阵乘方 运算符：运算规则：（1）当A为方阵，P为大于0的整数时，AP表示A的P次方，即A自乘P次；P为小于0的整数时，AP表示A的-|P|次方。（2）当A为方阵，p为非整数时，则,其中V为A的特征向量，为特征值对角矩阵。如果有重根，以上指令不成立。其中V为A的特征向量，（3）标量的数组乘方P.A，标量的数组乘方定义为 数组乘方：A.P：表示A的每个元素的P次乘方。,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.5 重要矩阵函数方阵指数，函数 expm矩阵的对数，函数 logm格式 Y=expm(A)Y=logm(X)%计算矩阵X的对数，它是expm(X)的反函数。例 A=1 1 0;0 0 2;0 0-1;Y=expm(A)Y=2.7183 1.7183 1.0862 0 1.0000 1.2642 0 0 0.3679 A=logm(Y)A=1.0000 1.0000 0.0000 0 0 2.0000 0 0-1.0000,3.3 矩阵运算,第3章 矩阵及线型方程组,方阵的函数函数 funm格式 F=funm(A,fun)%A为方阵，计算由fun指定的A的矩阵函数，fun可以是任意基本函数，如sin、cos等等，例如：funm(A,exp)=expm(A)与前面指数运算的方式是一样的。,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.7 矩阵转置运算符：运算规则：若矩阵A的元素为实数，则与线性代数中矩阵的转置相同。若A为复数矩阵，则A转置后的元素由A对应元素的共轭复数构成。,3.3 矩阵运算,第3章 矩阵及线型方程组,方阵的行列式函数 det格式 d=det(X)例：A=1 2 3;4 5 6;7 8 9A=1 2 3 4 5 6 7 8 9 D=det(A)D=0,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.9 逆与伪逆命令 逆,函数 inv格式 Y=inv(X)%求方阵X的逆矩阵。若X为奇异阵或 近似奇异阵，将给出警告信息。例:求的逆矩阵A=1 2 3;2 2 1;3 4 3;Y=inv(A)则结果显示为 Y=1.0000 3.0000-2.0000-1.5000-3.0000 2.5000 1.0000 1.0000-1.0000,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.9 逆与伪逆命令 伪逆,函数 pinv格式 B=pinv(A)%求矩阵A的伪逆 B=pinv(A,tol)%tol为误差：max(size(A)*norm(A)*eps说明:当矩阵为长方阵时，方程AX=I和XA=I至少有一个无解，这时A的伪逆能在某种程度上代表矩阵的逆，若A为非奇异矩阵，则pinv(A)=inv(A)。A=1 1 1;2 2 2;1 2 3，inv(A),pinv(A)关于求伪逆的方法有很多钟，感兴趣的可以参考相关资料，求伪逆是矩阵以及方程组中十分重要的内容。,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.10 矩阵的迹、秩、范数、条件数 矩阵的迹 函数 trace格式 b=trace(A)%返回矩阵A的迹，即A的对角线元素之和。矩阵的秩函数 rank格式 k=rank(A)%求矩阵A的秩,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.10 矩阵的迹、秩、范数、条件数 矩阵和向量的范数函数 norm矩阵的条件数函数 cond格式 c=cond(X)%X的最大奇异值和最小奇异值的商。说明 线性方程组AX=b的条件数是一个大于或者等于1的实数，用来衡量关于数据 中的扰动，也就是A/或b对解X的灵敏度，一个差条件的方程组的条件数很大。,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.11 符号矩阵运算 关于符号计算的问题在第四章来详细讲，这里记住Matlab 把符号矩阵的四则运算简化为与数值矩阵完全相同的运算方式，其运算符为：加（）、减（）、乘（）、除（/、）等。符号矩阵的其他一些基本运算包括转置（）、行列式（det）、逆（inv）、秩（rank）、幂（）和指数（exp和expm）等都与数值矩阵相同。例：A=sym(1/x,1/(x+1);1/(x+2),1/(x+3);B=sym(x,1;x+2,0);C=B-A则显示：C=x-1/x 1-1/(x+1)x+2-1/(x+2)-1/(x+3),3.2 矩阵运算,第3章 矩阵及线型方程组,3.3.12 矩阵关系运算M AT L A B有用于比较矩阵的六个关系运算符，也可以对矩阵与一个标量进行比较，即矩阵中的每个元素与标量进行比较。关系运算符如下：大于=大于等于=等于=不等于 在一个表达式中，算术运算符优先级最高，其次是关系运算符，最低级别是逻辑运算符，圆括号可以改变其顺序。,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.13 矩阵逻辑运算设矩阵A和B都是mn矩阵或其中之一为标量，在MATLAB中定义了如下的逻辑运算：（1）矩阵的与运算格式 A&B或and(A,B)说明 A与B对应元素进行与运算，若两个数均非0，则结果元素的值为1，否则为0。（2）或运算格式 A|B或or(A,B)说明 A与B对应元素进行或运算，若两个数均为0，则结果元素的值为0，否则为1。（3）非运算格式 A或not(A)说明 若A的元素为0，则结果元素为1，否则为0。（4）异或运算格式 xor(A,B)说明 A与B对应元素进行异或运算，若相应的两个数中一个为0，一个非0，则结果为0，否则为1。,3.3 矩阵运算,第3章 矩阵及线型方程组,3.3.12 矩阵逻辑运算例：A=0 2 3 4;1 3 5 0,B=1 0 5 3;1 5 0 5A=0 2 3 4 1 3 5 0B=1 0 5 3 1 5 0 5,3.3矩阵运算,C1=A&B,C2=A|B,C3=A,C4=xor(A,B)C1=0 0 1 1 1 1 0 0C2=1 1 1 1 1 1 1 1C3=1 0 0 0 0 0 0 1C4=1 1 0 0 0 0 1 1,第3章 矩阵及线型方程组,LU分解 矩阵的三角分解又称LU分解，它的目的是将一个矩阵分解成一个下三角矩阵L和一个上三角矩阵U的乘积，即A=LU。函数 lu（）格式 L,U=lu(X)%U为上三角阵，L为下三角阵 或其变换形式，满足LU=X。L,U,P=lu(X)%U为上三角阵，L为下三角 阵，P为单位矩阵的行变换矩阵，满 足LU=PX。,3.4 矩阵分解,第3章 矩阵及线型方程组,LU分解例：A=1 2 3;4 5 6;7 8 9L,U=lu(A)L=0.1429 1.0000 0 0.5714 0.5000 1.0000 1.0000 0 0U=7.0000 8.0000 9.0000 0 0.8571 1.7143 0 0 0.0000,3.4 矩阵分解,第3章 矩阵及线型方程组,Cholesky分解函数 chol格式 R=chol(X)%如果X为n阶对称正定矩阵，则存在一个实的非奇异上三角阵R，满足R*R=X；若X非正定，则产生错误信息。R,p=chol(X)%不产生任何错误信息，若X为正定阵，则p=0，R与上相同；若X非正定，则p为正整数，R是有序的上三角阵。,3.4 矩阵分解,第3章 矩阵及线型方程组,Cholesky分解X=pascal(4)%产生4阶pascal矩阵X=1 1 1 1 1 2 3 4 1 3 6 10 1 4 10 20 R,p=chol(X)R=1 1 1 1 0 1 2 3 0 0 1 3 0 0 0 1p=0,3.4 矩阵分解,第3章 矩阵及线型方程组,3.4.3 特征值分解*函数 eig格式 d=eig(A)%求矩阵A的特征值d，以向量形式存放d。d=eig(A,B)%A、B为方阵，求广义特征值d，以向量形式存放d。V,D=eig(A)%计算A的特征值对角阵D和特征向量V，使 AV=VD成立。V,D=eig(A,B)%计算广义特征值向量阵V和广义特征值阵D，满足AV=BVD。说明 一般特征值问题是求解方程：解的问题。广义特征值问题是求方程：解的问题。,3.4 矩阵分解,第3章 矩阵及线型方程组,3.4.4 奇异值分解函数 svd格式 s=svd(X)%返回矩阵X的奇异值向量U,S,V=svd(X)%返回一个与X同大小的对角矩阵S，两个酉矩阵U和V，且满足A=U*S*V。若A为mn阵，则U为 mm阵，V为nn阵。奇异值在S的对角线上，非负且按降序排列。U,S,V=svd(X,0)%得到一个“有效大小”的分解，只计算出矩阵U的前n列，矩阵S的大小为nn。,3.4 矩阵分解,第3章 矩阵及线型方程组,3.4.4 奇异值分解例:A=1 2;3 4;5 6;7 8;U,S,V=svd(A)U=-0.1525-0.8226-0.3945-0.3800-0.3499-0.4214 0.2428 0.8007-0.5474-0.0201 0.6979-0.4614-0.7448 0.3812-0.5462 0.0407,3.4 矩阵分解,S=14.2691 0 0 0.6268 0 0 0 0V=-0.6414 0.7672-0.7672-0.6414,第3章 矩阵及线型方程组,线性方程的求解分为两类：一类是方程组求唯一解或求特解，另一类是方程组求无穷解即通解。两类问题可以通过系数矩阵的秩来判断：若系数矩阵的秩r=n（n为方程组中未知变量的个数），则有唯一解；若系数矩阵的秩rn，则可能有无穷解；线性方程组的无穷解=对应齐次方程组的通解+非齐次方程组的一个特解；其中，特解的求法属于解的第一类问题，通解部分属第二类问题。,3.5 线性方程的组的求解,第3章 矩阵及线型方程组,3.5.1 求线性方程组的唯一解或特解（第一类问题）这类问题的求法分为两类：一类主要用于解低阶稠密矩阵 直接法；另一类是解大型稀疏矩阵 迭代法。1利用矩阵除法求线性方程组的特解（或一个解）方程：AX=b解法：X=Ab,3.5 线性方程的组的求解,第3章 矩阵及线型方程组,3.5.1 求线性方程组的唯一解或特解（第一类问题）求方程组 的解解：A=5 6 0 0 0；1 5 6 0 0；0 1 5 6 0；0 0 1 5 6；0 0 0 1 5;B=1 0 0 0 1;R_A=rank(A)%求秩X=AB%求解运行后结果如下R_A=5X=2.2662-1.7218 1.0571-0.5940 0.3188这就是方程组的解。,3.5 线性方程的组的求解,第3章 矩阵及线型方程组,3.5.1 求线性方程组的唯一解或特解（第一类问题）1利用矩阵除法求线性方程组的特解例 求方程组 的一个特解。解：A=1 1-3-1;3-1-3 4;1 5-9-8;B=1 4 0;X=AB%由于系数矩阵不满秩，该解法可能存在误差。X=0 0-0.5333 0.6000（一个特解近似值）。,3.5 线性方程的组的求解,第3章 矩阵及线型方程组,3.5.1 求线性方程组的唯一解或特解（第一类问题）2 利用矩阵的LU、QR和cholesky分解求方程组的解函数（1）LU分解：LU分解又称Gauss消去分解，可把任意方阵分解为下三角矩阵的基本变换形式（行交换）和上三角矩阵的乘积。即A=LU，L为下三角阵，U为上三角阵。则：A*X=b 变成L*U*X=b 所以X=U(Lb)这样可以大大提高运算速度。说明 这三种分解，在求解大型方程组时很有用。其优点是运算速度快、可以节省磁盘空间、节省内存。,3.5 线性方程的组的求解,第3章 矩阵及线型方程组,3.5.2 求线性齐次方程组的通解（第二类问题）在Matlab中，函数null用来求解零空间，即满足AX=0的解空间，实际上是求出解空间的一组基（基础解系）格式 z=null%z的列向量是方程AX=0的有理基例：求解方程组 的通解,3.5 线性方程的组的求解,A=1 2 2 1;2 1-2-2;1-1-4-3;format rat%指定有理式格式输出B=null(A,r)%求解空间的有理基运行后显示结果如下：B=2 5/3-2-4/3 1 0 0 1,第3章 矩阵及线型方程组,非齐次线性方程组需要先判断方程组是否有解，若有解，再去求通解。因此，步骤为：第一步：判断AX=b是否有解，若有解则进行第二步第二步：求AX=b的一个特解第三步：求AX=0的通解第四步：AX=b的通解=AX=0的通解+AX=b的一个特解。,3.5 求非齐次线性方程组的通解,第3章 矩阵及线型方程组,例 求解方程组 解：在Matlab中建立M文件如下：,3.5 求非齐次线性方程组的通解,A=1-2 3-1;3-1 5-3;2 1 2-2;b=1 2 3;B=A b;n=4;R_A=rank(A);R_B=rank(B)format ratif R_A=R_BR_B=3 X=equition no solve,第3章 矩阵及线型方程组,例 求解方程组：,3.5 求非齐次线性方程组的通解,A=1 1/2 1/2-1;1 1-1 1;1-1/4-1 1;-8-1 1 1b=0 10 0 1 B=A b;n=4;R_A=rank(A);R_B=rank(B)format ratif R_A=R_B&R_A=n%判断有唯一解 X=Abelseif R_A=R_B&R_An%判断有无穷解 X=Ab%求特解 C=null(A,r)%求AX=0的基础解系else X=equition no solve%判断无解end运行后结果显示：R_A=R_B=4,第3章 矩阵及线型方程组,线性方程组的LQ解法 函数 symmlq 双共轭梯度法解方程组 bicg 复共轭梯度平方法解方程组 cgs广义最小残差法 qmres 最小残差法解方程组 minres.,3.6 其它求解线性方程的组的方法,本章小结,MATLAB矩阵表示和运算方法。一些重要的矩阵运算方法和命令。线型方程组的一般解法.,