【问题标题】:Factoring matrices into elementary matrices将矩阵分解为基本矩阵
【发布时间】:2010-08-03 17:53:39
【问题描述】:

在 MATLAB、Maple 或 Mathematica 中是否有一个包可以做到这一点?

【问题讨论】:

    标签: matlab matrix wolfram-mathematica matrix-multiplication maple


    【解决方案1】:

    我想你所说的“基本”矩阵是指那些进行行交换、行乘法和行加法等基本操作的矩阵。

    您可能有兴趣知道这是 PLU 分解(因式分解)结果的一部分。 PLU分解得到的U是高斯消元的结果,而那个PLU分解只是GE的变相。 PLU 分解的 P 和 L 对完成 GE 所采取的基本操作进行编码。而且 Maple、Matlab 和 Mathematica 都有很好的 PLU 分解例程。所以你可以得到基本的因素。

    现在假设我们不需要进行任何行交换​​。因此,给定矩阵 M,我们可以得到下三角 ​​L 和上三角 M。位于主对角线下方的 L 的条目是构造基本行加法矩阵的值。

    最后是 Maple 代码来展示它是如何完成的。那里产生了三组基本矩阵。表 T1 中的第一组是由于将 M 变为行梯形的 GE 步骤,并且来自使用 l1 M 的 PLU 分解的 L。这就是完成的下三角形。接下来我们将M的PLU分解的U转置u1,从而处理M的上三角。

    表 T2 中的第二组基本行加法矩阵是由于得到 u1^%T(从 M 的 PLU 分解中的 U 的转置)到行梯形形式的 GE 步骤。它们是使用 u1^%T 的 PLU 分解的 L l2 中的条目构造的。

    只剩下u2 u1^%T 的 PLU 分解的 U。它是一个对角矩阵(如果没有执行行交换)。所以我们为u2的每一行构造基本的行缩放矩阵。

    最后,我们可以把它们按正确的顺序排列,然后将它们相乘。请注意,T2 矩阵以相反的顺序出现,即转置,因为它们必须相乘以形成 u1^%T。同样,T3 出现在 T1 和 T2 集之间,因为 T3 构造 u2

    作为以后的编辑,这里是作为 Maple 程序。现在它从排列结果生成行交换矩阵。而且它不会返回一些不必要的因素,这些因素恰好只是身份。

    请注意,这是针对精确矩阵,而不是浮点矩阵(您的里程可能会有所不同,这取决于如何按幅度选择枢轴以及它如何进行比较)。

    ElemDecomp:=proc(M::Matrix(square))
    local p1,u1,i,j,T1,T2,T3,p2,m,n,lu1,lu2,P1,P2;
    uses LinearAlgebra;
      (m,n):=Dimensions(M);
      p1,lu1:=LUDecomposition(M,output=[':-NAG']);
      for i from 1 to m-1 do
        for j from 1 to i do
          if lu1[i+1,j]<>0 then
            T1[i*j]:=IdentityMatrix(m,compact=false);
            T1[i*j][i+1,j]:=lu1[i+1,j];
          end if;
      end do; end do;
      for i from 1 to m do
        if p1[i]<>i then
          P1[i]:=IdentityMatrix(m,compact=false);
          P1[i][p1[i],i],P1[i][i,p1[i]]:=1,1;
          P1[i][p1[i],p1[i]],P1[i][i,i]:=0,0;
        end if;
      end do;
      u1:=Matrix(lu1,shape=triangular[upper]);
      p2,lu2:=LUDecomposition(u1^%T,output=[':-NAG']);
      for i from 1 to m-1 do
        for j from 1 to i do
          if lu2[i+1,j]<>0 then
            T2[i*j]:=IdentityMatrix(m,compact=false);
            T2[i*j][i+1,j]:=lu2[i+1,j];
          end if;
      end do; end do;
      for i from 1 to m do
        if lu2[i,i]<>1 then
          T3[i]:=IdentityMatrix(m,compact=false);
          T3[i][i,i]:=lu2[i,i];
        end if;
      end do;
      for i from 1 to m do
        if p2[i]<>i then
          P2[i]:=IdentityMatrix(m,compact=false);
          P2[i][p2[i],i],P2[i][i,p2[i]]:=1,1;
          P2[i][p2[i],p2[i]],P2[i][i,i]:=0,0;
        end if;
      end do;
      `if`(type(P1,table),entries(P1,':-nolist'),NULL),
      seq(seq(`if`(assigned(T1[i*j]),T1[i*j],NULL),j=1..i),i=1..m-1),
      seq(`if`(assigned(T3[i]),T3[i],NULL),i=1..min(m,n)),
      seq(seq(`if`(assigned(T2[i*j]),T2[i*j]^%T,NULL),j=i..1,-1),i=m-1..1,-1),
      `if`(type(P2,table),entries(P2,':-nolist'),NULL);
    end proc:
    
    A:=LinearAlgebra:-RandomMatrix(3,generator=1..4);
    
    ElemDecomp(A);
    
    LinearAlgebra:-Norm( `.`(%) - A);
    

    【讨论】:

      【解决方案2】:

      Mathematica 文档中的Matrix Decompositions 页面列出了所有内置的矩阵分解函数,如SingularValueDecompositionLUDecompositionCholeskyDecompositionSchurDecomposition 等。

      HTH!

      【讨论】:

        【解决方案3】:

        MATLAB 中有许多分解/分解函数(请参阅“特征值和奇异值”和“矩阵分解”下的列表 here),例如 LU factorizationorthogonal-triangular decompositionblock LDL' factorization 到仅举几例。

        【讨论】:

        • 别忘了 SVD(奇异值分解)!
        • 别忘了 Cholesky 分解!
        • 嗯,我确实链接到entire list of functions。你们真的要我在这里重新输入整个内容吗? =P
        猜你喜欢
        • 2014-02-09
        • 1970-01-01
        • 2015-03-16
        • 2015-07-29
        • 2013-08-03
        • 1970-01-01
        • 1970-01-01
        • 1970-01-01
        • 1970-01-01
        相关资源
        最近更新 更多