We experienced an issue with Matrix-Vector (or row-vector times column-vector) product. In the example the result should return zero: a = [1.3190624946305824, 1.3190624946305824] b = [-3.74 3.74] # Matrix-Vector-Product b * a # not always zero sum(b .* a) # zero On Machine 1, only the se…

[EDIT: FMA is (still) worrying, see my next post.] That seemed to be very worrying. And there seemed no limit to the error. I changed the code slightly to see if I can get larger err: err = 8.731107228996319e43 a = 6.208364721873767e29 b = 8.740381590017811e29 but actually the calculated fma is…

Sorry, should have looked at the relative error instead…

Actually real worrying, yes, value very far from 0.0: julia> fma(4.738285026766486e28, 4.0686784105686624e29, -(4.738285026766486e28 * 4.0686784105686624e29)) 8.475399346267984e41 So what does that mean? I was starting to think we might have discovered a hardware bug, but possibly this is the only…

It’s very worrying (based on my last post). Can you check: a = [4.738285026766486e28, 4.738285026766486e28] b = [-4.0686784105686624e29 4.0686784105686624e29] i.e. b * a that I get 0.0 for, on your affected machine (well both). And someone with alternative BLAS, also Elrod’s.

On Machine 1 we get julia> a = [4.738285026766486e28, 4.738285026766486e28] 2-element Vector{Float64}: 4.738285026766486e28 4.738285026766486e28 julia> b = [-4.0686784105686624e29 4.0686784105686624e29] 1×2 Matrix{Float64}: -4.06868e29 4.06868e29 julia> b*a 1-element Vector{Float64}: 8.475399346…

[image] Palli: Is there a way to fix BLAS to avoid FMA when dangerous? I don’t think that the conclusion here is that FMA is “dangerous”. Sometimes it will give more accurate results than non-fused operations, sometimes (probably less often) it will happen to give less accurate results, but…

Consider a decimal 2-significant-digit floating point type: 1.6*1.6 # = 2.56 --> 2.6 fma(1.6,1.6,-2.6) # = -0.04 --> -0.04 In the case of fma(-x,y,x*y), the result is the roundoff error of x*y. That is to say, it’ll never be zero unless x*y requires no rounding of the result. FMA is doing exactly…

Assuming your real use-case is ill-conditioned (*), but much larger than the toy example here, you could perhaps benefit from the compensated dot product implementation in AccurateArithmetic.jl : julia> using AccurateArithmetic julia> a = [4.738285026766486e28, 4.738285026766486e28] 2-element Vecto…

[image] Palli: 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…

[image] stevengj: [image] Palli: Is there a way to fix BLAS to avoid FMA when dangerous? I don’t think that the conclusion here is that FMA is “dangerous”. Sometimes it will give more accurate results than non-fused operations, sometimes (probably less often) it will happen to give less a…

Inaccurate Matrix Multiplication

General Usage

stevengj February 28, 2022, 1:58pm 39

Note that you can already use Julia’s generic linear-algebra backends with extended-precision floating-point types (double-double etc.). See also the discussion at: MPLAPACK (multiprecision blas+lapack)

Topic		Replies	Views
Mul! generates machine noise General Usage matrices	21	720	February 21, 2024
X * y + z does not automatically use FMA instruction Performance fast-math	31	2921	July 27, 2023
Benchmark MATLAB & Julia for Matrix Operations Performance	148	20808	October 15, 2019
Why mul! is so fast? General Usage question , linearalgebra	7	7054	November 26, 2019
Vector - Matrix - Vector multiplication Performance	19	4662	March 14, 2021

Inaccurate Matrix Multiplication

Related topics