MATLAB1：对任意一个常数矩阵，求其所有实可交换矩阵的MATLAB算法

（以下是matlab算法，我觉得这个算法的结果我觉得还是比较有美感的，但是中间有点效率不行，但我没时间改进，希望大家能指点一下。这是我的期末作业。虽然我也觉得这算法没多大用，但我还是先公开发布一波吧。数学证明不难，我也就不发布了。大家要的话可以私信我）
举个例子：

>> A=[0 1 0 0 0 0 0;0 0 1 0 0 0 0;0 0 0 1 0 0 0 ;0 0 0 0 1 0 0;0 0 0 0 0 1 0;0 0 0 0 0 0 1;0 0 0 0 0 0 0];
>> mat(A)
时间已过 18.127244 秒。
ans =
[ x7_7, x6_7, x5_7, x4_7, x3_7, x2_7, x1_7]
[    0, x7_7, x6_7, x5_7, x4_7, x3_7, x2_7]
[    0,    0, x7_7, x6_7, x5_7, x4_7, x3_7]
[    0,    0,    0, x7_7, x6_7, x5_7, x4_7]
[    0,    0,    0,    0, x7_7, x6_7, x5_7]
[    0,    0,    0,    0,    0, x7_7, x6_7]
[    0,    0,    0,    0,    0,    0, x7_7]

function [ result ] = mat( A )ticg=1;[n,n]=size(A);
X=sym('x',[n,n]) ;b=A*X-X*A; 
B= zeros(n^2,n^2);
[r,c] = size(b);    
for i = 1:r       for j = 1:cb(i,j);    for p=1:nfor q=1:nB(g,(p-1)*n+q)=diff(b(i,j),X(p,q));endendg=g+1; end
end
B=rref(B);
r=rank(B);
C=zeros(n^2,2);
s=1;t=1;f=1;
for i=f:n^2f=f+1;for j=1:n^2if B(i,j)==1B(i,j)=0;C(s,t)=i;C(s,t+1)=j;s=s+1;breakelsef=f+1;endend
end
X1=reshape(X',1,n^2);
K=(-1).*B*X1';
g=1;d=0;for i=1:rif C(i,1)<C(i,2)&&C(i,1)~=C(i,2)-d;d=d+1;endif mod(i+d,n)==0k=n;X(g,k)=K(i,1);else k=mod(i+d,n);g=ceil((i+d)/n);X(g,k)=K(i,1);end       
end
toc
result=X;
end

如果你喜欢我的作品，please，帮我点个赞吧。这会鼓励我在CSDN写出更多的作品哦！谢谢您耐心的观看哦！