[julia-users] Re: Examples of integrating Fortran code in Julia
Hi,
first of all, I am new to Julia (and Fortran). I tried to follow OP's
example and call Fortran from Julia. First, I was using the Intel Fortran
compiler and tried to compile the following .f90 file (saved as f90tojl.f90)
module m
contains
integer function five()
five = 5
end function five
end module m
as a shared library, by entering the following in the command line:
ifort f90tojl.f90 -O2 -dll -fPIC -o ifortlib.dll
the compiler ignores -fPIC since that's an unknown option apparently but
two files are created, one object file and the dll. I then try to call the
function from Julia and that's where the trouble starts.
The suggested command in this thread doesn't work because I guess the Intel
compiler has different name mangling than the gfortran compiler. So I tried
to find the name of my function, wikipedia suggests m_MP_five_ for ifort,
so I tried
julia: ccall( (:m_MP_five, "ifortlib.dll"), Int, () )
LoadError: ccall: could not find function m_MP_five in library ifortlib.dll
So my guess is that I am not using the correct name mangling. I couldn't find
anything online so I tried to view the function name via the object manager in
visual studio and an external dll viewer program.
In visual studio I got a meaningless error the external viewer just didn't do
anything (although it worked for other dll files). When I type
julia: Libdl.dlopen("ifortlib.dll")
Ptr{Void} @0x2a9d8fa0
fwiw. At this point I got so pissed that I decided to install the gfortran
compiler and just follow this thread step by step. So in the cmd window, I type:
gfortran -shared -O2 f90tojl.f90 -fPIC -o gfortlib.dll
(I get a warning that -fPIC is ignored, as written previously in this thread).
I use the dll viewer to determine the name of the function, its __m_MOD_five
indeed. Then
julia: ccall( (:__m_MOD_five, "gfortlib.dll"), Int, () )
LoadError: error compiling anonymous: could not load library "gfortlib.dll"
The specified module could not be found.
And
julia: Libdl.dlopen("gfortlib.dll")
LoadError: could not load library "gfortlib.dll"
The specified module could not be found.
And I have no clue what to do now.
[julia-users] Calling C from Julia with OpenMP (Windows)
Hi,
I am trying to call C code from Julia.It works for 'regular' C code but now
I have code that (potentially) uses OpenMP which causes trouble. I am
trying to interpolate a function using a spline approximation. In short, if
I type:
gcc -shared -O3 nspline_omp.c -o libcspline.dll -ffast-math
and then call one of the functions in the library from Julia using ccall,
it works. But judging from the speed and comparison with other languages,
only 1 thread is used. The code certainly supports multi-threading, for a
my friend of mine on his Mac, it worked. If I type:
gcc -shared -O3 nspline_omp.c -o libcspline2.dll -ffast-math -fopenmp
iand try to call it from Julia, it doesn't work anymore, in particular:
LoadError: error compiling spline: could not load library "libcspline2.dll"
The specified module could not be found.
Libdl.dlopen("libcspline2.dll")
gives the same error. As an example, consider:
#include
int banana(int hello) {
int nt;
#pragma omp parallel
{
nt = omp_get_num_threads();
}
return nt;
}
#include
int main() { printf("%d\n", banana(12)); }
which perfectly works when compiled as .exe (gcc helloworld.c -fopenmp). When I
try to compile
#include
#include
int banana(int hello) {
int nt;
#pragma omp parallel
{
nt = omp_get_num_threads();
}
return nt;
}
as a shared library using
gcc -shared -O3 helloworld2.c -o libcpsline3.dll -ffast-math -fopenmp
and then in Julia:
Libdl.dlopen("libcspline3.dll")
LoadError: could not load library "libcspline3.dll"
The specified module could not be found.
Finally, I am new to Julia and C and basically I have not much of an idea what
I am doing, so it could be an obvious mistake. Thanks for your help.
Re: [julia-users] Re: Examples of integrating Fortran code in Julia
Yes, sorry, next time I will open a new thread. I actually already have a new problem, so I am going to post that in a new thread. This particular problem here has been solved. I didn't manage to get it working with ifort (I guess there is no hope?) and with gfortran I noticed that the standard installation of this mingw was somehow 32bit. So I had to reinstall it using the 64bit version and then calling Fortran from Julia worked.
[julia-users] Re: Calling C from Julia with OpenMP (Windows)
copying the libgomp.dll into my working directory worked! Thanks!
Am Mittwoch, 20. Januar 2016 03:36:54 UTC+1 schrieb Tony Kelman:
>
> Take a look at the libcspline3.dll in Dependency Walker. You're likely
> missing a libgomp dll which you'll need to copy from your compiler
> installation and put it next to libcspline3.dll. Where did you install the
> gcc compiler from and are you calling it from inside an environment like
> cygwin or msys[2]?
>
>
> On Tuesday, January 19, 2016 at 4:15:45 PM UTC-8,
> pokerho...@googlemail.com wrote:
>>
>>
>>
>> Hi,
>>
>>
>> I am trying to call C code from Julia.It works for 'regular' C code but
>> now I have code that (potentially) uses OpenMP which causes trouble. I am
>> trying to interpolate a function using a spline approximation. In short, if
>> I type:
>>
>> gcc -shared -O3 nspline_omp.c -o libcspline.dll -ffast-math
>>
>>
>> and then call one of the functions in the library from Julia using ccall,
>> it works. But judging from the speed and comparison with other languages,
>> only 1 thread is used. The code certainly supports multi-threading, for a
>> my friend of mine on his Mac, it worked. If I type:
>>
>> gcc -shared -O3 nspline_omp.c -o libcspline2.dll -ffast-math -fopenmp
>>
>> iand try to call it from Julia, it doesn't work anymore, in particular:
>>
>> LoadError: error compiling spline: could not load library "libcspline2.dll"
>> The specified module could not be found.
>>
>>
>> Libdl.dlopen("libcspline2.dll")
>>
>> gives the same error. As an example, consider:
>>
>> #include
>>
>> int banana(int hello) {
>> int nt;
>> #pragma omp parallel
>> {
>> nt = omp_get_num_threads();
>> }
>> return nt;
>> }
>>
>>
>> #include
>> int main() { printf("%d\n", banana(12)); }
>>
>>
>> which perfectly works when compiled as .exe (gcc helloworld.c -fopenmp).
>> When I try to compile
>>
>> #include
>> #include
>>
>> int banana(int hello) {
>> int nt;
>> #pragma omp parallel
>> {
>> nt = omp_get_num_threads();
>> }
>> return nt;
>> }
>>
>> as a shared library using
>>
>> gcc -shared -O3 helloworld2.c -o libcpsline3.dll -ffast-math -fopenmp
>>
>> and then in Julia:
>>
>> Libdl.dlopen("libcspline3.dll")
>> LoadError: could not load library "libcspline3.dll"
>> The specified module could not be found.
>>
>>
>> Finally, I am new to Julia and C and basically I have not much of an idea
>> what I am doing, so it could be an obvious mistake. Thanks for your help.
>>
>>
>>
>>
[julia-users] Julia vs Matlab: interpolation and looping
Hi,
my original problem is a dynammic programming problem in which I need to
interpolate the value function on an irregular grid using a cubic spline
method. I was translating my MATLAB code into Julia and used the Dierckx
package in Julia to do the interpolation (there weren't some many
alternatives that did spline on an irregular grid as far as I recall). In
MATLAB I use interp1. It gave exactly the same result but it was about 50
times slower which puzzled me. So I made this
(http://stackoverflow.com/questions/34766029/julia-vs-matlab-why-is-my-julia-code-so-slow)
stackoverflow post.
The post is messy and you don't need to read through it I think. The bottom
line was that the Dierckx package apparently calls some Fortran code which
seems to pretty old (and slow, and doesn't use multiple cores. Nobody knows
what exactly the interp1 is doing. My guess is that it's coded in C?!).
So I asked a friend of mine who knows a little bit of C and he was so kind
to help me out. He translated the interpolation method into C and made it
such that it uses multiple threads (I am working with 12 threads). He also
put it on github here (https://github.com/nuffe/mnspline). Equipped with
that, I went back to my original problem and reran it. The Julia code was
still 3 times slower which left me puzzled again. The interpolation itself
was much faster now than MATLAB's interp1 but somewhere on the way that
advantage was lost. Consider the following minimal working example
preserving the irregular grid of the original problem which highlights the
point I think (the only action is in the loop, the other stuff is just
generating the irregular grid):
MATLAB:
spacing=1.5;
Nxx = 300 ;
Naa = 350;
Nalal = 200;
sigma = 10 ;
NoIter = 1;
xx=NaN(Nxx,1);
xmin = 0.01;
xmax = 400;
xx(1) = xmin;
for i=2:Nxx
xx(i) = xx(i-1) + (xmax-xx(i-1))/((Nxx-i+1)^spacing);
end
f_util = @(c) c.^(1-sigma)/(1-sigma);
W=NaN(Nxx,1);
W(:,1) = f_util(xx);
W_temp = NaN(Nalal,Naa);
xprime = NaN(Nalal,Naa);
tic
for banana=1:NoIter
% tic
xprime=ones(Nalal,Naa);
W_temp=interp1(xx,W(:,1),xprime,'spline');
% toc
end
toc
Julia:
include("mnspline.jl")
spacing=1.5
Nxx = 300
Naa = 350
Nalal = 200
sigma = 10
NoIter = 1
xx=Array(Float64,Nxx)
xmin = 0.01
xmax = 400
xx[1] = xmin
for i=2:Nxx
xx[i] = xx[i-1] + (xmax-xx[i-1])/((Nxx-i+1)^spacing)
end
f_util(c) = c.^(1-sigma)/(1-sigma)
W=Array(Float64,Nxx,1)
W[:,1] = f_util(xx)
spl = mnspline(xx,W[:,1])
function performance(NoIter::Int64)
W_temp = Array(Float64,Nalal*Naa)
W_temp2 = Array(Float64,Nalal,Naa)
xprime=Array(Float64,Nalal,Naa)
for banana=1:NoIter
xprime=ones(Nalal,Naa)
W_temp = spl(xprime[:])
end
W_temp2 = reshape(W_temp,Nalal,Naa)
end
@time performance()
In the end I want to have a matrix, that's why I do all this reshaping in
the Julia code. If I comment out the loop and just consider one iteration,
Julia does it in (I ran it several times, precompiling etc)
0.000565 seconds (13 allocations: 1.603 MB)
MATLAB on the other hand: Elapsed time is 0.007864 seconds.
The bottom line being that in Julia the code is much faster (more than 10 times
in this case), which should be since I use all my threads and the method is
written in C. However, if I don't comment out the loop and run the code as
posted above:
Julia:
3.205262 seconds (118.99 k allocations: 15.651 GB, 14.08% gc time)
MATLAB:
Elapsed time is 4.514449 seconds.
If I run the whole loop apparently MATLAB catches up a lot. It is still slower
in this case. In the original problem though Julia was only about 3 times
faster within the loop and once I looped through MATLAB turned out be 3 times
faster. I am stuck here, does anybody have an idea what might be going? / If I
am doing something wrong in my Julia code? A hint what I could do better?
Right now I am trying to parallelize also the loop. But that's obviously unfair
compared with MATLAB because you could also parallelize that (If you buy that
expensive toolbox).
[julia-users] Re: Julia vs Matlab: interpolation and looping
@Tomas: maybe check out Numerical Recipes in C: The Art of Scientific Computing, 2nd edition. There is also an edition for Fortran. The code that I use in C is basically from there. @Andrew: The xprime needs to be in the loop. I just made it ones to simplify but normally it changes every iteration. (In the DP problem, the loop is calculating an expecation and xprime is the possible future value of the state variable for each state of the world). Concerning the Dierckx package. I don't know about the general behaviour but for my particular problem (irregular grid + cubic spline) it is very slow. Run the following code: using Dierckx spacing=1.5 Nxx = 300 Naa = 350 Nalal = 200 sigma = 10 NoIter = 1 xx=Array(Float64,Nxx) xmin = 0.01 xmax = 400 xx[1] = xmin for i=2:Nxx xx[i] = xx[i-1] + (xmax-xx[i-1])/((Nxx-i+1)^spacing) end f_util(c) = c.^(1-sigma)/(1-sigma) W=Array(Float64,Nxx,1) W[:,1] = f_util(xx) spl = Spline1D(xx,W[:,1]) function performance2(NoIter::Int64) W_temp = Array(Float64,Nalal*Naa) W_temp2 = Array(Float64,Nalal,Naa) xprime=Array(Float64,Nalal,Naa) for banana=1:NoIter xprime=ones(Nalal,Naa) W_temp = spl(xprime[:]) end W_temp2 = reshape(W_temp,Nalal,Naa) end @time performance2(1) 30.878093 seconds (100.01 k allocations: 15.651 GB, 2.19% gc time) That's why I went on and asked my friend to help me out in the first place. I think the mnspline is really (not saying it's THE fastest) fast in doing the interpolation itself (magnitudes faster than MATLAB). But then I just don't understand how MATLAB can catch up by just looping through the same operation over and over. Intuitively (maybe I'm wrong) it should be somewhat proportional. If my code in Julia is 10 times faster whitin a loop, and then I just repeat the operation in that particular loop very often, how can it turn out to be only equally fast as MATLAB. Again, the mnspline uses all my threads maybe it has something to do with overhead, whatever. I don't know, hints appreciated.
[julia-users] Re: Julia vs Matlab: interpolation and looping
Ok. My original code certainly spends most of the time on looping the interpolation. In that sense, the example I post here is similar to the original code and hightlights the problem I am facing I think. Fwiw, I wrap the original code in a function. I also do it above, at least for the performance relevant part I guess. If I do it for the whole code above -- no difference. The first time I translated the code from MATLAB to Julia, I actually computed the expectation point by point with 2 double loops (or also with 1 loop and then reshaping) -- it was slower than vectorizing. (see also the stackoverflow post) . @Lutfullah Tomak: I don't really get that point, maybe I am misunderstanding something but the innermost part of the loop is wrapped in a function, so there should be no global variable?!
[julia-users] Re: Julia vs Matlab: interpolation and looping
Thanks for the research. I also did some testing with the original code and to me it seems like the problem has nothing to do with the interpolation method but with the memory allocation. Matlab is just faster because it doesn't reallocate memory in every iteration, as you said. It also explains my timings I guess. Matlab with 1 iteration, relatively slow but once it's looped, it doesn't reallocate memory. So I am going to look into the @inbounds and so on. I am pretty confident that the spl is faster than Matlabs interp1.
