Is this a correct python port of a fortran function call

Question

I'm trying to port fortran code to python (call me crazy), and am wondering whether I'm handling the input of a variable to a function in a correct way.

note that I've read:

Fortran77: what does the asterisk (*) means in subroutine argument?

How can a scalar be passed to a vector (1D array) to a Fortran subroutine?

What is * in FORTRAN 77

, but my question stands.

Focusing on the variable wrk(ia) used in curfit.f to input to subroutine fpcurf.f. wrk is a 1-d array.

curfit.f source

 c  ..array arguments..
      real*8 ... wrk(lwrk) ...
 c  ... no values assigned to wrk here, it is still empty initialized ...
 
 call fpcurf(...,
 * wrk(ia),...)
 ...

(from the links I understand the * is just a line continuation.

curfit.py

...
n, c, fp, ier = call_fpcurf(..., wrk[ia], ...)
...

fpcurf.f takes the variable wrk(ia) as variable a, which is defined as an array. Now I understand that this somehow becomes a "dummy" array (whatever that means).

fpcurf.f source

  recursive subroutine fpcurf(...,
 *   ...,a,...)
  implicit none
c  ..
c  ..scalar arguments..
      ...
c  ..array arguments..
      real*8 ...,
     * ...,a(nest,k1),...

Since wrk is empty initialized in curfit.f, and the element wrk(ia) is inputted to fpcurf.f as a, which has array dimensions (nest, k1), I basically ignore the inputted value of a in fpcurf.py and simply initialize a new empty / zeros array.

fpcurf.py

def fpcurfpy(... a, ...):
    ...
    a = np.zeros(shape=(nest, k1), dtype=np.float64)
    # alt.
    # a = np.empty(shape=(nest, k1), dtype=np.float64)
    ...

Is this python code a correct representation of the fortran code? Assume a is not returned from fpcurf.f , i.e. not used in curfit.f

More specifically: is a = np.zeros(shape=(nest, k1), dtype=np.float64) in my python file the best representation of what's going on in the fortran code, especially how a is treated in fpcurf.f

It will be helpful for you to show complete code (see [mre]), otherwise there's a fair amount of guessing and noise. For example, you say `wrk` is "empty initialized" but you could show more convincingly what that actually means. No assignment happening is far from the same thing as "all values initially zero". (The code you do show is not Fortran 77.) — francescalus, Aug 07 '22 at 23:58
You left everything out you didn't think would matter - but you also have trouble understanding why it's not working as you expect. A likely cause is that something you left out does in fact matter. This is why you need to provide *all* of an example that shows the problem you're having, not just the parts you think matter - if you weren't wrong about something, you wouldn't be having this problem after all. — Grismar, Aug 08 '22 at 00:08
@francescalus it was my understanding that this is fortran 77.. if not, what is it? Also I said "empty initialized" because `real*8 ... wrk(lwrk)` is all that happens to the variable. Also added links to source code. — niCk cAMel, Aug 08 '22 at 06:30
@Grismar added links to source codes. I don't have trouble understanding because it's not working, the work is not finished and I can't know whether this is right or not. I do have trouble understanding whether this particular variable is ported correctly. Hope the links help. — niCk cAMel, Aug 08 '22 at 06:34
Source code in a link is the same as no source code. The code uses features of Fortran 90 or later (`recursive`, `implicit none`) and features that are not standard Fortran at all (`real*8`). There is no such thing as "empty initialized" in Fortran. It is undefined unless given some value. Is the array actually used anywhere at all? Where? How? — Vladimir F Героям слава, Aug 08 '22 at 07:09
@VladimirFГероямслава .. Fortran 90, got it, thanks. Ok, I should've said "undefined" instead of "empty". Oh yes, the the array is used inside `fpcurf.f` as a local variable. What should I do? Paste in 360 lines of code for `fpcurf.f` in the question? — niCk cAMel, Aug 08 '22 at 08:15
Then just make it a local variable. This way of reusing one big work array thatis passed around to save memory belongs to the 1980s, no reason to replicate it in a Python code. No, one should not dump hundreds of lines of code but it should be at least clear what happens to that variable. — Vladimir F Героям слава, Aug 08 '22 at 09:17
@VladimirFГероямслава it's not a question about locality. The question is whether `a = np.zeros(shape=(nest, k1), dtype=np.float64)` in my python file is the best representation of what's going on in the fortran codes, especially how `a` is treated in `fpcurf.f` — niCk cAMel, Aug 08 '22 at 09:22
You do not show how it is treated there! I trusted your description that it is used as a local variable. — Vladimir F Героям слава, Aug 08 '22 at 09:23
@VladimirFГероямслава by "treated" I meant, how it is initialized, given that it is inputted as an element of an array, but then initialized to a `(nest, k1)` array — niCk cAMel, Aug 08 '22 at 09:26
@VladimirFГероямслава imagine that was the entire code. `fpcurf.f` just does `real*8 a(nest,k1)` , makes a "hello world" print and exits. — niCk cAMel, Aug 08 '22 at 09:29
@VladimirFГероямслава I'm obviously simplifying the situation in order to get my question (found at end of post) answered.. Why you avoid it and have opinions about the usage, I can't understand. — niCk cAMel, Aug 08 '22 at 11:48
The question makes no sense at all. I posted my answer as the best description of the situation I could find. If there is nothing done to the array, there is no reason to do any specific `np.zeros`. If you need a work array inside your function, create it there. But in that case do not pass it as an argument. **Either** create it somewhere else and pass it as an argument **or** create it locally using `a=`. Doing both makes no sense. You will just redefine the reference. — Vladimir F Героям слава, Aug 08 '22 at 12:22
@VladimirFГероямслава since I know nothing about Fortran and the question was about representing a piece of simple fortran code in python, for me, the question does make sense. Of course I understand, if a variabel was passed back-and-forth between the functions, that would have an impact on the answer, but I clearly stated that it's not being returned, anyway.. your standpoint is clear, and I appreciate the time you've spent on this. — niCk cAMel, Aug 08 '22 at 14:24

Vladimir F Героям слава · Accepted Answer · 2022-08-08T20:27:05.167

In the previous century (except its last decade) it was common to re-use one large array for various local working spaces that various subroutines needed. One just created one array and passed it around so that the working space was shared and memory was saved.

In the current century you can use local working spaces and allocate it dynamically or have some object that contains various data structures that the algorithms needs.

Therefore, there is no need to pass around a. It should be created locally or in some class that represents the library or algorithm.

If you do pass a global working array around, do not redefine it locally, it will just redefine the reference and ignore what it received from above:

 def fun(a): 
     a = 5 
     print(a)

 a = []

 fun(a)                                                                                                                                                                                                                                                                                                          
5

In the original code the element wrk(ia) is passed. That is a legal way how to pass a part of the array that begins at element wrk(ia). In modern Fortran we would more likely pass wrk(ia:) instead. It is not clear why this is done. Perhaps to separate the section that is used in fpcurf from other sections used elsewhere.

The function interprets the passed array as a two-dimensional one a(nest, k1). This is again legal. Both of these actions are allowed by the sequence association rules which allow re-interpreting sequences of array elements. The receiving functions just see a part of the wrk array as a 2D array that is used for some internal work.

This way of passing around one large working array is archaic and should not be transferred to new code.

It would be great if you provided the python example in reference to the fortran code in OP, where `a` is an K-darray, an element of `a` is inputted to `fun`, and `a` is somehow (fortran-magic) turned into an NxM darray inside `fun` — niCk cAMel, Aug 08 '22 at 15:07
Great... So basically the 1d array is reshaped to (nest, k1) d array. Got it thanks. Then the answer is "yes, the python code is a good representation of what's happening in the Fortran code" (my interpretation). Thanks buddy — niCk cAMel, Aug 09 '22 at 05:39
@niCkcAMel Well, yes, it is a reshape. The Python code is will ultimately do the same computation, but it works differently with the argument. It actually just uses some local array (named `a`) which completely replaces the argument `a` which becomes unused. Basically, if you completely remove `a` from the `def fpcurfpy(... a, ...):`, the code will remain doing the same. — Vladimir F Героям слава, Aug 09 '22 at 06:17
"it is a reshape".. this is the answer! I know python all too well to realize I can skip having `a` as an input. But I had to keep it "stupid" since I don't know ANY fortran at all. — niCk cAMel, Aug 09 '22 at 08:31

Is this a correct python port of a fortran function call

1 Answers1