Symmetric Matrix Fails Symmetry Check

[Ben Goodrich]

Hi,

The following code:

Lambda = Matrix(5, 2, lambda i, j: Symbol("Lambda_%d%d" % (i, j)))
Upsilon = eye(2)
Upsilon[0,1] = Upsilon[1,0] = Symbol("Upsilon_01")
Sigma = Lambda * Upsilon * Lambda.transpose()
Sigma == Sigma.transpose() # WTF
False

disputes the fact that Sigma is necessarily symmetric, and thus an
exception would be thrown by any function that first checks the
symmetry of a matrix.

It appears as if this is some (lack of) simplification issue

(Sigma[0,1] - Sigma[1,0]) == 0
False
Sigma[0,1] - Sigma[1,0]
Λ₁₀⋅(Λ₀₀ + Λ₀₁⋅Υ₀₁) + Λ₁₁⋅(Λ₀₁ + Λ₀₀⋅Υ₀₁) - Λ₀₀⋅(Λ₁₀ + Λ₁₁⋅Υ₀₁) - Λ₀₁⋅
(Λ₁₁ + Λ₁₀⋅Υ₀₁)

Does anyone have a suggestion about how best to convince Sympy that
Sigma is symmetric?

Thanks,
Ben

[smichr]

Think of sympy as a tool not a mathematician ;-) *You* know they are
equal, but sympy doesn't go through contortions to figure out if two
things are equal or not when they don't have the same structure. Your
question is answered as FAQ 1 at http://tinyurl.com/nwa87j where it
reminds us that

"The equality operator (==) tests whether expressions have identical
form, not whether they are mathematically equivalent."

So, for example, consider the following mathematical failure

>>> var('x y')
(x, y)
>>> x*(1+y)==x+x*y
False

To sympy, it looks like (at the highest level) that a Mul is being
compared to an Add. So they aren't in the same form.

If you generate expresions which are taking on different forms for the
same quantities, you might get by with a simple expansion to show that
they are equal:

>>> a,b = x*(1+y),x+x*y
>>> a==b
False
>>> a.expand() == b.expand()
True
>>> (a-b)==0
False
>>> (a-b).expand()==0
True
>>>

Looking at your matrix entries you will see:

>>> sa=list(iter(Sigma))
>>> st=list(iter(Sigma.transpose()))
>>> sa[0]
Lambda_00*(Lambda_00 + Lambda_01*Upsilon_01) + Lambda_01*(Lambda_01 +
Lambda_00*Upsilon_01)
>>> st[0]
Lambda_00*(Lambda_00 + Lambda_01*Upsilon_01) + Lambda_01*(Lambda_01 +
Lambda_00*Upsilon_01)
>>> sa[0]==st[0]
True

Well that make sense...it's the corner which doesn't get moved in the
transpose. It's the non-diagonals that are the problem:

>>> print sa[1]
Lambda_10*(Lambda_00 + Lambda_01*Upsilon_01) + Lambda_11*(Lambda_01 +
Lambda_00*Upsilon_01)
>>> print st[1]
Lambda_00*(Lambda_10 + Lambda_11*Upsilon_01) + Lambda_01*(Lambda_11 +
Lambda_10*Upsilon_01)
>>> sa[1]==st[1]
False

Solution? Same as the simple example above:

>>> (sa[1]-st[1]).expand()==0
True

And this suggests a little helper function knowing that a simple
expansion will do:

>>> def mat_eq(a, b):
... z = zip(iter(a), iter(b))
... for zi in z:
... if (zi[0]-zi[1]).expand():
... return False
... else:
... return True
...
>>> mat_eq(Sigma, Sigma.transpose())
True

HTH,
/chris

[Mateusz Paprocki]

Hi,

On Wed, Feb 03, 2010 at 03:40:31AM -0800, smichr wrote:
> And this suggests a little helper function knowing that a simple
> expansion will do:
>
> >>> def mat_eq(a, b):
> ... z = zip(iter(a), iter(b))
> ... for zi in z:
> ... if (zi[0]-zi[1]).expand():
> ... return False
> ... else:
> ... return True
> ...
> >>> mat_eq(Sigma, Sigma.transpose())
> True
>

Or you can use a one-liner, like this:

In [5]: Sigma.applyfunc(simplify) == Sigma.transpose().applyfunc(simplify)
Out[5]: True

--
Mateusz

[Ben Goodrich]

Thanks Chris and Mateusz. I reopened issue 1472

http://code.google.com/p/sympy/issues/detail?id=1472

with a patch to substitute the one-liner and add a test.

Ben