I am attempting to wrap a fortran package in Julia. In many of the functions, there are 2D arrays that are specified as ALLOCATABLE. Does Julia have a way to handle this? And what is the proper way to pass these values in?
For example, if in the fortran code we have
REAL(dp),DIMENSION(:,:),ALLOCATABLE:: P !dp has 64-bits precision
and then in Julia, would we do something like
P=zeros(n,m)
ccall((fun_name_, "path/to/so"), Cvoid, (Ref{Cdouble}), P)
I don’t think there is a portable way to do this. It’s better to first modify/wrap your Fortran package to use an ISO_C_BINDING API. See e.g. c++ - Allocating memory in C for a Fortran allocatable - Stack Overflow
Do we need to wrap this because of the name mangling only, or are there issues with the data types as well?
As I understand it, the problem is the data type. Fortran compilers represent allocatable arrays with some internal data structure that they don’t document.
I think I’m almost getting it here.
So I settled on fortran code that looks like
module tester
!
use iso_c_binding
!
contains
!
subroutine test1(n,arr) BIND(C)
!
implicit none
integer(c_int):: n
type(c_ptr),intent(inout):: arr
real(SELECTED_REAL_KIND(15,307)),pointer:: farr(:)
integer:: i
!
allocate(farr(4))
n=4
do i=1,n
farr(i)=i
end do
arr = c_loc(farr)
!
end subroutine test1
!
end module tester
I’m able to call it in Julia with
n = Ref(Int32(2))
arr=pointer_from_objref([2.0,3.0])
a=ccall((:test1,"/home/jmezier2/test/test.so"),Cvoid,(Ptr{Int32},Ptr{Float64}),n,arr)
println(unsafe_pointer_to_objref(arr))
but this only prints the first 2 elements of my array, instead of the 4 expected after allocating. Is there a way to recover the other 2 elements?
Don’t do use pointer_from_objref. That gives you a pointer to an internal Julia datastructure, not raw array data. To pass an array, you would normally just do:
arr = [2.0, 3.0]
ccall((:test1,"/home/jmezier2/test/test.so"),Cvoid,(Ref{Int32},Ptr{Float64}),n,arr)
Except here, that’s not right either, because arr is a pointer argument passed by reference in Fortran, so it should be a Ref{Ptr{Float64}} (corresponding to a double** in C).
So I think you need something like:
arr = Ref{Ptr{Float64}}() # uninitialized pointer
ccall((:test1,"/home/jmezier2/test/test.so"),Cvoid,(Ref{Int32},Ref{Ptr{Float64}}),n,arr)
a = unsafe_wrap(Array, arr[], Int(n[]))
to get the pointer back and wrap it in an Array.
However, you’ll have to be careful about de-allocating the array. You’ll eventually need to deallocate it on the Fortran side, but make sure you don’t reference the Julia array a after that.
Note also that you declared arr as inout but you are only using it as an output argument? Similarly with n?
Ok, wow that works directly as advertised. But why can’t we do something like
arr = Ptr{Float64}} # uninitialized pointer
ccall((:test1,"/home/jmezier2/test/test.so"),Cvoid,(Ref{Int32},Ref{Ptr{Float64}}),n,arr)
Isn’t convert supposed to handle this? Or does it get confused with the Ref{Ptr}?
Besides the fact there are unbalanced braces so what you wrote is a syntactic error, there’s no such a thing as an uninitialised pointer in Julia. Side note, I personally find the @ccall macro much more readable than plain ccall
No. Even if you corrected the code to arr = Ptr{Float64}(0), it won’t work because arr is bound to an immutable object.
The problem is that passing this arr for a Ref{Ptr{Float64}} argument doesn’t pass a “pointer to arr” — it’s not analogous to passing &arr in C. What it does is to construct a new Ref(arr) mutable object that wraps a copy of arr’s value, and the contents of this new object are what get mutated inside the Fortran subroutine. However, you can’t access the mutated contents after the ccall returns, because you didn’t give a name to this new Ref object.
To give a simpler example, suppose you have a C function:
void getRandomNumber(int *x) {
*x = 4; // chosen by fair dice roll, guaranteed to be random
}
Suppose we call it in Julia by:
ccall(:getRandomNumber, Cvoid, (Ref{Cint},), 17)
How would you get the value 4 out? Clearly, you can’t (and don’t want) to “change the value of 17” … 17 is an immutable value. What the ccall does is to construct a Ref{Cint}(17) object, which gets passed as a pointer to the C function, and the contents of that object get mutated.
But if you do
x = 17
ccall(:getRandomNumber, Cvoid, (Ref{Cint},), x)
it is exactly equivalent to passing 17 directly. The function call cannot change the value of x, for the same reason that
y = x
y = 4
does not change the value of x to 4. See also the manual section on assignment vs. mutation.
Instead, you do:
xref = Ref{Cint}(17)
ccall(:getRandomNumber, Cvoid, (Ref{Cint},), xref)
after which xref[] (the contents of the Ref object) should be 4.
This is analogous in C to:
int xref[] = {17}; // a mutable container for 17, analogous to Ref{Cint}(17)
getRandomNumber(xref);
after which *xref should be 4.
Ah, I see, so Julia will automatically convert everything for you, but then you lose access to the reference. So if I just passed a pointer in, I would have access to the old pointer. But I changed to a new pointer in the fortran function and the old pointer would have no knowledge of this. But a reference to that pointer would.
This makes sense. Thank you.