MATLAB中矩阵方程求解的实现

MATLAB中矩阵方程求解的实现

一、矩阵方程

1、定义：

 

 

2、分类

http://naotu.baidu.com/file/14d36860667d356a54490320cdab2950?token=56a299fe0e0815fe

二、M代码实现

1、M代码

  1 function d=CDBH_for_sov_JZFC(a,b)
  2 [m1,n1]=size(a);
  3 [m2,n2]=size(b);
  4 c=[a,b];
  5 ra=rank(a);             %矩阵a的秩
  6 rb=rank(b);             %矩阵b的秩
  7 rc=rank(c);             %矩阵[a,b]的秩
  8 zero=zeros(m2,n2);      %构造与b规格相同的零矩阵
  9 pj=zero~=b;             %确定b中非零元素的个数
 10 pj=sum(pj);
 11 pj=sum(abs(b));         %更新，把判据改为b中所有值的绝对值之和
 12 global i1;              %用于阶梯型的计算
 13 global i0;              %用于阶梯型的计算，其值为当前列于当前行的差值
 14 global cx;              %用于记录阶梯型的首个元素的位置
 15 i1=0;
 16 i0=0;
 17 x_fqc=[];                  %非齐次计算中，用于记录阶梯型的首个元素的位置
 18 x_qc=[];                  %在齐次计算中，用于记录阶梯型的首个元素的位置
 19 cx=[];
 20 if m1~=m2
 21     error('输入有误，无法计算');
 22     return;
 23 end
 24 switch pj
 25     %这种情况为其次方程组
 26     case 0
 27         switch rc
 28             %这种情况只有零解
 29             case n1
 30                 
 31                 d=zeros(n1,1);
 32                 %disp(d);
 33                 
 34                 %这种情况有基础解系
 35             case num2cell([0:n1-1])
 36                 %求解思路：
 37                 %一.矩阵变换
 38                 %1）化为行阶梯型
 39                 %2）化为标准型
 40                 %二.求基础解系
 41                 %1）分别取非线性向量
 42                 %2）则线性向量的值为，上述向量对应元素的相反数
 43                 %三.输出结果
 44                 
 45                 %一、矩阵变换
 46                 for i=1:m1-1;    %需要对行数-1行进行行变换
 47                     %选中了第i行
 48                     %如果都为非线性相关的向量，则阶梯的行列数相等，若存在线性相
 49                     %的向量，则需要构造一个i0来表示阶梯对应的列，并用i1表示列于
 50                     %行数的差值。
 51                     i0=i1+i;
 52                     %因为i1的值根据本行元素具体情况确定，因此需要注意，应当在这
 53                     %个for循环内完成下三角、化为1、上三角的所有操作。
 54                     %Step01 检查阶梯元素是否等于零
 55                     %       若等于零，则需要与第一个不为零的行对调，若全为零，
 56                     %       则i1+1，并再次循环。
 57                     pj_01=0;%初始化以下循环的判据
 58                     while pj_01==0 %利用判据pj_01来判断是否需要再次循环，在循环中，通过修改pj的值来跳出循环。
 59                         pj_01=1;    %没有特殊情况执行完跳出循环
 60                         i0=i+i1;
 61                         if c(i,i0)==0
 62                             %01 找到第i0列，i行及以下第一个非零元素的位置k
 63                             k=find(c([i+1:end],i0));
 64                             %k=k(1);
 65                             %02 将k行与i行对调（若k为空集，应当i1+1并再循环）
 66                             if k~=[];           %k不为空集时，将k行与i行对调
 67                                 k=k(1);
 68                                 e=