There's a few valid options/work-rounds that you can use.
1. Probably the simplest would be to arrange for python to call the function n_bins
, and then use the result from that to call (a slightly modified version of) the create_hist
function with non-allocatable output arrays of the correct sizes.
i.e. in Python:
n_b = n_bins(a,dr) # don't need to pass n - inferred from a
bins,hist = create_hist(a,dr,n_b) # bins and hist are auto-allocated
with the Fortran interface to create_hist
now defined as
subroutine create_hist(a, n, dr, n_b, bins, hist)
integer, intent(in) :: n
real(8), intent(in) :: a(n)
real(8), intent(in) :: dr
integer, intent(in) :: n_b
integer, intent(out) :: hist(n_b)
real(8), intent(out) :: bins(n_b+1)
! code follows
! you also need to specify the function n_bins...
That only works in cases where you can call n_bins
cheaply and from outside create_hist
. I know of 2 methods for simulating allocatable arrays for cases where that doesn't apply (i.e. the code to work out the array sizes is expensive and can't easily be separated).
2. The first is to use module level allocatable arrays (described in the documentation here). This is essentially an allocatable global variable - you call your function, it saves the data into the global variable and then you access that from Python. The disadvantage is that it isn't thread-safe (i.e. it's bad if you call create_hist
simultaneously in parallel)
module something
real(8), allocatable :: bins(:)
integer, allocatable :: hist(:)
contains
subroutine create_hist(a,n,dr)
integer, intent(in) :: n
real(8), intent(in) :: a(n)
real(8), intent(in) :: dr
integer :: n_b
n_b = n_bins(a,n,dr)
allocate(bins(n_b+1))
allocate(hist(n_b))
! code follows
end subroutine
end module
The Python call then looks like
something.create_hist(a,n,dr)
bins = something.bins # or possible bins = copy.copy(something.bins)
hist = something.hist # or possible hist = copy.copy(something.hist)
3. The other way which I quite like is to only allocate the arrays within the function (i.e. do not pass them in/out as parameters). However, what you do pass is a Python callback function that is called at the end and saves the arrays. This is thread-safe (I believe).
The fortran code then looks something like
subroutine create_hist(a,n,dr,callback)
integer, intent(in) :: n
real(8), intent(in) :: a(n)
real(8), intent(in) :: dr
external callable ! note: better to specify the type with an interface block (see http://www.fortran90.org/src/best-practices.html#callbacks)
integer :: n_b
real(8), allocatable :: bins(:)
integer, allocatable :: hist(:)
n_b = n_bins(a,n,dr)
allocate(bins(n_b+1))
allocate(hist(n_b))
! code follows
call callable(bins,hist,n_b)
end subroutine
unfortunately it's then slightly more involved. You need to create a signature file with the command f2py -m fortran_module -h fortran_module.pyf my_fortran_file.f90
(this is to build a Python module called fortran_module
- change the names as appropriate), and then modify the relevant lines of it to clarify the dimensions of the callback function:
python module create_hist__user__routines
interface create_hist_user_interface
subroutine callable(bins,hist,n_b) ! in :f:my_fortran_file.f90:stuff:create_hist:unknown_interface
real(kind=8), dimension(n_b+1) :: bins
integer, dimension(n_b) :: hist
integer :: n_b
end subroutine callable
end interface create_hist_user_interface
end python module create_hist__user__routines
compile that with f2py -c fortran_module.pyf my_fortran_file.f90
and then a short Python wrapper (to give you an easy interface) that looks like
def create_hist(a,dr):
class SaveArraysCallable(object):
def __call__(self,bins,hist):
self.bins = bins.copy()
self.hist = hist.copy()
f = SaveArrayCallable()
fortran_module.create_hist(a,dr,f)
return f.bins, f.hist
Summary
For many cases option 1 is probably the best. Otherwise option 2 or option 3. I prefer option 3 because there aren't multithreading pitfalls to avoid (but really this is a corner-case that you'll likely never see). Option 2 is easier to implement.