Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Calling Fortran subroutines with array arguments from Julia

I can call compile this fortran code 'test.f90'

subroutine test(g,o)
double precision, intent(in):: g
double precision, intent(out):: o
o=g*g
end subroutine

with

gfortran -shared -fPIC test.f90 -o test.so

and create this wraper function test.jl for Julia:

function test(s)
  res=Float64[1]
  ccall((:test_, "./test.so"), Ptr{Float64}, (Ptr{Float64}, Ptr{Float64}), &s,res);
  return res[1]
end

and run these command with the sought for output:

julia> include("./test.jl")    
julia> test(3.4)
11.559999999999999

But I want to return an array instead of a scalar. I think I have tried everyting, including using the iso_c_binding, in this answer. But everything I try throws errors at me looking like this:

ERROR: MethodError: `convert` has no method matching convert(::Type{Ptr{Array{Int32,2}}}, ::Array{Int32,2})
This may have arisen from a call to the constructor Ptr{Array{Int32,2}}(...),
since type constructors fall back to convert methods.
Closest candidates are:
  call{T}(::Type{T}, ::Any)
  convert{T}(::Type{Ptr{T}}, ::UInt64)
  convert{T}(::Type{Ptr{T}}, ::Int64)
  ...
 [inlined code] from ./deprecated.jl:417
 in unsafe_convert at ./no file:429496729

As an example, I'd like to call the following code from julia:

subroutine arr(array) ! or  arr(n,array)
implicit none
integer*8, intent(inout) :: array(:,:)
!integer*8, intent(in) :: n
!integer*8, intent(out) :: array(n,n)
integer :: i, j

do i=1,size(array,2) !n
    do j=1,size(array,1) !n
        array(i,j)= j+i
    enddo
enddo
end subroutine

Using the commented out variant is also an alternative as only changing the argument doesn't seem to be useful when calling from julia.

So how do I call fortran subroutines with arrays from Julia?

like image 648
Jonatan Öström Avatar asked Oct 03 '16 22:10

Jonatan Öström


2 Answers

When using ccall, Julia arrays should be passed as pointers of the element type, along with an extra argument describing the size.

Your example test.f90 should be:

subroutine arr(n,array)
implicit none
integer*8, intent(in) :: n
integer*8, intent(out) :: array(n,n)
integer :: i, j

do i=1,size(array,2) !n
    do j=1,size(array,1) !n
        array(i,j)= j+i
    enddo
enddo
end subroutine

Compiled as before with

gfortran -shared -fPIC test.f90 -o test.so

Then in Julia:

n = 10
X = zeros(Int64,n,n) # 8-byte integers
ccall((:arr_, "./test.so"), Void, (Ptr{Int64}, Ptr{Int64}), &n, X)
like image 162
Simon Byrne Avatar answered Oct 14 '22 00:10

Simon Byrne


Suggestion: '&n' will cause invalid syntax error for newer versions. And, using Ref is safer than Ptr.

You may try:

ccall((:arr_, "./test.so"), Cvoid, (Ref{Int64}, Ref{Int64}), Ref(n), X)
like image 2
Charlie Chu Avatar answered Oct 13 '22 23:10

Charlie Chu