There are several valid options / work rounds that you can use.
1. Probably the easiest way would be to arrange for python to call the n_bins function, and then use the result to call the (slightly modified version) create_hist function with unallocated output arrays of the correct size.
i.e. in Python:
n_b = n_bins(a,dr)
with the Fortran interface up to create_hist , which is 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...
This only works when you can call n_bins cheaply and create_hist outside. I know two methods for modeling distributed arrays for cases when this is not applied (i.e., code for developing array sizes is expensive and cannot be easily split).
2. First of all, you should use distributed module level arrays (described in the documentation here ). This is essentially a distributed global variable - you call your function, it saves the data in a global variable, and then you access it from Python. The downside is that it is not thread safe (that is, bad if you call create_hist at the same time)
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
Then the Python call 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. Another way that I really like is to select only arrays inside the function (that is, do not pass them to / from as parameters). However, what you are doing is a Python callback function that is called at the end and stores arrays. It is thread safe (I believe).
The fortran code looks something like this:
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
Unfortunately, this is a bit more active. You need to create a signature file using the command f2py -m fortran_module -h fortran_module.pyf my_fortran_file.f90 (this is to build a Python module called fortran_module - how you can change the names), and then change the corresponding lines to specify the size 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 this with f2py -c fortran_module.pyf my_fortran_file.f90
and then a short Python shell (to give you a simple 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 you do not need to break multithreading to avoid (but this is actually an angular case, which you probably will never see). Option 2 is easier to implement.