[Numpy-discussion] rewriting NumPy code in C or C++ or similar
ondrej at certik.cz
Mon Mar 14 16:24:35 EDT 2011
On Tue, Mar 8, 2011 at 6:25 AM, Sturla Molden <sturla at molden.no> wrote:
> Den 08.03.2011 05:05, skrev Dan Halbert:
>> Thanks, that's a good suggestion. I have not written Fortran since 1971,
>> but it's come a long way. I was a little worried about the row-major vs
>> column-major issue, but perhaps that can be handled just by remembering
>> to reverse the subscript order between C and Fortran.
> In practice this is not a problem. Most numerical libraries for C assume
> Fortran-ordering, even OpenGL assumes Fortran-ordering. People program
> MEX files for Matlab in C all the time. Fortran-ordering is assumed in
> MEX files too.
> In ANSI C, array bounds must be known at compile time, so a Fortran
> routine with the interface
> subroutine foobar( lda, A )
> integer lda
> double precision A(lda,*)
> end subroutine
> will usually be written like
> void foobar( int lda, double A);
> in C, ignoring different calling convention for lda.
> Now we would index A(row,col) in Fortran and A[row + col*lda] in C. Is
> that too difficult to remember?
> In ANSI C the issue actually only arises with small "array of arrays"
> having static shape, or convoluted contructs like "pointer to an array
> of pointers to arrays". Just avoid those and stay with 1D arrays in C --
> do the 1D to 2D mapping in you mind.
> In C99 arrays are allowed to have dynamic size, which mean we can do
> void foobar( int lda, double *pA )
> typedef double row_t [lda];
> vector_t *A = (vector_t*)((void*)&pA);
> Here we must index A[k][i] to match A(i,k) in Fortran. I still have not
> seen anyone use C99 like this, so I think it is merely theoretical.
> Chances are if you know how to do this with C99, you also know how to
> get the subscripts right. If you are afraid to forget "to reverse the
> subscript order between C and Fortran", it just tells me you don't
> really know what you are doing when using C, and should probably use
> something else.
> Why not Cython? It has "native support" for NumPy arrays.
> Personally I prefer Fortran 95, but that is merely a matter of taste.
+1 to all that you wrote about Fortran. I am pretty much switching to
it from C/C++ for all my numerical work, and then I use Cython to call
it from Python, as well as cmake for the build system.
More information about the NumPy-Discussion