# Symv producing incorrect results

**URL:** <https://discourse.julialang.org/t/symv-producing-incorrect-results/7548>\
**Category:** Numerics\
**Tags:** blas\
**Created:** [December 5, 2017, 8:10pm UTC](https://discourse.julialang.org/t/symv-producing-incorrect-results/7548 "2017-12-05T20:10:36Z")\
**Posts on this page:** 3\
**Page:** 1

<div class="post-metadata">

**Author:** ![Elrod](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/elrod/32/22461_2.png) [@Elrod](https://discourse.julialang.org/u/Elrod)\
**Post date:** [December 5, 2017, 8:10pm UTC](https://discourse.julialang.org/t/symv-producing-incorrect-results/7548/1 "2017-12-05T20:10:36Z")

</div>

Can anyone reproduce the following bug?

```julia
julia> A = randn(10,20) |> X -> X * X';

julia> x = rand(10);

julia> sA = Symmetric(A);

julia> all( A .== sA )
true

julia> (A * x)'
1×10 RowVector{Float64,Array{Float64,1}}:
 11.538 7.90415 -1.03592 -3.84546 15.5033 3.37954 7.7948 4.26935 6.29426 5.57295

julia> (sA * x)'
1×10 RowVector{Float64,Array{Float64,1}}:
 21.5629 15.7696 -4.14023 -10.6005 28.2271 5.50057 17.9761 6.60855 6.29426 5.57295

julia> Base.LinAlg.BLAS.symv('U', A, x)'
1×10 RowVector{Float64,Array{Float64,1}}:
 21.5629 15.7696 -4.14023 -10.6005 28.2271 5.50057 17.9761 6.60855 6.29426 5.57295

```

I do on both v0.6.1 and v0.7.0. Both are linked to a local BLAS build. Additionally, `symv` is slow:

```julia
julia> y = similar(x);

julia> @benchmark Base.LinAlg.BLAS.symv!('U', 1.0, $A, $x, 0.0, $y)
BenchmarkTools.Trial: 
  memory estimate: 0 bytes
  allocs estimate: 0
  --------------
  minimum time: 1.014 μs (0.00% GC)
  median time: 1.278 μs (0.00% GC)
  mean time: 1.426 μs (0.00% GC)
  maximum time: 187.736 μs (0.00% GC)
  --------------
  samples: 10000
  evals/sample: 10

julia> @benchmark A_mul_B!($y, $A, $x)
BenchmarkTools.Trial: 
  memory estimate: 0 bytes
  allocs estimate: 0
  --------------
  minimum time: 75.414 ns (0.00% GC)
  median time: 75.662 ns (0.00% GC)
  mean time: 77.367 ns (0.00% GC)
  maximum time: 144.184 ns (0.00% GC)
  --------------
  samples: 10000
  evals/sample: 971

```

On a laptop, with the JuliaPro v0.6.0.1 MKL version (Julia v0.6), I cannot reproduce the bug. Besides getting the correct answer, symv was also much faster – while `A_mul_B` took 200 ns, symv took only 150 ns.  
(EDIT: v0.7.0 on my laptop also gets the correct answer / can’t reproduce the bug.)

Perhaps this is a clue? Is it calling the wrong function somehow?

Also, the last two elements of the vector match.

---

<div class="post-metadata">

**Author:** ![fredrikekre](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/fredrikekre/32/1688_2.png) [@fredrikekre](https://discourse.julialang.org/u/fredrikekre)\
**Post date:** [December 5, 2017, 8:52pm UTC](https://discourse.julialang.org/t/symv-producing-incorrect-results/7548/2 "2017-12-05T20:52:23Z")

</div>

> [@Elrod](#):
>
> Can anyone reproduce the following bug?

No.

> [@Elrod](#):
>
> Additionally, symv is slow:

Yes, `symv` is generally slower than `gemv`.

---

<div class="post-metadata">

**Author:** ![Elrod](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/elrod/32/22461_2.png) [@Elrod](https://discourse.julialang.org/u/Elrod)\
**Post date:** [December 5, 2017, 10:02pm UTC](https://discourse.julialang.org/t/symv-producing-incorrect-results/7548/3 "2017-12-05T22:02:56Z")

</div>

> No

Sigh… Seems `OpenBLAS`, not Julia is at fault:

```fortran
program symvtest
implicit none
integer, parameter :: n = 20, p = 12
real, dimension(n, p) :: Xb
real, dimension(p, p) :: A
real, dimension(p) :: x, y, ys

call random_number(Xb)
call random_number(x)

A = matmul( transpose(Xb), Xb )
y = matmul( A, x )
call ssymv('U', p, 1.0, A, p, x, 1, 1.0, ys, 1)

print *, 'Correct:', y
print *, 'ssymv:', ys

end program symvtest

```

```bash
$ gfortran-7.2 symvtest.f08 -o symvt2 -pthread /opt/OpenBLAS/lib/libopenblas.a

$ ./symvt2
 Correct: 21.2189522 18.4703217 23.1896820 17.4168549 24.2416267 19.9038563 19.7339249 18.0757847 21.3725090 25.3771706 24.4452839 22.0691833    
 ssymv: 34.6876221 30.6038113 -1435453.88 28.6646252 40.5542870 32.9165344 -1692865.25 30.3047485 21.3725128 25.3771706 24.4452820 22.0691833 

```

Reproducing the results, down to the first 8 answers being wrong (I get the same thing when changing 10 → 12 in the Julia example). For reference, OpenBLAS [issue](https://github.com/xianyi/OpenBLAS/issues/1388).

> Yes, symv is generally slower than gemv.

`symv` is faster on my laptop at `n=10` and `n=2000`, but perhaps that’s MKL.
