# 1 \* X' \* X is not symmetric

**URL:** <https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444>\
**Category:** General Usage\
**Tags:** linearalgebra\
**Created:** [June 5, 2018, 1:40pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444 "2018-06-05T13:40:42Z")\
**Posts on this page:** 13\
**Page:** 1

<div class="post-metadata">

**Author:** ![e3c6](https://avatars.discourse-cdn.com/v4/letter/e/e79b87/32.png) [@e3c6](https://discourse.julialang.org/u/e3c6)\
**Post date:** [June 5, 2018, 1:40pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/1 "2018-06-05T13:40:42Z")

</div>

I find the following snippet unexpected.

```julia
srand(1);
X = rand(1000, 1000) - 0.5;
(X' *X) == (X' *X)' #true
(1 * X' *X) == (1 * X' *X)' #false

```

Multiplying by `1` changes the symmetry, probably due to floating-point arithmetic issues. But I fail to understand where the failure happens.

Let `A = X' * X`. Naively I expect that:

```julia
A[i,j] == sum(X[k,i] * X[k,j] for k = 1 : N)
A[j,i] == sum(X[k,j] * X[k,i] for k = 1 : N)

```

Since multiplication is commutative (this holds in floating-point arithmetic as well), it follows that `A[i,j] == A[j,i]`, so `A` is symmetric. I don’t see how multiplying by 1 changes this.

---

<div class="post-metadata">

**Author:** ![e3c6](https://avatars.discourse-cdn.com/v4/letter/e/e79b87/32.png) [@e3c6](https://discourse.julialang.org/u/e3c6)\
**Post date:** [June 5, 2018, 1:47pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/2 "2018-06-05T13:47:23Z")

</div>

After some discussion in the slack chat (see [Slack](https://julialang.slack.com/archives/C67910KEH/p1528205759000199)), it seems that the culprit code is here:

[https://github.com/JuliaLang/julia/blob/2c5aec90c9f9b58158eb61a516f6c2565c25eb9a/stdlib/LinearAlgebra/src/matmul.jl#L216](https://github.com/JuliaLang/julia/blob/2c5aec90c9f9b58158eb61a516f6c2565c25eb9a/stdlib/LinearAlgebra/src/matmul.jl#L216)

`1 * A` creates a copy of `A`, which the fails the identity test `A === B`, which results in a call to `gemm_wrapper!`…

But I don’t understand why that affects the result and breaks the symmetry.

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [June 5, 2018, 1:49pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/3 "2018-06-05T13:49:49Z")

</div>

I think that `1 * X' * X` is evaluated as `(1 * X') * X`, and `1 * X'` cannot be recognized as coming from the same matrix `X`.

You should use `1 * (X' * X)`, or even better, `Symmetric(X' * X) * c` for some generic real `c` (I am assuming the `1` is just for the MWE).

---

<div class="post-metadata">

**Author:** ![kristoffer.carlsson](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/kristoffer.carlsson/32/22_2.png) [@kristoffer.carlsson](https://discourse.julialang.org/u/kristoffer.carlsson)\
**Post date:** [June 5, 2018, 1:51pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/4 "2018-06-05T13:51:47Z")

</div>

> [@e3c6](#):
>
> But I don’t understand why that affects the result and breaks the symmetry.

Different algorithms computing the same thing will give different results in floating point arithmetic. Don’t use `==` to compare floating point values.

---

<div class="post-metadata">

**Author:** ![e3c6](https://avatars.discourse-cdn.com/v4/letter/e/e79b87/32.png) [@e3c6](https://discourse.julialang.org/u/e3c6)\
**Post date:** [June 5, 2018, 1:51pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/5 "2018-06-05T13:51:49Z")

</div>

Thanks. However I am not trying to fix the issue. I want to understand where it comes from to understand the kind of floating point errors I might be committing here.

---

<div class="post-metadata">

**Author:** ![e3c6](https://avatars.discourse-cdn.com/v4/letter/e/e79b87/32.png) [@e3c6](https://discourse.julialang.org/u/e3c6)\
**Post date:** [June 5, 2018, 1:53pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/6 "2018-06-05T13:53:34Z")

</div>

What exactly is the difference between `syrk_wrapper!` and `gemm_wrapper!`?

---

<div class="post-metadata">

**Author:** ![Nosferican](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/nosferican/32/9275_2.png) [@Nosferican](https://discourse.julialang.org/u/Nosferican)\
**Post date:** [June 5, 2018, 1:57pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/7 "2018-06-05T13:57:03Z")

</div>

Those are different BLAS routines. If not mistaken `syrk` (symmetric rank-k update) takes advantage of crossproducts while `gemm` (General Matrix Multiply) is general.

---

<div class="post-metadata">

**Author:** ![e3c6](https://avatars.discourse-cdn.com/v4/letter/e/e79b87/32.png) [@e3c6](https://discourse.julialang.org/u/e3c6)\
**Post date:** [June 5, 2018, 1:59pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/8 "2018-06-05T13:59:16Z")

</div>

Is there a reason for `===` instead of `==` in the snippet of code of `mul!`?

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [June 5, 2018, 1:59pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/9 "2018-06-05T13:59:30Z")

</div>

> [@e3c6](#):
>
> I want to understand where it comes from to understand the kind of floating point errors I might be committing here.

A good starting point would be

> [@PSA: floating-point arithmetic](https://discourse.julialang.org/t/psa-floating-point-arithmetic/8678):
>
> Sometimes people are surprised by the results of floating-point calculations such as julia\> 5/6 0.8333333333333334 # shouldn't the last digit be 3? julia\> 2.6 - 0.7 - 1.9 2.220446049250313e-16 # shouldn't the answer be 0? These are not bugs in Julia. They’re consequences of the IEEE-standard 64-bit binary representation of floating-point numbers that is burned into computer hardware, which Julia and many other languages use by default. Brief explanation You can t…

TL;DR: order matters, if you know (from theory) that two things are equal and want to use it like that, don’t calculate it twice, or discard one of them. Hence `Symmetric`.

---

<div class="post-metadata">

**Author:** ![StefanKarpinski](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stefankarpinski/32/24_2.png) [@StefanKarpinski](https://discourse.julialang.org/u/StefanKarpinski)\
**Post date:** [June 5, 2018, 2:06pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/10 "2018-06-05T14:06:16Z")

</div>

> [@e3c6](#):
>
> Is there a reason for `===` instead of `==` in the snippet of code of `mul!` ?

Yes, `A === B` checks if `A` and `B` are the same matrix rather than checking if they are equal.

---

<div class="post-metadata">

**Author:** ![e3c6](https://avatars.discourse-cdn.com/v4/letter/e/e79b87/32.png) [@e3c6](https://discourse.julialang.org/u/e3c6)\
**Post date:** [June 5, 2018, 2:12pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/12 "2018-06-05T14:12:27Z")

</div>

Wouldn’t the same thing work with `==` instead of `===`? Why do you need them to be the same matrix, instead of just having the same numerical values?

---

<div class="post-metadata">

**Author:** ![StefanKarpinski](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stefankarpinski/32/24_2.png) [@StefanKarpinski](https://discourse.julialang.org/u/StefanKarpinski)\
**Post date:** [June 5, 2018, 2:14pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/13 "2018-06-05T14:14:37Z")

</div>

Because checking `===` is very fast: it only requires comparing two pointers. Checking `==` is expensive: in general, it requires looking at all elements of both matrices (which may be very large).

---

<div class="post-metadata">

**Author:** ![e3c6](https://avatars.discourse-cdn.com/v4/letter/e/e79b87/32.png) [@e3c6](https://discourse.julialang.org/u/e3c6)\
**Post date:** [June 5, 2018, 2:20pm UTC](https://discourse.julialang.org/t/1-x-x-is-not-symmetric/11444/14 "2018-06-05T14:20:41Z")

</div>

I see. Thanks!
