allocate on the fly vectors into fortran from python using f2py

Question

I have the following start of a subroutine called by python:

Called from python using: libomp.cic(pos,mpart,boxsize,dim,npart)

subroutine cic(pos,mpart,boxsize,dim,npart,mesh)
use omp_lib
implicit none 
real*8 mpart
!f2py intent(in) :: mpart
real*8 boxsize
!f2py intent(in) :: boxsize
integer dim
!f2py intent(in) :: dim
integer*8 npart
!f2py intent(in) :: npart
real*8, dimension(0:dim-1,0:dim-1,0:dim-1) :: mesh
!f2py intent(out) :: mesh
real*8, dimension(npart,3) :: pos
!f2py intent(in) :: pos

Error when running python script:

error: (shape(pos,0)==npart) failed for 1st keyword npart: cic:npart=47099622
WARNING: Failure executing file: <test.py>

For some reason it works when I set npart as an integer but it can't receive npart from python and create the array pos on the fly. It works for mesh using dim but not pos, why? Does it have something to do with the ordering in which they are called?

Answer 1

pos is declared to be intent(in) which means that it is supposed to be allocated when you enter the routine. Note that this is completely different than using an automatic array in fortran since pos is a dummy argument. The problem is that since pos is intent(in), f2py is checking it's declared dimensions with the dimensions that you're passing in (mpart) because npart is also supposed to be intent(in). Since they don't match, you're getting an error.

In other words, you have 2 intent(in) variables and the values that you're passing in for them are incompatible. The actual solution really depends on what the subroutine is supposed to be doing, and what the intent of the dummy arguments actually should be.

If pos is really intent(in), then you could in principle change npart to intent(hide) in which case, numpy will pass the npart to your subroutine based on the allocated dimensions of pos. e.g., I believe it will be called as:

libomp.cic(pos,mpart,boxsize,dim)