DLAMI OpenBLAS Behaviors Different from Build from Source OpenBLAS

[Zha0q1]

OpenMathLib/OpenBLAS#2779
Please refer to ^ this issue. Float point behaviors are different in ILP64 OpenBLAS and the vanilla 32-bit LP64 OpenBLAS due to difference in the usage of FMA/blocking. Essentially we will lose precision when two float numbers are more than significand apart. This is better handled by vanilla openblas and hence we may get different numerical results on same inputs migrating to ILP64 openblas when using gemm, gemv, and others.

We should discuss the implication of this issue before preceding to migrate to ILP64 openblas

@sandeep-krishnamurthy @access2rohit @josephevans
@szha @leezu

Update: I am using dlami p2.x16large. This issue might vary on different hardware according to Martin at openblas

[sandeep-krishnamurthy]

Thanks Zhaoqi. This is applicable for tensors with greater than 2pow31 right? It doesn't have any implications for regular tensors?

[Zha0q1]

Thanks Zhaoqi. This is applicable for tensors with greater than 2pow31 right? It doesn't have any implications for regular tensors?

After more discussion with openblas people my understanding is that loss of precision will happen whenever two numbers to be added are certain magnitude apart from each other. (in fact loss of precision always happen but is especially bad under such condition)

Direct quotes from them:

Expected precision is as of standard netlib blas with single addition per loop. If you need to sum vaues more than significand+1 apart you need bigger floating point variable with obviously bigger significand range.

What the naive loop reaches is limit of "significand"
https://en.wikipedia.org/wiki/Floating-point_arithmetic#IEEE_754:_floating_point_in_modern_computers
i.e 2^24 (it rounds properly even near the end) i.e close to 1,6777216 ^ 10^7

With that said much less error will be accumulated with regularly sized tensors. Regular openblas is able to avoid this error accumulation by FMA/blocking so hence the discrepancy we saw between ilp64 openblas and vanilla openblas.

My demo program calculates (1, 1, 1, ....) dot (1, 1, 1, ...).T. The error will start to accumulate when the sum reaches significand. Here are some result recaps

ILP64 openblas
With shape = (2^31, 1)

2147483647 // 2^31
1.93274e+09 // ilp64 gemm answer
2.14748e+09 // accurate value

With shape = (2^29)

536870911 // 2^29
5.36871e+08 // ilp64 openblas gemm answer
5.36871e+08 // accurate value

compared with vanilla openblas
With shape = (2^31, 1)

2147483647 // 2^31
2.14748e+09 // vanilla gemm answer
2.14748e+09 // accurate value

[Zha0q1]

Update: The difference in numerical value is caused by the fact that DLAMI (what i referred to earlier as vanilla openblas) uses a special (possibly optimized) openblas lib.

I did more tests and found that:

DLAMI openblas (LP64, int32): no loss of precision on my example 2.14748e+09
openblas dev (ILP64, int64): loss of precision 1.93274e+09
openblas dev (LP64, int32): loss of precision 1.93274e+09

In that sense the results are consistent between the 32 and 64 openblas if we build from the same source with the same configs (except for int size).

Now the problem seems to be that whoever uses DLAMI will potentially get different numerical results from the ones who build openblas dev from source.

The DLAMI openblas version is:

ubuntu@ip-172-31-11-14:~$ ll /usr/local/lib/ |grep blas
lrwxrwxrwx  1 root root         30 Aug  7 01:38 libopenblas.a -> libopenblas_haswellp-r0.2.20.a
lrwxrwxrwx  1 root root         31 Aug  7 01:38 libopenblas.so -> libopenblas_haswellp-r0.2.20.so*
lrwxrwxrwx  1 root root         31 Aug  7 01:38 libopenblas.so.0 -> libopenblas_haswellp-r0.2.20.so*
-rw-r--r--  1 root root   28498048 Aug  7 01:12 libopenblas_haswellp-r0.2.20.a
-rwxr-xr-x  1 root root   14344848 Aug  7 01:13 libopenblas_haswellp-r0.2.20.so*

My openblas dev build is:

-rw-r--r--  1 root root   65650108 Aug 15 21:38 libopenblasp-r0.3.10.dev.a
-rwxr-xr-x  1 root root   39439000 Aug 15 21:38 libopenblasp-r0.3.10.dev.so*

Whether the root cause comes from different version or different config needs further investigation

---

update:
openblas team said the kernel for gemm was changed earlier this year and that might introduce the different rounding due to FMA operations

[sandeep-krishnamurthy]

Thanks Zhaoqi. If the issue turns out to be DLAMI open blas we should be able create an issue with AWS DLAMI team and get it resolved by bringing this issue to their notice.

Outside DLAMI what is your summary ? Are you saying there will be no loss of precision? You Github iasue with openblas team says it is expected to have different behavior?

[Zha0q1]

Thanks Zhaoqi. If the issue turns out to be DLAMI open blas we should be able create an issue with AWS DLAMI team and get it resolved by bringing this issue to their notice.

Will do by this week.

Outside DLAMI what is your summary ? Are you saying there will be no loss of precision? You Github iasue with openblas team says it is expected to have different behavior?

The summary is that under the same version and configuration 32 and 64 openblas will give consistent results. This is what we and the openblas team expect.

Both my 32 and 64 openblas builds have the same loss of precision. Openblas team explains that this is caused by the nature of float number computation and is not a bug. However, it appears that DLAMI ships with an openblas binary which does not suffer from the same loss of precision.

In other words, migrating from 32 to 64 openblas with other settings and version being identical should not break things. But users can get different numerical results if 1. they use opeblas binaries with different optimization levels and 2. they have (large) computations that accumulates numerical error.

---

update:
openblas team said the kernel for gemm was changed earlier this year and that might introduce the different rounding due to FMA operations. DLAMI does ship with a pretty old version of openblas (0.2.20 vs. latest = 0.3.10)

[leezu]

openblas team said the kernel for gemm was changed earlier this year and that might introduce the different rounding due to FMA operations. DLAMI does ship with a pretty old version of openblas (0.2.20 vs. latest = 0.3.10)

Should it be reverted / fixed?

[Zha0q1]

openblas team said the kernel for gemm was changed earlier this year and that might introduce the different rounding due to FMA operations. DLAMI does ship with a pretty old version of openblas (0.2.20 vs. latest = 0.3.10)

Should it be reverted / fixed?

So according to them this was done to boost performance. They probably won't revert