13

I'm porting some MATLAB code to Numpy. This task includes stripping the MEX out of some C++ code and replacing it with equivalent calls to Numpy's C-API. One problem is that the MEX code treats the incoming data as Fortran-ordered because that's how MATLAB orders its arrays. Numpy, on the other hand, uses a C ordering by default.

Short of completely re-writing the MEX code for C ordering, I can:

  • (A) Reorder the arrays that go into the C code with .copy('F') and reorder those that come out with .copy('C')
  • (B) Figure out how to get numpy to "emulate" MATLAB by doing everything in Fortran order from the get-go.

Option A -- currently implemented -- works just fine but is terribly inefficient. Anybody know how to make option B work?

BrianTheLion
  • 2,618
  • 2
  • 29
  • 46
  • Possible duplicate: http://stackoverflow.com/questions/7688304/how-to-force-numpy-array-order-to-fortran-style – cyborg Oct 27 '11 at 11:08
  • @cyborg - Ah! You're right. I think mine is a little clearer, though. Hopefully the gods will be kind to me. – BrianTheLion Oct 27 '11 at 15:33

1 Answers1

5

My approach to this problem (when I wrap fortran code with f2py) is to explicitly declare all the relevant numpy arrays in fortran order, because numpy can happily work with them transparently, and it even works well combining fortran and C order arrays. Unfortunately, it seems that numpy operations do not preserve fortran order. Therefore you should pre-allocate your destination arrays which will be passed to MEX in fortran order, for example:

A = np.empty((10, 10))
B = np.empty((10,2))
# fill A and B with the data you want

C = np.empty((10, 2), order='F')
C[:] = np.dot(A, B) # note that this index is essential
# C is then passed to your MEX routine

I'm not sure if this is much more efficient than your solution A, since the assignment has an implicit copy.

However, there should be no need to reorder the fortran arrays coming out of your MEX routine - numpy will deal with them quite transparently, as long as it knows what order they are in.

DaveP
  • 6,952
  • 1
  • 24
  • 37
  • PS: I would also like to be able to set numpy to default to fortran order, but I don't think it's possible. – DaveP Oct 27 '11 at 03:44