In [31]:

    

import sympy as sp

A = []
for deriv in ['x', 'y', 'z']:
    A.append([])
    for field in ['x', 'y', 'z']:
        A[-1].append(sp.Symbol('A' + deriv + field))

A = sp.Matrix(A)

A2 = A**2
A3 = A**3
Q = -sp.horner(sp.simplify(sum(A2[i, i] for i in range(3))/2))
R = -sp.horner(sp.simplify(sum(A3[i, i] for i in range(3))/3))
print(Q)
print(R)

S = (A + A.T)/2
S2 = S**2
trS2 = sp.horner(sp.simplify(sum(S2[i, i] for i in range(3))))
print(trS2)


    

    




-Axx**2/2 - Axy*Ayx - Axz*Azx - Ayy**2/2 - Ayz*Azy - Azz**2/2
-Axx*(Axx**2/3 + Axy*Ayx + Axz*Azx) - Axy*(Ayx*Ayy + Ayz*Azx) - Axz*(Ayx*Azy + Azx*Azz) - Ayy*(Ayy**2/3 + Ayz*Azy) - Ayz*Azy*Azz - Azz**3/3
Axx**2 + Axy*(Axy/2 + Ayx) + Axz*(Axz/2 + Azx) + Ayx**2/2 + Ayy**2 + Ayz*(Ayz/2 + Azy) + Azx**2/2 + Azy**2/2 + Azz**2

In [15]:

    

def alt_measure(expr):
    POW = sp.Symbol('POW')
    count = sp.count_ops(expr, visual=True).subs(POW, 10)
    count = count.replace(sp.Symbol, type(sp.S.One))
    return count

Qalt = sp.simplify(Q, measure = alt_measure)
print(Qalt)
print(sp.count_ops(Qalt, visual=True))
Ralt = sp.simplify(R, measure = alt_measure)
print(Ralt)
print(sp.count_ops(Ralt, visual=True))


    

    




-Axx**2/2 - Axy*Ayx - Axz*Azx - Ayy**2/2 - Ayz*Azy - Azz**2/2
3*DIV + 3*MUL + NEG + 3*POW + 5*SUB
-Axx**3/3 - Axx*Axy*Ayx - Axx*Axz*Azx - Axy*Ayx*Ayy - Axy*Ayz*Azx - Axz*Ayx*Azy - Axz*Azx*Azz - Ayy**3/3 - Ayy*Ayz*Azy - Ayz*Azy*Azz - Azz**3/3
3*DIV + 16*MUL + NEG + 3*POW + 10*SUB

In [23]:

    

Qalt = - sum(A[i, (i+1)%3] * A[(i+1)%3, i] for i in range(3)) - sum(A[i, i]**2 for i in range(3))/2
print Qalt
print(sp.simplify(Qalt - Q))


    

    




-Axx**2/2 - Axy*Ayx - Axz*Azx - Ayy**2/2 - Ayz*Azy - Azz**2/2
0

In [24]:

    

Ralt = - (sum(A[i, i]*(A[i, i]**2/3 + sum(A[i, (i+j)%3]*A[(i+j)%3, i]
                                          for j in range(1, 3)))
              for i in range(3)) +
          A[0, 1]*A[1, 2]*A[2, 0] + A[0, 2]*A[1, 0]*A[2, 1])
print Ralt
print(sp.count_ops(Ralt, visual=True))
print(sp.simplify(Ralt - R))


    

    




-Axx*(Axx**2/3 + Axy*Ayx + Axz*Azx) - Axy*Ayz*Azx - Axz*Ayx*Azy - Ayy*(Axy*Ayx + Ayy**2/3 + Ayz*Azy) - Azz*(Axz*Azx + Ayz*Azy + Azz**2/3)
6*ADD + 3*DIV + 13*MUL + NEG + 3*POW + 4*SUB
0

In [37]:

    

AxxAxx = A[0, 0]**2
AyyAyy = A[1, 1]**2
AzzAzz = A[2, 2]**2
AxyAyx = A[0, 1]*A[1, 0]
AyzAzy = A[1, 2]*A[2, 1]
AzxAxz = A[2, 0]*A[0, 2]

Qalt = - (AxxAxx + AyyAyy + AzzAzz)/2 - AxyAyx -AyzAzy - AzxAxz
print(sp.simplify(Qalt - Q))
Ralt = - (A[0, 0]*(AxxAxx/3 + AxyAyx + AzxAxz) +
          A[1, 1]*(AyyAyy/3 + AxyAyx + AyzAzy) +
          A[2, 2]*(AzzAzz/3 + AzxAxz + AyzAzy) +
          A[0, 1]*A[1, 2]*A[2, 0] +
          A[0, 2]*A[1, 0]*A[2, 1])
print sp.simplify(Ralt - R)
#print sp.simplify(Ralt + A.det())
trS2alt = (AxxAxx + AyyAyy + AzzAzz +
           ((A[0, 1] + A[1, 0])**2 +
            (A[1, 2] + A[2, 1])**2 +
            (A[2, 0] + A[0, 2])**2)/2)
print sp.simplify(trS2alt - trS2)

print Ralt
print A.det()


    

    




0
0
0
-Axx*(Axx**2/3 + Axy*Ayx + Axz*Azx) - Axy*Ayz*Azx - Axz*Ayx*Azy - Ayy*(Axy*Ayx + Ayy**2/3 + Ayz*Azy) - Azz*(Axz*Azx + Ayz*Azy + Azz**2/3)
Axx*Ayy*Azz - Axx*Ayz*Azy - Axy*Ayx*Azz + Axy*Ayz*Azx + Axz*Ayx*Azy - Axz*Ayy*Azx

In [ ]: