2

Is I mentioned in this post: "Is there locally saved data (like Fortran's) in Julia?" I am attempting to re-write an entire fortran77 code into Julia. This code is more than 2000 lines, so I won't be able to put it in here, but I will summarize the core of my problem.

In this old program, there are extensive use of functions from NUERICAL RECIPIES. In particular, ran1.f and ran3.f (both are random generator functions). The declaration lines in both functions are as follow:

FUNCTION ran1(idum)
      INTEGER idum,IA,IM,IQ,IR,NTAB,NDIV
      REAL ran1,AM,EPS,RNMX
      PARAMETER (IA=16807,IM=2147483647,AM=1./IM,IQ=127773,IR=2836,
     *NTAB=32,NDIV=1+(IM-1)/NTAB,EPS=1.2e-7,RNMX=1.-EPS)
      INTEGER j,k,iv(NTAB),iy
      SAVE iv,iy
      DATA iv /NTAB*0/, iy /0/

and

FUNCTION ran3(idum)
    INTEGER idum
    INTEGER MBIG,MSEED,MZ
C   REAL MBIG,MSEED,MZ
    REAL ran3,FAC
    PARAMETER (MBIG=1000000000,MSEED=161803398,MZ=0,FAC=1./MBIG)
C   PARAMETER (MBIG=4000000.,MSEED=1618033.,MZ=0.,FAC=1./MBIG)
    INTEGER i,iff,ii,inext,inextp,k
    INTEGER mj,mk,ma(55) 
C   REAL mj,mk,ma(55) Knuth.
    SAVE iff,inext,inextp,ma
    DATA iff /0/

In this mentioned post(above), I learned about the mutable structures in Julia that "apparently" let me simulate the expected behavior of the variables inside the SAVE and DATA statements. Below, my mutable structures

mutable struct ran1_init     # or simply "type" in older versions
    NTAB :: Int64
    idum :: Int64
    iv :: Array{Int64,1}
    iy :: Int64
    # constructors if necessary
end
# Initialize mutable structure for ran1
ran1_in = ran1_init(32,ISEED,zeros(Int64,32),0)

mutable struct ran3_init     # or simply "type" in older versions
    idum :: Int64
    iff :: Int64
    inext :: Int64
    inextp :: Int64
    ma :: Array{Int64,1}
    # constructors if necessary
end
# Initialize mutable structure for ran3
ran3_in = ran3_init(ISEED,0,0,0,zeros(Int64,55))

The fortran CODE uses ran1 until some point, where started to implement ran3 (I guess the program was assembled throughout the years, and in some point ran3 was implemented to generate random numbers, instead ran1). The first time ran3 is used in the program is in this operation:

    do while (iii.le.N_emean)
               t0=XA(1,1)+dble(ran3(ISEED))*(XB(1,1)-XA(1,1))
                d=XA(1,2)+dble(ran3(ISEED))*(XB(1,2)-XA(1,2))
               z0=XA(1,3)+dble(ran3(ISEED))*(XB(1,3)-XA(1,3))
                v=XA(1,4)+dble(ran3(ISEED))*(XB(1,4)-XA(1,4))
           iflag=0
stop

(Stopped for the purpose of showing this matter)

In Julia:

while (iii <= N_emean)
    for i = 1:4
       randomnumber = ran3(ran3_in)
       XX[i] = XA[1,i] + randomnumber*(XB[1,i]-XA[1,i]);
    end
    t0, d, z0, v = XX;
stop

At first, I thought it was working fine. But then I saw the following outputs in both Julia and FORTRAN CODE... the experiment bellow:

In Julia

