[julia-users] Re: Euler Maruyama for neural networks, slower than python

[Pablo Zubieta]

I just stumble upon this comment 
<https://github.com/JuliaLang/julia/issues/13157#issuecomment-169452236> on 
this GitHub issue <https://github.com/JuliaLang/julia/issues/13157>. It 
seems that the current situation in Julia won't let you match the numpy 
code yet (until that issue is solved). But from the comment in the link, 
the following should be faster (it is at least on my machine).

### Julia code
import Base.LinAlg, Base.LinAlg.BlasReal, Base.LinAlg.BlasComplex

Base.disable_threaded_libs()

function sig(x::Union{Float64,Vector})
  return 1.0 ./ (1. + exp(-x))
end

function mf_loop(Ndt::Int64,V0::Vector{Float64},V::Vector{Float64},dt::
Float64,W::Matrix{Float64},J::Matrix{Float64})
    sv = copy(V0)
    V  = copy(V0)
    for i=1:Ndt
        #sv = sig(V)
        #V = V + ( -V + J*sv ) * dt + W[:,i]
        BLAS.gemv!('N',dt,J,sv,1.0-dt,V)
        BLAS.axpy!(1., W[:,i], V)
    end
    nothing 
end

N = 100
dt = 0.1
Ndt = 10000

sigma  = 0.1
W = (randn(N,Ndt)) * sigma * sqrt(dt)
J = (rand(N,N))/N
  
V0 = rand(N)
V = copy(V0)

println( "calcul\n")
mf_loop(1,V0,V,dt,W,J)
@time mf_loop(Ndt,V0,V,dt,W,J)

Cheers.