0

I am trying to make a wrapper for a Fortran function with f2py from this code: http://arxiv.org/e-print/1601.07458v1 . It is from a paper concerning reduced density matrices, (http://arxiv.org/abs/1601.07458v1, and the code works fine when running it in fortran.

Now trying to make a wrapper with f2py on the Fortran file "partial_trace.f90", results in

/tmp/tmpBwIuHh/src.linux-x86_64-2.7/partial_tracemodule.c:325:15: error: invalid type argument of unary ‘*’ (have ‘int’)
w_Dims[0]=2**nqb;
            ^

I am pretty unfamiliar with both Fortran and C, and couldn't find a solution on google. What could be the reason for the error? I suspect its the use of a Fortran derived type, but I'm not so sure about it.

Any hints would be appreciated!

edit: - the fortran code (tar archive): arxiv.org/e-print/1601.07458v1 - the command I ran: f2py -c partial_trace.f90 -m partial_trace - its output: pastebin.com/g7QNnaCR

Felix Huber
  • 133
  • 3
  • Welcome on Stack Overflow. Use tag [tag:fortran] and add the version only if you need to restrict the question to that particular version. You should add more information. You must show the command which produced the error. You should also show the code, namely "partial_trace.f90". Be sure to read http://stackoverflow.com/help/how-to-ask I stress that the code must be in the question itself, not as a link to a paper. – Vladimir F Героям слава Feb 02 '16 at 17:33
  • 1
    I wonder if you're confusing `**` and `^`. The former is how Fortran expresses exponentiation, the latter is used by many other, lesser, languages. – High Performance Mark Feb 02 '16 at 18:03
  • Yes, but it seems that the `**` appears in the compiled C based f2py module for Python. It shouldn't appear there and it is not a place where users put their input. That's why I asked for the exact command line and the code. – Vladimir F Героям слава Feb 02 '16 at 18:38
  • Please, next time, as I already wrote: "Questions seeking debugging help (*"why isn't this code working?"*) must include the desired behavior, a specific problem or error and the shortest code necessary to reproduce it **in the question itself**. " In the question itself. This is en exact quote from a reason to close a question as off-topic. Including the emphasis. The code must be in the question, not as a link to somewhere. The same holds for the output. Do not use pastebin! Include the output **into the question itself**. – Vladimir F Героям слава Feb 04 '16 at 11:00
  • @VladimirF The compile error happens in a temporary file generated by f2py. – jofel Feb 04 '16 at 11:07
  • @VladimirF: The code is more than 300 lines long. How am I supposed to include it in the post? Recall that my knowledge of Fortran is close to zero, so "shortest code necessary to reproduce it in the question itself" is in this case 300 lines of code. Same goes for output, which is roughly 120 lines long. I posted above the (to me) most telling error message. – Felix Huber Feb 04 '16 at 12:54
  • You could have easily included only the subroutine which causes the error. – Vladimir F Героям слава Feb 04 '16 at 12:58
  • Seriously, I had pretty much no clue about how to construct a minimal example. Just because it seems to be easy in your eyes, it may not be for others, not knowing basic fortran syntax. In hindsight: yes, I could have, if I had known how to. – Felix Huber Feb 04 '16 at 13:28

2 Answers2

1

There seems to be a bug in f2py. It does not translate powers correctly.

You can fix this by creating the c file first (do not use the -c flag of f2py) and changing

W_Dims[0]=2**nqb;

to

W_Dims[0]=1<<nqb;

by hand and compile the source code manually. You could also use an ipow function from The most efficient way to implement an integer based power function pow(int, int) and use

W_Dims[0]=ipow(2,nqb);
Community
  • 1
  • 1
jofel
  • 3,297
  • 17
  • 31
  • Version of f2py, I use Fortran daily. – Vladimir F Героям слава Feb 04 '16 at 11:10
  • I can reproduce the bug with f2py from numpy 1.8.2 – jofel Feb 04 '16 at 11:23
  • I can it reproduce now too. I couldn't with a simple program I wrote myself. It was horrible to do it from OP's link. Download the gz. Deflate. Realise it is not a .gz, but a tar.gz, even if it is not named as such, untar... Not a MCVE indeed. – Vladimir F Героям слава Feb 04 '16 at 11:28
  • I can reproduce the bug also with f2py from numpy 1.10.4. – jofel Feb 04 '16 at 11:29
  • Sorry guys. I my kubuntu I downloaded the file, it recognised the filetype and uncompressed it with a double click.. I'm very new to fortran and never used f2py before, so I can't even make the code minimal, because I basically don't even know fortran syntax. That's why I provided the full original file. I have f2py Version: 2, numpy Version: 1.9.2 – Felix Huber Feb 04 '16 at 12:28
  • Thanks Jofel! The problem was indeed the power. As compiling it by hand raises the difficulty of the task to yet another, I simply circumvent the whole exponentiation to "solve" the problem. Thanks a lot for your help guys! – Felix Huber Feb 04 '16 at 12:58
0

The problem in Numpy happens when the power operator is used in array dimension declaration. This is a Minimal, Complete, and Verifiable example

subroutine s(n, a)
  integer :: n, a(1:2**n)
end

Basically, you delete (or comment) unneeded lines until you stop getting the error. That's why it is no way necessary to include all 300 lines of code, but 1 subroutine suffices. Even the original form of the subroutine is OK, it is small.

As a workaround one can explicitly pass the array dimension:

subroutine s(n, na, a)
  integer :: n, na, a(1:na)
end

or use assumed shape array:

subroutine s(n, a)
  integer :: n, a(1:)
end

but be careful, this requires an explicit interface (module).

If you use the power operator for normal computation, it is OK. I reported the issue at https://github.com/numpy/numpy/issues/7184

Community
  • 1
  • 1
Vladimir F Героям слава
  • 57,977
  • 4
  • 76
  • 119