Requiring comprehension of equivalence to enable removal

Help
1

I have the following lines in my code…

            double precision, dimension(mj) :: temp, xl, r, pres
            equivalence     (x(1, 1), pres), (x(1, 4), temp)
            equivalence     (x(1, 2), r), (x(1, 3), xl)

I understand that equivalence causes the two variables to share common memory allocation, thus as one variable’s contents change, the other changes as well (since they have the same memory allocation).

However, I note that pres, temp, r and xl are arrays with a dimension of mj (set to 4000).

I am reading this as that every single part of those 4 arrays are linked to a single value in the corresponding index of x.

For example, all 4000 elements of pres are set to reflect the contents of x(1, 1). If x(1, 1)'s contents change, then every part of pres will change accordingly. Equally, if any part of pres changes, then all parts of pres change, as does x(1, 1)

  1. Am I correct here?
  2. If so, does this not appear to be overly complex, wasting of space and time?
  3. If not, where am I going wrong in my thinking?

Thank you for any help received.

2

@garynewport , that is not the case with EQUIVALENCE. The code you show is equivalent (pun intended!) to the following:

double precision, target :: x(mj,4)
double precision, pointer :: pres(:)
double precision, pointer :: r(:)
double precision, pointer :: xl(:)
double precision, pointer :: temp(:)
..
pres(1:mj) => x(1:mj,1)
r(1:mj) => x(1:mj,2)
xl(1:mj) => x(1:mj,3)
temp(1:mj) => x(1:mj,4)
..

So basically the pres vector of mj elements is mapped to the mj elements in the 1st column of your data array x, r to the next column, and so forth.

3

I think you should extract this code into a separate program and play around with it, I have an old copy of Monro’s book, and I have converted similar F77 code, but my understanding/memory is:

  1. Suppose x was declared x(4000,4) then
    x(1,1) and pres(1) are the same location in memory, say address ‘1’, then
    x(2,1) and pres(2) are both address ‘2’
    x(4000,1) and pres(4000) are both address ‘4000’ etc

But things get more complicated if, say, x is declared x(1000, 4)
x(1,1) and pres(1) are both at address 1
x(1000,1) and pres(1000) are at address 1000, but
x(1,2) and pres(1001) are at address 1001 which is now also r(1)

(This also means that completely impossible equivalence statements can be written. e.g.

equivalence(r, pres)

happy days!)

My advice is therefore to draw a diagram of memory and check it with some print statements

Good luck, it looks like you might need it depending on how x is declared.

4

Good heavens, no! Such questions are best resolved by reading language manuals and standards, not by letting one’s imagination run wild.

1 Like
5

@garynewport , you can create "minimal working example"s of sections of your legacy code and try them out to understand the behavior per FORTRAN 77 and other dialects/extensions upon which such legacy code may be based. So for example,

! legacy
   integer, parameter :: mj = 2
   double precision x(mj,4)
   double precision, dimension(mj) :: temp, xl, r, pres
   equivalence     (x(1, 1), pres), (x(1, 4), temp)
   equivalence     (x(1, 2), r), (x(1, 3), xl)
   integer :: i
   do i = 1, mj
      x(i,1) = 1D5*i
      x(i,2) = 1D0*i
      x(i,3) = 1D1*i
      x(i,4) = 300D0*i
   end do
   print *, "pres = ", pres
   print *, "r = ", r
   print *, "x1 = ", xl
   print *, "temp = ", temp
end

vis-a-vis

! EQUIVALENCE replaced by TARGET<->POINTER paired semantics starting Fortran 90 standard
   integer, parameter :: mj = 2
   double precision, target :: x(mj,4)
   double precision, pointer, dimension(:) :: temp, xl, r, pres
   integer :: i
   pres(1:mj) => x(1:mj,1)
   r(1:mj) => x(1:mj,2)
   xl(1:mj) => x(1:mj,3)
   temp(1:mj) => x(1:mj,4)
   do i = 1, mj
      x(i,1) = 1D5*i
      x(i,2) = 1D0*i
      x(i,3) = 1D1*i
      x(i,4) = 300D0*i
   end do
   print *, "pres = ", pres
   print *, "r = ", r
   print *, "x1 = ", xl
   print *, "temp = ", temp
end

You can then run both the programs and check whether the output is equivalent - it should be the same.

6

Yes, of course. When I read the passages I have on Fortran I remembered that the arrays can be of different dimensions but the size equates to the same; I just forgot this until I read your reply (scary, since I only read that part of the book about 2 days ago!).

Thank you for the clarity.

7

Creating multiple copies of the code, reducing it down and experimenting with those parts is exactly what I am doing (it might not appear so, from my queries, but I am constantly exploring and testing).

Thank you for the explanation and, yes, as I have said in another reply, this echoes the very thing I read - just forgot I had read it until reading through these. :slight_smile:

8

I’m more of the opinion in this case that I can remove the equivalence altogether. It appears that these are used within the local routine, rather than directly working with the globals (as defined within the common blocks; something I am also replacing with modules).

9

When I was reading through the code tonight, I came across this equivalence and could not, for the life of me, remember what I had read. I felt I needed the re-assurance from this community to allow me to move forward (something that has worked wonders so far).

