On 7/2/2012 7:17 PM, Casey W. Stark wrote:
Hi numpy.
Does anyone know if f2py supports allocatable arrays, allocated inside
fortran subroutines? The old f2py docs seem to indicate that the
allocatable array must be created with numpy, and dropped in the
module. Here's more background to explain...
I have a fortran subroutine that returns allocatable positions and
velocities arrays. I wish I could get rid of the allocatable part, but
you don't know how many particles it will create until the subroutine
does some work (it checks if each particle it perturbs ends up in the
domain).
module zp
implicit none
contains
subroutine ics(..., num_particles, particle_mass, positions, velocities)
use data_types, only : dp
implicit none
... inputs ...
integer, intent(out) :: num_particles
real (kind=dp), intent(out) :: particle_mass
real (kind=dp), intent(out), dimension(:, :), allocatable ::
positions, velocities
...
end subroutine
end module
I tested this with a fortran driver program and it looks good, but
when I try with f2py, it cannot compile. It throws the error "Error:
Actual argument for 'positions' must be ALLOCATABLE at (1)". I figure
this has something to do with the auto-generated "*-f2pywrappers2.f90"
file, but the build deletes the file.
If anyone knows an f2py friendly way to rework this, I would be happy
to try. I'm also fine with using ctypes if it can handle this case.
Best,
Casey
_______________________________________________
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion
Hi,
... not sure what version of numpy you are using but it can be done with
the f2py included with numpy 1.6.x. Here is a (hopefully intelligible)
code fragment from a working application:
module Grid_
! a representation of a two-dimensional, rectangular, (image) sensor
grid
! -- this component is used in the "smoothing" portion of the
Thematic MaP solution algorithm (smoothness prior probabilities)
use runtime_
implicit none
integer, private, parameter :: begin_ = 1, end_ =
2, row_ = 1, column_ = 2
integer, private, dimension(begin_:end_) :: rows__ = (/0,0/),
columns__ = (/0,0/) ! grid rows and columns extent
integer, allocatable, dimension(:,:) :: pixel_neighbors
! the neighbors of a specified pixel -- ((row,column), ...,
(row,column))
! + this array is managed (allocated, deallocated, populated)
by the neighbors_of() procedure
! + when allocated, the first dimension is always 2 -- pixel
(row,column)
! + it would be preferable for this array to be the return
result, of procedure neighbors_of(),
! but F2PY does not seem to support an allocatable array as
a function result
contains
[snip]
function neighbors_of (pixel) result(completed)
! determines all second-order (grid) neighbors of *pixel*
! -- returns .true. if successful
! -- use the error_diagnostic() procedure if .false. is
returned
! -- upon successful completion of this procedure,
Grid_.pixel_neighbors contains the neighbors of *pixel*
! -- *pixel* may not be on the grid; in this case,
Grid_.pixel_neighbors is not allocated and .false. is returned
! integer, dimension(row_:column_), intent(in) :: pixel !
pixel to be considered ... f2py does not support this
integer, dimension(2), intent(in) :: pixel ! pixel under
consideration (row,column)
logical :: completed
! integer, dimension(begin_:end_) :: neighbor_rows,
neighbor_columns ... f2py does not support this
! integer, dimension(row_:column_) :: neighbor ... f2py does
not support this
integer, dimension(2) :: neighbor_rows,
neighbor_columns ! each is: (begin,end)
integer, dimension(2) :: neighbor ! (row,column)
integer :: row, column, code, count
character (len=100) :: diagnostic
completed = .false.
count = 0 ! *pixel* has no neighbors
if (allocated (pixel_neighbors)) deallocate (pixel_neighbors)
if (is_interior (pixel)) then
count = 8 ! interior pixels have eight,
second-order neighbors
neighbor_rows (begin_) = pixel(row_)-1
neighbor_rows (end_) = pixel(row_)+1
neighbor_columns (begin_) = pixel(column_)-1
neighbor_columns (end_) = pixel(column_)+1
else if (is_border (pixel)) then
count = 5 ! non-corner, border pixels have five,
second-order neighbors, but ...
if (is_corner (pixel)) count = 3 ! corner pixels have
three, second-order neighbors
neighbor_rows (begin_) = max (pixel(row_)-1, rows__(begin_))
neighbor_rows (end_) = min (pixel(row_)+1, rows__(end_))
neighbor_columns (begin_) = max (pixel(column_)-1,
columns__(begin_))
neighbor_columns (end_) = min (pixel(column_)+1,
columns__(end_))
end if
if (count > 0) then
allocate (pixel_neighbors(row_:column_,count), stat=code,
errmsg=diagnostic)
if (code /= 0) then
call set_error (code,diagnostic)
return
end if
count = 0
do row = neighbor_rows(begin_), neighbor_rows(end_)
do column = neighbor_columns(begin_), neighbor_columns(end_)
neighbor(row_) = row; neighbor(column_) = column
if (neighbor(row_) == pixel(row_) .and.
neighbor(column_) == pixel(column_)) cycle ! neighbor is pixel
count = count + 1; pixel_neighbors(:,count) = neighbor
end do ! columns
end do ! rows
completed = .true.
end if
end function neighbors_of
end module Grid_
_______________________________________________
NumPy-Discussion mailing list
NumPy-Discussion@scipy.org
http://mail.scipy.org/mailman/listinfo/numpy-discussion
