Wrap dotu from BLAS and rename dot to dotc for complex.

[andreasnoack]

No description provided.

[kmsquire]

Would it worthwhile to leave dotc as dot? Alternatively, perhaps an
additional parameter could be added for dispatch.

I can imagine users trying to use dot on complex arrays and being confused.

Cheers, Kevin

On Thursday, January 2, 2014, Andreas Noack Jensen wrote:

---

You can merge this Pull Request by running

git pull https://github.com/JuliaLang/julia anj/dotu

Or view, comment on, or merge it at:

#5283
Commit Summary

• Wrap dotu from BLAS and rename dot to dotc for complex.

File Changes

• M base/linalg/blas.jlhttps://github.com/Wrap dotu from BLAS and rename dot to dotc for complex. #5283/files#diff-0(46)
• M base/linalg/matmul.jlhttps://github.com/Wrap dotu from BLAS and rename dot to dotc for complex. #5283/files#diff-1(11)
• M test/blas.jlhttps://github.com/Wrap dotu from BLAS and rename dot to dotc for complex. #5283/files#diff-2(101)

Patch Links:

• https://github.com/JuliaLang/julia/pull/5283.patch
• https://github.com/JuliaLang/julia/pull/5283.diff

[andreasnoack]

I thought about that but decided that the most transparent thing to do was to consistently use the BLAS names inside the BLAS module.

[kmsquire]

For users familiar with blas, that's fine. For those that aren't, how about
making both dot and dotc equivalent for complex arrays?

On Thursday, January 2, 2014, Andreas Noack Jensen wrote:

I thought about that but decided that the most transparent thing to do was
to consistently use the BLAS names inside the BLAS module.

—
Reply to this email directly or view it on GitHubhttps://github.com//pull/5283#issuecomment-31457088
.

[goszlanyi]

I would also prefer making dot and dotc equivalent for complex arrays.

[johnmyleswhite]

+1 for having dot act as dotc.

[andreasnoack]

Just to avoid misunderstandings. Outside the BLAS module dot is still defined for complex types and does conjugate dot product. Inside the BLAS module there is dotc and dotu but no dot for complex types.

[johnmyleswhite]

Ah, I was confused and thought the proposal was to change things in a way that would affect code outside the BLAS module. Using BLAS names in the BLAS module definitely seems like the right way to go.

[jiahao]

lgtm

[lindahua]

This caused a test failure on my machine:

dhlin@julia:make testall
    JULIA test/all
    From worker 5:       * strings
    From worker 3:       * keywordargs
    From worker 6:       * unicode
    From worker 7:       * collections
    From worker 4:       * numbers
    From worker 8:       * hashing
    From worker 9:       * remote
    From worker 2:       * core
    From worker 9:       * iobuffer
    From worker 9:       * arrayops
    From worker 3:       * linalg
    From worker 8:       * blas
    From worker 6:       * fft
    From worker 2:       * dsp
exception on 8: ERROR: assertion failed: |:(BLAS.dotc(z1,z2)) - :(sum(.*(conj(z1),z2)))| <= 0.0023841858
  :(BLAS.dotc(z1,z2)) = 0.0f0 + 0.0f0im
  :(sum(.*(conj(z1),z2))) = -2.3032098f0 + 0.3048067f0im
  difference = 2.3232913 > 0.0023841858

[lindahua]

The dotc function also caused segfault when I tried to call it:

julia> z1 = complex(randn(3), randn(3))
3-element Array{Complex{Float64},1}:
  0.534293+0.334639im
 -0.584508-0.784822im
   0.787467+1.59046im

julia> z2 = complex(randn(3), randn(3))
3-element Array{Complex{Float64},1}:
  0.206515-1.04789im
 0.339515-0.997745im
  0.889693-1.06187im

julia> BLAS.dotc(z1, z2)
Segmentation fault: 11

Here is my version info:

Julia Version 0.3.0-prerelease+800
Commit eda220b (2014-01-02 18:25 UTC)
Platform Info:
  System: Darwin (x86_64-apple-darwin13.0.0)
  CPU: Intel(R) Core(TM) i7-3720QM CPU @ 2.60GHz
  WORD_SIZE: 64
  BLAS: libmkl_rt
  LAPACK: libmkl_rt
  LIBM: libopenlibm

I am using MKL.

[andreasnoack]

I can confirm the problem with MKL. My guess is that the problem has been there all the time, but before this pr the exported dot did not call BLAS for complex arguments and the BLAS tests did not test dot. I could find this

http://software.intel.com/en-us/forums/topic/487184

and will try the Fortran test program now.

[andreasnoack]

I have tried the test program from the link to the Intel forum and I get a segfault with gfortran but not ifort. Hence I think this must be a problem with gfortran and MKL.

[andreasnoack]

There is some clash between the calling conventions here. If I call zdotc as a subroutine when using MKL it works

module myBLAS
    import Base.LinAlg.BlasInt
    function dotc(x::Vector{Complex{Float64}}, y::Vector{Complex{Float64}})
        result = Array(Complex{Float64}, 1)
        ccall(("zdotc_", Base.libblas_name), Void,
            (Ptr{Complex{Float64}}, Ptr{BlasInt}, Ptr{Complex{Float64}}, Ptr{BlasInt},
             Ptr{Complex{Float64}}, Ptr{BlasInt}),
            result, &length(x), x, &1, y, &1)
        return result[1]
    end
end

I don't know how we should deal with this.

[lindahua]

Probably we should use vendor-specific code here?

[JeffBezanson]

There are still problems with returning complex values, especially on 32-bit where it doesn't work at all.

[kmsquire]

Apologies! I should have looked closer at the pull request.

On Thursday, January 2, 2014, Andreas Noack Jensen wrote:

Just to avoid misunderstandings. Outside the BLAS module dot is still
defined for complex types and does conjugate dot product. Inside the BLASmodule there is
dotc and dotu but no dot for complex types.

—
Reply to this email directly or view it on GitHubhttps://github.com//pull/5283#issuecomment-31470326
.

[jiahao]

Apologies for pulling the trigger early