Alongside this I have bookmarked sites and several books that I have at differing points, depending upon style of writing and my requirements at the time; this is not my only resource.

I can assure you that I am not allowing my imagination run wild but believe that sometimes it is worth asking the stupid questions, since the answers help secure understanding. I do not care if others think me stupid; that’s their issue, not mine.

10

In your original post, you had the rank-2 array x equivalenced to various variables, but did not show the declaration of x itself. I think that @FortranFan has guessed that declaration correctly.

The Equivalence statement is marked obsolescent, and its use in modern Fortran code can lead to several problems. It was useful in the early days of Fortran, when memory was a scarce resource, with a big mainframe coming equipped with as little as 32 or 64 kWords. On such a machine, Equivalence allowed memory to be shared by different variables at different times, as long as no overlapping usage patterns were needed.

Another use for Equivalence, the one that I believe is relevant in your example code, is for packing multiple arrays into a single container array. Thus, x(:,1:4) could be regarded as a composite rank-2 array whose first column holds pres, the last column holds temp, etc. The use of pres(i) rather than x(i,1) makes the code more readable.

There are several restrictions on the use of this obsolescent feature, and it is difficult to work with code that uses the feature unless the conventions, intentions and usage patterns are well-documented by the original author.

1 Like
11

Yes, @FortranFan was correct in that assumption, and I agree with your assumption that it’s use in all of this is to make the code a little easier to read/code.

I do have one more equivalence left that I am struggling to see how I should resolve this.

The variable that the equivalence appears to be resolving is found in an included file, called xvar.h…

      double precision :: ha(mh, 2 * mh + 1)
     
      common/  hy  /ha

mh is defined a little earlier, in another included file called parm.h…

        integer, parameter :: mh          = 4

So, this gives ha a dimension of ha(4, 9).

In the main code (stelcor.f), along with other subroutines, is a subroutine called henyey.

In henyey appears the following declarations…

            double precision :: d(mh, mh + 1), ca(mh)
            
            equivalence (ha(1, mh + 1), d(1, 1)), (d(1, mh + 1), ca(1))

Now, this makes for the following dimensions; d(4, 5) and ca(4).

If I am reading this correctly, then when I get d(i, j) = -he(i, j), this could be replaced with ha(j, mh + j). Am I correct there?

If so (or even if not), what would I replace this with… d(i, mh + 1) = -g(i)?

And how does ha() relate to ca()?

Basically, I need to lose d() and ca(), lose the equivalence and use the included ha() instead.

Sorry and thanks for any help. :slight_smile:

12

Your description is incomplete, and probably contains some typographical errors. The arrays He and g have not been declared, so I have no idea what they are. Without being able to see the entire code, one cannot answer your questions, but one may well ask, “why don’t you leave the code as it is, if it is working correctly?”

Conceptually, it seems to me that you have a number of matrices and vectors, all with MH rows. The composite matrix Ha = [H0 D] = [H0 [D0 Ca]], where H0 and D0 are square matrices and Ca is a column vector. Given all this, the assignment statement

d(i,j) = -He(i,j)

could be rewritten as

Ha(i,j+MH) = -He(i,j)

The question, however, is “for what purpose?”

13

Part of the reasoning is that I need to understand how the code works, because later I have got to adapt it. Secondly, there is a lot of time lost in the current format. Thirdly, some of the existing code has been deprecated and will, at some point, become unsupported; meaning that the task I am completing will have to be done anyway.

The alterations I have to make later are substantial and require that I get this code into a structure that I can not only easily follow but adapt.

In terms of loss, my working through the code has allowed me to identify subroutines that do not work properly (or not even included), loops that aren’t required, file calls (the thing that is slowing it down the most) that aren’t necessary, etc

The equivalences are using up space; another aspect that needs resolving. The use of a lot of variables runs throughout the code, with very little re-use of existing variables (so, you end up with x1, x2, x3, x4, etc; when x was perfectly capable of being re-used each time). I need to identify each variable (currently 1064 of them) and recognise what each one is used for. Of those 1064, only a few are likely to be related to the actual purpose of the code (the astrophysical element); whichare the parts I need to identify the most.

Another point is that the code works, but it uses understanding (astrophysics understanding) from before 2012. I need to bring this up to date AND be able to adapt it further once we begin to get new data.

Yes, sorry, I should have recognised that I had not defined the other variables I gave.

g is defined as double precision :: g(mh), whilst he is defined as double precision :: he(mh, mh).

14

A task of the type that you described requires a lot of resources and experience, in addition to subject knowledge. One would need a large number of test cases with certified results, a rather complete set of computer tools to transform, debug and verify modifications, and sufficient personnel allocated for the task.

How big is the code, and is the source code freely available? Is it similar to Garstec?

15

You overplay the complexity; the actual astrophysics part is relatively straight forward; it’s the implied ration that this code offers.
No, it is not open-source, hence why I haven’t shared too much of the code.
The data being outputted is detailed and easy to check against changes. Equally, the results of that data has already been checked against actual data and maintains consistency.

16

The code is around 7000 lines and appears similar to your example. I can’t see any obvious disconnects between the two, though mine was finalised in 2012 and this paper is from 2007.
It might simply come down to agreement between the original author and the university.