with(linalg): with(combinat): Help:=proc() print(`A collection of MAPLE procedures to compute the standard polynomial in k symbolic matrices efficiently by making use of symmetry.`): print(`List of procedures:`): print(`perm_sign(A), FSP(mraw,k,vars), FSPL(Li1,vars), FSP1(mraw,k,vars), SFSPL (Li1,vars)`): end: #perm_sign(A): Given a permutation A, compute its sign. #Example: perm_sign([2,1,3])=-1 perm_sign:= proc(A::list) local B,i,j,n,count: B:=A: n:=nops(B): count:=0: for i from 1 to n-1 do if B[i] <> i then j:=i+1: while j<=n and B[j]<>i do j:=j+1: od: if (B[j]<>i) then ERROR("The input is not a permutation"): fi: B[j]:=B[i]: count:=count+1: fi: od: if B[n] <> n then ERROR("The input is not a permutation"): fi: (-1)^count: end proc: #FSP(mraw,k,vars): Given the first term 'mraw'=m_1*m_2*...*m_k, #compute the standard polynomial in 'k' matrices: m_1,m_2,...m_k. Assume #the input matrices have the same general form, whose elements #are variables in set 'vars', and index of a variable signifies #the matrix it belongs to. Note the result is in 'compact form', #meaning a monomial in an element actually means a sum of monomials #generated by permutating the k matrices in that monomial # (times the sign of the permutation). # Thus if the standard polynomial in 'k' matrices: m_1,m_2,...m_k # is equal to 0, a zero matrix will be output. FSP:=proc(mraw,k,vars) local i,j,d,r,c,mono,m,newm,shadow,Li_mono,s,ele,newele,vars1: m:=simplify(mraw): d:=rowdim(m): if coldim(m)<>d then ERROR("The input matrix must be a square matrix."): fi: newm:=matrix(d,d,[0$(d*d)]): for r from 1 to d do for c from 1 to d do ele:=m[r,c]: newele:=0: #Break a matrix element into a list of monomials. if type(ele,`*`) or type(ele,indexed) then Li_mono:=[ele]: else Li_mono:=[op(ele)]: fi: for s from 1 to nops(Li_mono) do mono:=Li_mono[s]: shadow:=subs({seq(seq(vars[i][j]=vars[i]^j,j=1..k), i=1..nops(vars))}, mono); vars1:=sort([op(indets(shadow))]): # First, we remove(by skipping) the monomials with repeated letters . # Justification: for example, if a1a2 is part of the monomial, # then the sum of monomials generated by permutating the k # matrices in that monomial (times the sign of the permutation) # will vanish. if nops(vars1)=k then #Now we normalize the monomial #print(mono,shadow,vars1,map2(degree,shadow,vars1)): newele:=newele+perm_sign(map2(degree,shadow,vars1)) *lcoeff(mono)*convert(vars1,`*`): fi: od: newm[r,c]:=simplify(newele): od: od: return(op(newm)): end: #FSPL(Li1,vars): compute mraw=Li1[1]*Li1[2]*...*Li1[k] and call FSP(mraw,k,vars ) FSPL:=proc(Li1,vars) local k,mraw,i: k:=nops(Li1): mraw:=Li1[1]: for i from 2 to k do mraw:=multiply(mraw,Li1[i]): od: return(FSP(mraw,k,vars)): end: #A reduced version of FSP, it only removes monomials with repeated letters FSP1:=proc(mraw,k,vars) local i,j,d,r,c,mono,m,newm,shadow,Li_mono,s,ele,newele,vars1: m:=simplify(mraw): d:=rowdim(m): if coldim(m)<>d then ERROR("The input matrix must be a square matrix."): fi: newm:=matrix(d,d,[0$(d*d)]): for r from 1 to d do for c from 1 to d do ele:=m[r,c]: newele:=0: if type(ele,`*`) or type(ele,indexed) then Li_mono:=[ele]: else Li_mono:=[op(ele)]: fi: for s from 1 to nops(Li_mono) do mono:=Li_mono[s]: shadow:=subs({seq(seq(vars[i][j]=vars[i],j=1..k), i=1..nops(vars))}, mono); # we remove(skip) the monomials with repeated letters. if nops(indets(shadow))=k then newele:=newele+mono: fi: od: newm[r,c]:=simplify(newele): od: od: return(op(newm)): end: #SFSPL(Li1,vars) is a faster version of FSPL(Li1,vars). It prunes the #monomials with repeated letters while the product Li1[1]*Li1[2]*...*Li1[k] #is being formed. SFSPL:=proc(Li1,vars) local k,mraw,i: k:=nops(Li1): mraw:=Li1[1]: for i from 2 to k do if i=k then return(FSP(multiply(mraw,Li1[i]),i,vars)): fi: mraw:=FSP1(multiply(mraw,Li1[i]),i,vars): od: end: