[SciPy-dev] Cython, f2py and GSoC

[Dag Sverre Seljebotn]

Sturla Molden wrote:
> On 3/26/2009 8:43 PM, Dag Sverre Seljebotn wrote:
>
>> b) While there's a clear technical seperation between Fortran <-> C and
>> C <-> Cython, the main challenge is the same both places (dealing with
>> strided arrays on a C level) and so it's benefitial to have them within
>> the same project.
>
> Do you need to pass strided arrays to Fortran? Would not a local copy
> suffice?

Thanks for the feedback. You know a lot more about this than me so I'm
very grateful for any input.

Is that necesarrily easier? It's a different approach which I think would
involve about the same amount of work (shifting some logic from the new
tool to the Cython compiler instead...) There's two tools generating
respectively C and Fortran code, and one of those tools need to emit the
logic to do the copying in either C or Fortran.

(Note that Cython doesn't want to directly depend on NumPy as such, but
PEP 3118, and would need to reimplement copying to contiguous memory
anyway. But it's not much code either way.)

Some sort custom struct (or change of signature in some other way)
containing the bounds of the array are necesarry anyway.

But even if I think the work amount is the same, contiguous memory might
simpler and therefore more robust to maintain etc. Hmm..

>
> Also consider that the Fortran 2003 function C_F_POINTER (used to
> convert a C pointer to a Fortran array pointer) assumes the C pointer
> references a continuous array. You would have to pass in the original
> unstrided array, create one Fortran pointer with C_F_POINTER and then
> index this with appropriate strides. Most likely the Fortran compiler
> would react to these strides by making its own local copy, which would
> make these efforts a bit futile. It is easier to always pass contiguous
> arrays to Fortran, and in most cases just as efficient.

Yes, that was exactly the plan. And if the dimension ordering is different
I couldn't find a set of slicing operations to do it, so some Fortran
function (reshape?) must be used which copies anyway (which is probably
just as well).

Since it seems doing the copying either side is about the same amount of
work, my idea was that copying in Fortran would have the advantage that it
moves the decision to the Fortran compiler,  which is smarter and knows
more.

Is it possible to write Fortran code in a way which doesn't make the
compiler copy, in case you know that memory use is a bigger problem than
speed? If so I'd say that settles the issue for my part, no point "turning
a feature off".

If not, I suppose this should ideally be resolved by benchmarks. Perhaps
in some obscure cases an important stride is already in register C side
making the copying more efficient there :-)

Dag Sverre