Julia Fortran function call returning an array: unknown function, segfault?

Question

I want to call functions in my Julia Fortran library. In this case, I have an eye function that takes an Integer and returns a two-dimensional array of integers.

The Fortran module is compiled into a shared library using

 $ gfortran -shared -fPIC -o matrix_routines.so matrix_routines.f90 

And after that I try to call it from Julia’s interactive interpreter (the name is derived from nm ):

 julia> n=5 5 julia> ccall( (:__matrix_routines_MOD_eye, "/path/to/library/matrix_routines.so"), Array{Int64,2} , (Ptr{Int64},), &n ) 

This, however, immediately leads Julia to give me a segfault:

 signal (11): Segmentation fault __matrix_routines_MOD_eye at /path/to/library/matrix_routines.so (unknown line) anonymous at no file:0 unknown function (ip: -1137818532) jl_f_top_eval at /usr/bin/../lib/julia/libjulia.so (unknown line) eval_user_input at REPL.jl:53 jlcall_eval_user_input_19998 at (unknown line) jl_apply_generic at /usr/bin/../lib/julia/libjulia.so (unknown line) anonymous at task.jl:95 jl_handle_stack_switch at /usr/bin/../lib/julia/libjulia.so (unknown line) julia_trampoline at /usr/bin/../lib/julia/libjulia.so (unknown line) unknown function (ip: 4199613) __libc_start_main at /usr/bin/../lib/libc.so.6 (unknown line) unknown function (ip: 4199667) unknown function (ip: 0) zsh: segmentation fault (core dumped) julia 

Am I calling the function wrong? What is the correct function name? (This does not seem to be just an eye , as it does not work either.)

As an additional question: does Julia do anything with the memory orientation of the resulting arrays? Fortran and Julia are both the main columns, but I wonder if, due to ccall (), Julia might think that this should transform them?

 module matrix_routines implicit none private public :: eye contains pure function eye(n,offset) result(um) !{{{ integer, intent(in) :: n integer, intent(in), optional :: offset integer, dimension(n,n) :: um integer :: i, l, u, os um = 0 l = 1 u = n os = 0 if (present(offset)) then os = offset end if if (abs(os) < n) then if (os > 0) then u = n - os else if (os < 0) then l = 1 - os end if do i=l, u um(i, i+os) = 1 end do end if end function eye !}}} end module matrix_routines 

Answer 1

There are a couple of questions with your approach. Returning the array directly to julia is problematic because Fortran arrays are not C compatible unless special conditions are met. When you make the array interoperable (add bind(C) to your procedure and give the array type C), the compiler ( gfortran ) will complain:

 Error: Return type of BIND(C) function 'um' at (1) cannot be an array 

To get around this, we can return the array through a dummy argument. You will want to make this an intent(inout) argument and build an array in julia to avoid memory / scope problems when creating an array in Fortran.

Secondly, the optional argument is problematic and reduces Julia docs, I'm not sure if it is even supported. Note that even Fortran cannot call Fortran with optional arguments without an explicit interface, and since Julia does not interact with the .mod file and seems to expect C to do something, it probably won't work (and Fortran 2008 15.3 .7 p2.6 seems to say that it is not supported). However, there are workarounds - you can create several Fortran procedures with a different number of arguments, and then call the procedure with optional arguments from them.

First, consider this Fortran module, which started with your example, but compared to what is needed to demonstrate the interaction:

 module matrix_routines implicit none private public :: eye contains pure subroutine eye(n,um) bind(C,name="eye") !{{{ use, intrinsic :: iso_c_binding, only: c_int implicit none integer(c_int), intent(in) :: n integer(c_int), intent(inout), dimension(n,n) :: um integer :: i, j do j=1,n do i=1,n um(i,j) = i+j end do end do end subroutine eye !}}} end module matrix_routines 

Note that I moved um to a dummy argument inout , and since we are not returning a value, this changed the procedure to a subroutine. I also removed the optional argument. I also used the C interaction types and associated the C name with the procedure. You can compile this as in your question.

In Julia you can now do the following:

 julia> n = 2 2 julia> um = zeros(Int32, n, n) 2x2 Array{Int32,2}: 0 0 0 0 julia> ccall((:eye, "matrix_routines.so"), Void, (Ptr{Int32}, Ptr{Array{Int32,2}}), &n, um) julia> um 2x2 Array{Int32,2}: 2 3 3 4 julia> n = 4 4 julia> um = zeros(Int32, n, n) 4x4 Array{Int32,2}: 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 julia> ccall((:eye, "matrix_routines.so"), Void, (Ptr{Int32}, Ptr{Array{Int32,2}}), &n, um) julia> um 4x4 Array{Int32,2}: 2 3 4 5 3 4 5 6 4 5 6 7 5 6 7 8 

Note that we can call the function simply :eye , since we used inter bind(C,name="eye") C in our Fortran.

And finally, if we change the do loop in my Fortran example to um(i,j) = i*10+j , we will see that transposition does not occur in the array:

 julia> um 3x3 Array{Int32,2}: 11 12 13 21 22 23 31 32 33 

---

There could be several things that could cause your segfault to be private — inconsistent data types, problems with the return type, problems with an optional argument, or inconsistencies in the actual call and the variables passed.