println("initial parameters")
        println()
        println("ISEED = ",ran3_in.idum)
        println("iff = ",ran3_in.iff)
        println("inext = ",ran3_in.inext)
        println("inextp = ",ran3_in.inextp)
        while (iii <= N_emean)  

            for i = 1:4
            println()
            println("OUTSIDE: just before entering ran3")
            println("ISEED = ",ran3_in.idum)
            println("iff = ",ran3_in.iff)
            println("inext = ",ran3_in.inext)
            println("inextp = ",ran3_in.inextp)
            println()
            randomnumber = ran3(ran3_in)
            XX[i] = XA[1,i] + randomnumber*(XB[1,i]-XA[1,i]);
            println("OUTSIDE: just exiting ran3")
            println("ISEED = ",ran3_in.idum)
            println("iff = ",ran3_in.iff)
            println("inext = ",ran3_in.inext)
            println("inextp = ",ran3_in.inextp)
            println("random number : ", randomnumber)

                error("JUST PLAYING AROUND")

In Fortran:

write(*,*)"initial parameters"
write(*,*)
write(*,*)"ISEED = ",ISEED
write(*,*)"iff = ",iff
write(*,*)"inext = ",inext
write(*,*)"inextp = ",inextp
    do while (iii.le.N_emean)
        write(*,*)
        write(*,*)"OUTSIDE: just before entering ran3"
        write(*,*)"ISEED = ",ISEED
        write(*,*)"iff = ",iff
        write(*,*)"inext = ",inext
        write(*,*)"inextp = ",inextp
        write(*,*)
        a1ran = dble(ran3(ISEED))
        write(*,*)"OUTSIDE: just exiting ran3"
        write(*,*)"ISEED = ",ISEED
        write(*,*)"iff = ",iff
        write(*,*)"inext = ",inext
        write(*,*)"inextp = ",inextp
        write(*,*)"random number : ",a1ran
    stop

Getting: In Julia code

initial parameters

ISEED = 557527639
iff = 0
inext = 0
inextp = 0

OUTSIDE: just before entering ran3
ISEED = 557527639
iff = 0
inext = 0
inextp = 0

INSIDE ran3: just entering ran3
ISEED = 557527639
iff = 0
inext = 0
inextp = 0

INSIDE ran3: just exiting ran3
ISEED = 1
iff = 1
inext = 1
inextp = 32

OUTSIDE: just exiting ran3
ISEED = 1
iff = 1
inext = 1
inextp = 32
random number : 0.253305072

In fortran code

initial parameters

 ISEED =    557527639
 iff =            0
 inext =            0
 inextp =            0

 OUTSIDE: just before entering ran3
 ISEED =    557527639
 iff =            0
 inext =            0
 inextp =            0

 INSIDE ran3: just entering ran3
 ISEED =    557527639
 iff =            0
 inext =            0
 inextp =            0

 INSIDE ran3: just exiting ran3
 ISEED =            1
 iff =            1
 inext =            1
 inextp =           32
 OUTSIDE: just exiting ran3
 ISEED =            1
 iff =            0
 inext =            0
 inextp =            0
 random number :   0.25330507755279541

Seems like the SAVEed and DATAed BLOCKED variables are not exiting with the correct values in Fortran, while they do in Julia, 1) What is going on?. 2) Is there a problem between ran1 and ran3 for having the same named DATA and SAVE blocks. 3) Shouldn't be the Julia outputs the expected or corrected output?

I really need help understanding this...thanks

francescalus
  • 30,576
  • 16
  • 61
  • 96
