# Inaccurate Matrix Multiplication

**URL:** <https://discourse.julialang.org/t/inaccurate-matrix-multiplication/76903>\
**Category:** General Usage\
**Tags:** linearalgebra\
**Created:** [February 22, 2022, 11:04am UTC](https://discourse.julialang.org/t/inaccurate-matrix-multiplication/76903 "2022-02-22T11:04:45Z")\
**Posts on this page:** 1\
**Showing post:** 4

<div class="post-metadata">

**Author:** ![Palli](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/palli/32/3380_2.png) [@Palli](https://discourse.julialang.org/u/Palli)\
**Post date:** [February 22, 2022, 12:19pm UTC](https://discourse.julialang.org/t/inaccurate-matrix-multiplication/76903/4 "2022-02-22T12:19:50Z")

</div>

> [@judober](#):
>
> Given that matlab uses MKL we expected an issue there but using MKL.jl does not change the results.

I think MKL.jl should give the same result as MKL in MATLAB (assuming it does use MKL, and the exact same version on the same processor). I think you’re not actually getting 0.0 there, just MATLAB showing the value as such? Some languages, including it, show more user-friendly (well that’s debatable) version of your number by default.

> [@judober](#):
>
> Is the difference the expected inaccuracy due to Floating Point calculation?

The answer you get is because of floating-point (and more…), but the difference isn’t strictly because of it. The answer “yes” would be a simplification.

As you can see at [PSA: floating-point arithmetic - #8 by simonbyrne](https://discourse.julialang.org/t/psa-floating-point-arithmetic/8678/8)

Across machines, you “might get non-deterministic results”, for “BLAS or LAPACK operations” (also for e.g. “muladd”), and not just the default BLAS, i.e. openBLAS, also for MKL. And not just in Julia. That’s too be expected across languages using these dependencies, such as Python and MATLAB.

The reason is, e.g. “skylake” vs “tigerlake” even though the floating point-math is deterministic, and should I believe be the same across those (and all x86) architectures, that would only apply for individual scalar operations (or exact same SIMD instructions), but for the full code using them, when it’s not the same (or run in same sequence). It’s just that OpenBLAS has separate code for each architecture, different assembly code. You can use @Elrod’s pure Julia code, but I’m not sure it’s NOT special-cased too (or if not yet, then might be later). Another thing I could think of is, if you have many threads then that might, or likely will, add non-determinism, but I doubt for such small matrices. You could try disabling it.

So in short, this is to be expected, and isn’t strictly considered a bug. But:

If you need your answer to **contain** the right value, you could look into interval arithmetic (Julia has packages for), and use with Chris Elrod’s code. I’m not sure the end-points will be deterministic, actually. Another thing to look into are Posits (Unum) instead of traditional floating point.

Note while -3.583133801120832e-16 isn’t exactly 0.0 it is close enough, while not:

```julia
julia> -3.583133801120832e-16 ≈ 0.0
false

```

It’s been enough for me to use that operator (when porting from MATLAB), an alias for isapprox with the defaults.

Note, however:

```julia
julia> isapprox(-3.583133801120832e-16, 0.0; atol=0.05)
true

```

I’ve never had to use absolute tolerance (atol) before, or relative (rtol), so I’m not sure which values are most appropriate, I just took that one from the docs.

---

_[View the full topic](https://discourse.julialang.org/t/inaccurate-matrix-multiplication/76903)._