Basic Matrix Operations一、实验目的1、掌握向量和矩阵的创建方法；2、掌握向量和矩阵元素的索引方法；3、掌握向量和矩阵的基本操作；4、利用MATLAB编写程序进行矩阵运算。

二、基础知识1、常见数学函数函数名数学计算功能函数名数学计算功能Abs(x) 实数的绝对值或复数的幅值floor(x) 对x朝-∞方向取整Acos(x) 反余弦arcsin x gcd(m，n）求正整数m和n的最大公约数acosh(x) 反双曲余弦arccosh x imag(x) 求复数x的虚部angle(x) 在四象限内求复数 x 的相角lcm(m，n) 求正整数m和n的最小公倍数asin(x) 反正弦arcsin x log(x) 自然对数（以e为底数）asinh(x) 反双曲正弦arcsinh x log10(x) 常用对数（以10为底数）atan(x) 反正切arctan x real(x) 求复数x的实部atan2(x,y) 在四象限内求反正切Rem(m，n) 求正整数m和n的m/n之余数atanh(x) 反双曲正切arctanh x round(x) 对x四舍五入到最接近的整数ceil(x) 对x朝+∞方向取整sign(x) 符号函数：求出x的符号conj(x) 求复数x的共轭复数sin(x) 正弦sin xcos(x) 余弦cos x sinh(x) 反双曲正弦sinh xcosh(x) 双曲余弦cosh x sqrt(x) 求实数x的平方根：xexp(x) 指数函数xe tan(x) 正切tan xfix(x) 对x朝原点方向取整tanh(x) 双曲正切tanh x2、常量与变量系统的变量命名规则：变量名区分字母大小写；变量名必须以字母打头，其后可以是任意字母，数字，或下划线的组合。

此外，系统内部预先定义了几个有特殊意义和用途的变量，见下表：特殊的变量、常量取值ans 用于结果的缺省变量名pi 圆周率π的近似值（3.1416）eps 数学中无穷小（epsilon）的近似值（2.2204e - 016）inf 无穷大，如 1/0 = inf （infinity）NaN 非数，如 0/0 = NaN （Not a Number），inf / inf = NaNi，j 虚数单位：i = j =1−数值型向量（矩阵）的输入①任何矩阵（向量），可以直接按行方式……输入每个元素：同一行中的元素用逗号（，）或者用空格符来分隔；行与行之间用分号（；）分隔。

所有元素处于一方括号（[ ]）内；例1：>> Time = [11 12 1 2 3 4 5 6 7 8 9 10]>> X_Data = [2.32 3.43；4.37 5.98]②系统中提供了多个命令用于输入特殊的矩阵：函数功能函数功能compan 伴随阵toeplitz Toeplitz矩阵diag 对角阵vander Vandermonde矩阵hadamard Hadamard矩阵zeros 元素全为0的矩阵hankel Hankel矩阵ones 元素全为1的矩阵invhilb Hilbert矩阵的逆阵 rand 元素服从均匀分布的随机矩阵kron Kronercker张量积randn 元素服从正态分布的随机矩阵magic 魔方矩阵eye 对角线上元素为1的矩阵pascal Pascal矩阵meshgrid 由两个向量生成的矩阵上面函数的具体用法，可以用帮助命令help得到。

如：meshgrid(x,y)输入 x=[1 2 3 4]; y=[1 0 5]; [X,Y]=meshgrid(x, y)，则X = Y =1 2 3 4 1 1 1 11 2 3 4 0 0 0 01 2 3 4 5 5 5 5目的是将原始数据x，y转化为矩阵数据X，Y。

3、数组（矩阵）的点运算运算符：+（加）、-（减）、./（右除）、.\（左除）、.^（乘方），例3：>> g = [1 2 3 4]；h = [4 3 2 1]；>> s1 = g + h, s2 = g.*h, s3 = g.^h, s4 = g.^2, s5 = 2.^h4、矩阵的运算运算符：+（加）、-（减）、*（乘）、/（右除）、\（左除）、^（乘方）、’（转置）等；常用函数：det（行列式）、inv（逆矩阵）、rank（秩）、eig（特征值、特征向量）、rref（化矩阵为行最简形）例4：>> A= [2 0 –1; 1 3 2]; B= [1 7 –1; 4 2 3; 2 0 1];>> M = A*B % 矩阵A与B按矩阵运算相乘M =0 14 -317 13 10>> det_B = det(B) % 矩阵A的行列式det_B =20>> rank_A = rank(A) % 矩阵A的秩rank_A =2>> inv_B = inv(B) % 矩阵B的逆矩阵inv_B =0.1000 -0.3500 1.15000.1000 0.1500 -0.3500-0.2000 0.7000 -1.3000>> [V,D] = eig(B) % 矩阵B的特征值矩阵V与特征向量构成的矩阵D =7.2680 0 00 -1.6340 + 0.2861i 00 0 -1.6340 - 0.2861i>> X = A/B % A/B = A*B-1，即XB=A，求X0.4000 -1.4000 3.60000.0000 1.5000 -2.5000>> Y = B\A % B\A = B-1*A，即BY=A，求Y三、实验内容1、练习数据和符号的输入方式，将前面的命令在命令窗口中执行通过。

2、键入常数矩阵输入命令：a = [1 2 3] 与 a = [1;2;3]记录结果，比较显示结果有何不同；b = [1 2 5] 与 b = [1 2 5];记录结果，比较显示结果有何不同；a a’b b’记录结果，比较变量加“’”后的区别；c = a * b 1 2 52 4 103 6 15c = a* b’记录显示结果与出错原因；a = [1 2 3; 4 5 6; 7 8 0]，求a^2a^0.5。

30 36 1566 81 4239 54 690.5977 + 0.7678i 0.7519 + 0.0979i 0.5200 - 0.4680i1.4102 + 0.1013i 1.7741 + 0.6326i 1.2271 - 0.7467i1.2757 - 1.0289i 1.6049 - 1.0272i 1.1100 + 1.6175i3、使用冒号选出指定元素：已知 A=[1 2 3;4 5 6;7 8 9]，求 A 中第 3 列前 2 个元素，A 中所有列第 2，3 行的元素。

ns =36ns =4 5 67 8 04、输入 A=[7 1 5; 2 5 6; 3 1 5]，B=[1 1 1; 2 2 2; 3 3 3]，在命令窗口中执行下列表达式，掌握其含义：A(2,3)ans =6A(:,2)A(3,:) A(:,1:2:3)ans =7 52 63 5A(:,3).*B(:,2)ns =51215A(:,3)*B(2,:) A*B A.*B A^2 A.^2 B/A B./A5、建立M 文件，求123221343A⎛⎫⎪= ⎪⎪⎝⎭的逆矩阵。

nv_M =1.0000 3.0000 -2.0000 -1.5000 -3.0000 2.5000 1.0000 1.0000 -1.00006、设123221343A⎛⎫⎪= ⎪⎪⎝⎭,2153B⎛⎫= ⎪⎝⎭,132031C⎛⎫⎪= ⎪⎪⎝⎭，建立M文件，求矩阵X ，使满足：AXB = C。

四，详细设计（1）存储要点position[col]=position[col-1]+num[col-1];三元组表(row，col，v)稀疏矩阵（（行数m，列数n，非零元素个数t），三元组，...，三元组）row col vA->data[0] 0 4 81 1 0 12 1 2 33 2 1 24 3 0 6………max-1（2）乘法运算要点已知稀疏矩阵A(m1× n1)和B(m2× n2)，求乘积C(m1× n2)。

稀疏矩阵A、B、C及它们对应的三元组表A.data、B.data、C.data如图6所示。

由矩阵乘法规则知：C(i，j)=A(i,1)×B(1,j)+A(i,2)×B(2,j)+…+A(i,n)×B(n,j)=这就是说只有A(i,k)与B(k,p)(即A元素的列与B元素的行相等的两项)才有相乘的机会，且当两项都不为零时，乘积中的这一项才不为零。

矩阵用二维数组表示时，a11只有可能和B中第1行的非零元素相乘，a12只有可能和B中第2行的非零元素相乘，…，而同一行的非零元是相邻存放的，所以求c11和c12同时进行：求a11*b11累加到c11，求a11*b12累加到c12，再求a12*b21累加到c11，再求a12*b22累加到c22.，…，当然只有aik和bkj(列号与行号相等)且均不为零(三元组存在)时才相乘，并且累加到cij当中去。

（3）稀疏矩阵的快速转置要点：矩阵A中三元组的存放顺序是先行后列，对同一行来说，必定先遇到列号小的元素，这样只需扫描一遍A.data 。

所以需引入两个向量来实现：num[n+1]和position [n+1]，num[col]表示矩阵A中第col列的非零元素的个数(为了方便均从1单元用起)，position [col]初始值表示矩阵A中的第col列的第一个非零元素在B.data中的位置。

于是position的初始值为：position [1]=1;position [col]= position [col-1]+num[col-1]; 2≤col≤n依次扫描A.data，当扫描到一个col列元素时，直接将其存放在B.data的position [col]位置上，position [col]加1，position [col]中始终是下一个col列元素在B.data中的位置。

A->row[0] 01 12 33 4（4）逆矩阵⒈判断矩阵是否为方阵⒉逆矩阵的算法：①求行列式的值②求矩阵的伴随矩阵③用伴随矩阵除以行列式⒊求逆矩阵的流程：五、源程序（测试结果）#include<iostream>#include<string>#include<iomanip>#include<cmath>using namespace std;#define max 100#define datatype int typedef struct{int row,col;//行，列datatype v;//非0数值}Node;typedef struct{Node data[max];//稀疏矩阵int m,n,t;//m行，n列，t非0数个数}Matrix;/*求逆矩阵存储*/typedef struct{ //存储结构int m, n; //行、列数double *p; //矩阵基址}nMatrix;void In(nMatrix a){ //求逆输入cout<<"请将矩阵a的行、列数再次输入：";cin >>a.m>>a.n;int m = a.m, n = a.n;int i, j;double *p = a.p = new double[m * n]; //p是行指针cout<<"请按行优先输入矩阵a的全部数值：\n";for (i = 0; i < m; p += n, i++){for (j = 0; j < n; j++)cin>>p[j]; //即a.p[i*n+j]}}void Out(nMatrix a){ //求逆输出int m = a.m, n = a.n;int i, j;double *p = a.p;for (i = 0; i < m; p += n, i++){for (j = 0; j < n; j++)cout<<p[j]<<'\t';cout<<endl;}}istream& operator >>(istream& input,Matrix &A){ int i;cout<<"请输入行数:";input>>A.m;cout<<"请输入列数:";input>>A.n;cout<<"请输入非0值个数:";input>>A.t;for(i=1;i<=A.t;i++){cout<<"请输入行数，列数，非0值:"<<"("<<i<<")\n";input>>A.data[i].row>>A.data[i].col>>A.data[i].v;}return input;}ostream& operator <<(ostream& output,Matrix &A){int i,j,t=1,k=0;for(i=1;i<=A.m;i++){for(j=1;j<=A.n;j++){if(A.data[t].row==i&&A.data[t].col==j){output<<A.data[t].v<<'\t';t++;}elseoutput<<k<<'\t';}cout<<"\n";}return output;}Matrix operator +(Matrix A,Matrix B){//加法int i,j,k;Matrix C;if(A.m!=B.m||A.n!=B.n){cout<<"这两个矩阵不能相加"<<endl;exit(0);}C.m=A.m;C.n=A.n;C.t=0;if(A.t==0&&B.t==0) exit(0);i=j=k=1;while(i<=A.t&&j<=B.t){if(A.data[i].row<B.data[j].row){C.data[k]=A.data[i];i++;k++;}else{if(A.data[i].row>B.data[j].row){C.data[k]=B.data[j];j++;k++;}else{if(A.data[i].col<B.data[j].col){C.data[k]=A.data[i];i++;k++;}else{if(A.data[i].col>B.data[j].col){C.data[k]=B.data[j];j++;k++;}else{if(A.data[i].v+B.data[j].v!=0){C.data[k].row=A.data[i].row;C.data[k].col=A.data[i].col;C.data[k].v=A.data[i].v+B.data[j].v;k++;}i++;j++;}}}}}while(i<A.t){C.data[k]=A.data[i];i++;k++;}while(j<B.t){C.data[k]=B.data[j];j++;k++;}C.t=k;return C;}Matrix operator -(Matrix A,Matrix B){ Matrix C;int i,j,k;if(A.m!=B.m||A.n!=B.n){cout<<"这两个矩阵不能相减";exit(0);} C.m=A.m;C.n=A.n;C.t=0;if(A.t==0&&B.t==0) exit(0);i=j=k=1;while(i<=A.t&&j<=B.t){ if(A.data[i].row<B.data[j].row){ C.data[k]=A.data[i];i++;k++;}else{ if(A.data[i].row>B.data[j].row){ C.data[k]=B.data[j];j++;k++;}else{ if(A.data[i].col<B.data[j].col) { C.data[k]=A.data[i];i++;k++;}else{ if(A.data[i].col>B.data[j].col) { C.data[k]=B.data[j];j++;k++;}else{ if(A.data[i].v-B.data[j].v!=0){ C.data[k].row=A.data[i].row;C.data[k].col=A.data[i].col;C.data[k].v=A.data[i].v-B.data[j].v;k++;}i++;j++;}}}}}while(i<A.t){ C.data[k]=A.data[i];i++;k++;}while(j<B.t){ C.data[k]=B.data[j];j++;k++;}C.t=k;return C;}Matrix operator *(Matrix A,Matrix B){Matrix C;int k,p,crow,brow,q,ccol;int num[max],pos[max],ctemp[max];if (A.n==B.m){for(k=1;k<=B.m;k++)num[k]=0;for(k=1;k<=B.t;k++)num[B.data[k].row]++;pos[1]=1;for(k=2;k<=B.t;k++)pos[k]=pos[k-1]+num[k-1];pos[1+B.t]=pos[B.t]+1;C.m=A.m; C.n=B.n; C.t=0; p=1;while(p<=A.t){crow=A.data[p].row;for(k=1;k<=C.n;k++)ctemp[k]=0;while (p<=A.t&&A.data[p].row==crow){brow=A.data[p].col;for(q=pos[brow];q<=pos[brow+1]-1;q++){ccol=B.data[q].col;ctemp[ccol]=ctemp[ccol]+A.data[p].v*B.data[q].v;}p=p+1;}for(ccol=1;ccol<=B.n;ccol++)if(ctemp[ccol]!=0){C.t=C.t+1;C.data[C.t].row=crow;C.data[C.t].col=ccol;C.data[C.t].v=ctemp[ccol];}}}elsecout<<"这两个矩阵不能相乘";return C;}Matrix operator ~(Matrix A){Matrix B;int col,i,p,q;int num[max],position[max];B.t=A.t;B.m=A.n; B.n=A.m;if(B.t){for(col=1;col<=A.n;col++)num[col]=0;for(i=1;i<=A.t;i++)num[A.data[i].col]++;position[1]=1;for(col=2;col<=A.n;col++)position[col]=position[col-1]+num[col-1];for(p=1;p<=A.t;p++){col=A.data[p].col;q=position[col];B.data[q].row=A.data[p].col;B.data[q].col=A.data[p].row;B.data[q].v=A.data[p].v;position[col]++;}}return B;}nMatrix Trs(nMatrix a){ //求逆矩阵的先转置nMatrix trs;trs.m = a.n;trs.n = a.m;trs.p = new double[a.m * a.n];for (int i = 0; i < a.m; i++){for (int j = 0; j < a.n; j++){trs.p[j * a.m + i] = a.p[i * a.n + j];}}return trs;}nMatrix Adjunct(nMatrix a, int indexm, int indexn){ //求第indexm行indexn列元素的代数余子式nMatrix adj;adj.m=a.m - 1;adj.n=a.n - 1;adj.p = new double[(a.n - 1) * (a.n - 1)];for (int i = 0; i < indexm; i++){for (int j = 0; j < indexn; j++){adj.p[i * (a.n - 1) + j] = a.p[i * a.n + j];}for (int k = indexn + 1; k < a.n; k++){adj.p[i *(a.n - 1) + k -1] = a.p[i * a.n + k];}}for (int m = indexm + 1; m < a.n; m++){for (int j = 0; j < a.n - 1; j++){adj.p[(m - 1) * (a.n - 1) + j] = a.p[m * a.n + j];}for (int k = indexn + 1; k < a.n; k++){adj.p[(m - 1) * (a.n - 1) + k - 1] = a.p[m * a.n + k];}}return adj;}double Det(nMatrix a) //递归求行列式{double det = 0;if (a.m != a.n){cout<<"不是方阵，没有行列式!"<<endl;cout<<"求行列式退出"<<endl;}if (a.n == 1){det = a.p[0];return det;}else{for (int i = 0; i < a.n; i++){if (i % 2 == 0)det += a.p[i * a.n] * Det(Adjunct(a, i, 0));elsedet -= a.p[i * a.n] * Det(Adjunct(a, i, 0));}}return det;}nMatrix operator !(nMatrix a) { //求矩阵的逆nMatrix temp;temp.m=a.n;temp.n=a.m;temp.p = new double[a.m * a.n];//矩阵的逆 = 伴随矩阵 / 行列式double det = Det(a);if (det == 0) //如果行列式的值为0，则没有逆{cout<<"此矩阵没有逆!"<<endl;cout<<"求矩阵逆退出!"<<endl;}for (int i = 0; i < temp.m; i++){for (int j = 0; j < temp.n; j++){if ((i +j) % 2 == 0)temp.p[i * temp.m + j] = Det(Adjunct(a, i, j)) / det;elsetemp.p[i * temp.m + j] = -Det(Adjunct(a, i, j)) / det;}}return Trs(temp);}void main(){Matrix C={0},A={0},B={0};nMatrix a={0}, b={0},temp={0};cin>>A; cout<<"a=\n"<<A;cin>>B; cout<<"b=\n"<<B;C=A+B;cout<<"a+b=\n"<<C;C=A-B;cout<<"a-b=\n"<<C;C=A*B;cout<<"a*b=\n"<<C;C=~A;cout<<"~a=\n"<<C;In(a);cout<<"矩阵的行列式为:"<<endl;cout<<Det(a)<<endl;cout<<"矩阵的逆为:"<<endl;temp=!a;cout<<"!a=\n"<<Trs<<endl;六、总结这次数据结构课程设计的制作使我对matlab的理解更加深刻，也使我认识到了自己很多不足之处。