GEBRU
  • 495
  • 2
  • 11
  • 1
    You seem to be expecting that (in the Fortran code) `iseed` outside `ran3` is the same thing as `iseed` local to `ran3`? Is that what you are expecting? And if so, why do you think that? – francescalus Feb 21 '18 at 20:37
  • 1
    What is the whole purpose of that? Why would you want to replicate the crappy Numerical Recipes random number generator? Doesn't Julia have a decent built-in random number generator or a PRNG library? – Vladimir F Героям слава Feb 21 '18 at 20:40
  • @francescalus On the contrary, ISEED is the only one that is behaving like I would expect, that is, it has the same value at the begining, and when modified inside ran3, it modifies in the ouputs. Instead, the variables "iff","inext" and "inextp" change inside ran3, but outside they appear to be allways " 0 ". I would expect that variables to be "shared" outside the program, but apparently, SAVE and DATA statement, only save and declare variables within ran3. So, I am misundestaning the meaning of SAVE and DATA...right? – GEBRU Feb 22 '18 at 13:59
  • @VladimirF, Indeed, Julia has very nice random intrinsic functions, BUT...I am traking some errors at the final output of the main program and, attempting to achieve the same results, I need to reproduce the same exact random numbers. Otherwise, I guess I can not tell if the error happens to be inside some particular subroutine, the random number subroutines, or a mixture of those! – GEBRU Feb 22 '18 at 13:59
  • Yes, I should have said "`iseed` _and the other variables_". `save` and `data` have nothing at all to do with "sharing". If you explain what you think `save` and `data` mean then we can perhaps address your misunderstandings. – francescalus Feb 22 '18 at 14:01
  • @francescalus, Agreed.. So, I understand that the DATA statement has the purpose of initializing the value of some variable, for instance, DATA iff /0/, declare that iff will be 0 "that very first time", but(as far as I understood) if it changes, it will no longer be that value when entering the subroutine again... Am I right? In the case of SAVE, I understand that those values are stored in some place in momory, and can be used throughtout the program at my will, with the last values saved in those variables... I am completly wrong I guess... – GEBRU Feb 22 '18 at 14:11

2 Answers2

3

Here I address the misconception around the Fortran code. I won't attempt to detail the implications for a Julia program - I'll leave that to other answers, such as this one.

Important to Fortran is scope. What I'll say about that is closely related to the implicit statement.

Consider two procedures:

subroutine x
  integer i
end subroutine x

subroutine y
  integer i
end subroutine y

Here there are two scopes and in each a variable i. These are two totally distinct variables. There is no overlap between them at all.

Now, neither the data statement nor the save statement change that. The save statement such as

subroutine x
  integer i
  save i
end subroutine

says that the variable has the save attribute. This means that the variable retains its value when the subroutine completes execution (without it the variable would become undefined).

The data statement, such as

subroutine y
  integer i
  data i/4/
end subroutine y

explicitly initializes the variable with a particular value. On first entering the subroutine y the variable i has value 4. Crucially, the data statement also gives the initialized variables the save attribute. So, if y's i has its value changed from 4 that new value will be retained on next entering the subroutine.

That's all there is, in this context, to the data and save statements. Nothing affecting scope.

francescalus
  • 30,576
  • 16
  • 61
  • 96
  • Common indeed likely to be a point (especially given the first question on this topic). I'm slightly worried that my "long comment" (rather than a finely tuned answer) is overbearing even before bringing in common blocks - especially when storage association doesn't happen by "magic" (would need to be reciprocated in the other scoping unit). – francescalus Feb 22 '18 at 15:01
  • The example above was really useful. I think I have now a pretty good idea of what the problem was!. Thanks a lot! – GEBRU Feb 23 '18 at 13:18
2

The Fortran save attribute, applies to local variables. Therefore they cannot be visible outside of the function. I do not have the full Fortran code, but to me the "outside variable" in the Fortran code, are just different variables, that happen to have the same names than the local variables. I would recommend to use implicit none in Fortran.

To have "saved" variables in julia (or static variables as they are called in C), I prefer the following approach:

let
      my_saved_variable = 0 # the value corresponds to the Fortran DATA statement 

      global function do_something(inc)
           my_saved_variable += inc
           return my_saved_variable
      end
 end

The parameter inc gets added to the saved variable:

You get the following output:

julia> do_something(2)
2

julia> do_something(5)
7

julia> do_something(10)
17
Alex338207
  • 1,825
  • 12
  • 16
  • Apparently, there are no IMPLICIT NONE declarations (but one in some non critical subroutine with very little variables) in all the main program and subroutines, and there is onother subroutine named "cost" that has the following declaration lines: SAVE iff, ax, ay, az, and DATA iff/0/ ... I am wondering if this can be an issue ... I will try to pinpoint this variables... Thanks a lot! – GEBRU Feb 22 '18 at 14:16