Can't find declaration of local array.

Hi, I'm new to FORTRAN and struggling with the following subroutine. There's an array 'INDU' used in the subroutine 'INVERS'. I can't find the declaration of INDU in the code. I searched the whole programm but INDU only appears in subroutine 'INVERS', so it can't be a global variable. The FTN95 debugger says 'Variable Unknown' when try to read the elements of the array. Nevertheless, the array INDU is used in the subroutine. Please help me to understand whats going on here.

SUBROUTINE INVERS(NGL,NUGL,MAX,MAXUGL,ME0,DIAG,SPG,LUGL,MAT,
     *                 VEK,KOL,END,X,MAT0,VEK0,END0)

      IMPLICIT REAL*8(A-H,O-Z)
C
      COMMON /WINDOW/ WIN3
      COMMON /IODEF / VERSIO,NTEST,IKART,NLIST,ND1,INET,NNET,NPLOT,
     ,                NUEPL
C
      REAL*8    DIAG(1),MAT(1),VEK(1),X(1)
      REAL*8    MAT0(1),VEK0(1),SPG(1)
      LOGICAL   LUGL(1)
C
      INTEGER   KOL(1),END(1),END0(1)
      CHARACTER*10 TE

C Where is INDU declared??? 
      INDU(IZ,IS)=(IZ-1)*NUGL+IS

C .... additional code where INDU is used

Indu is not an array; it is what is called a 'statement function'. There is less confusion when statement functions take other than integer arguments.

INDU seems intended to map a subscript pair to a single offset into an array.

I'm not sure whether INDU is a function...

See the line IMPLICIT REAL8(A-H,O-Z). This means that all variables beginning with 'A', 'B', ... 'H' or 'O' ... 'Z' are declared as real8. So I think that INDU is declared as a real*8 variable.

But I don't like IMPLICIT declarations and always use IMPLICIT NONE instead. This forces the programmer to explicitely declare all of his variables and makes the code clearer.

Besides: IZ seems to be the line index, IS the column index in a matrix, NUGL the number of unknows / equations in an equation system (= number of lines of the matrix). So, INDU calculates a unique index for each position in the matrix. INDU itself is a matrix. Find out the maximum dimensions of a matrix in your programme (say 30 x 30), then declare INDU in this way:

integer*4  INDU(30,30)

The code is written by a German-speaking author:

Z = Zeile (row) S = Spalte (column) GL = Gleichung (equation)

Wilfried

But INDU begins with 'I' and I [u:d5f047af2a]isn't[/u:d5f047af2a] in the range A-H or O-Z.

So INDU isn't real*8.

INDU is Integer.

It is pretty clear that INDU is a statement function in the code posted.

David, you're right, of course INDU was implicitely declared as integer. Sorry.

But - Tenobaal88 wrote that INDU appears [u:cb41286b36]only[/u:cb41286b36] in this subroutine but nowhere else in the programme. Therefore I guess that it is not a function but a matrix 😉

Have a nice day Wilfried

Go back to mecej4's post

Indu is what is called a 'statement function'.

These are obsolete but should still work. It uses the 2 arguments IZ and IS and also the variable NUGL.

I would change the declarations to as below, to indicate they are arrays of unknown size:

  REAL*8    DIAG(*),MAT(*),VEK(*),X(*) 
  REAL*8    MAT0(*),VEK0(*),SPG(*) 
  LOGICAL   LUGL(*) 

C INTEGER KOL(),END(),END0(*)

FTN95 provides a way of clarifying whether an identifier represents an array or a function, with the /XREF and /LIS options. If the source is compiled with the /XREF option, a file with the suffix .XRF is produced, which will contain a line such as INTEGER :: FUNCTION INDU

As John Campbell has noted, the use of '1' as a placeholder in dimension declarations for array type dummy arguments can be problematic. If array subscript checking is desired for the code in the subroutine, one has to put in actual limits, which can be constants or expressions involving integer dummy arguments, instead of '*'.

Obsolete or not, statement functions still work, and that is definitely what INDU is. If you haven't read your code for a long time, or you are unaccustomed to statement functions, they look like something else, and cause the puzzlement that this has caused. My understanding is that they are compiled inline, and don't create the overheads that conventional subroutines and functions do - but then FTNxx is a law unto itself. In Fortran77 and earlier, the scope of any variable name is limited to the subprogram it is declared in and subprogram names have a global scope, but statement functions have a local scope and also have access to variables such as NUGL that for a conventional function would need to be passed to it as a parameter or via COMMON.

A more modern approach would be to declare that INVERS 'contains' INDU (and then write a conventional but 'contained' function), although the statement function does have the advantage of brevity even if it suffers the disadvantage of incomprehensibility.

Passing arrays of indeterminate length used to be done with the '1', because only the address of the first element is passed, and not its length, in most early Fortran variants and probably the same now. This has been addressed in Fortran 90 and onwards.

There will be divided opinions on IMPLICIT typing. I like it because I instantly know the types of all my variables. The two times in nearly 50 years I strayed from this, it caused me huge trouble. I know it marks me out as a dinosaur, or a museum piece, but it works for me. Others may prefer IMPLICIT NONE, and that is their prerogative. But I can tell you what will cause indigestion, and that is to mix implicit and explicit typing (as here)!

Eddie

You're right Wilfried Linder. The Subroutine INVERS solves a system of linear equations (Gauss algorithm). INDU seems to be an implicit declaration of an array of integer. When I delete the line 'INDU(IZ,IS)=(IZ-1)*NUGL+IS' the debugger terminates the programm with the message 'Call to missing routine - INDU'.

I used the Notepad++ and other stuff to find the declaration of 'INDU' but there's nothing outside of subroutine 'INVERS'.

Edit: I turned on the XREF and LIST Option. I found a Function INDU in the XREF File.

Output from XREF File:

INTEGER :: FUNCTION INDU
  6485*,   6539,    6544,    6550,    6598,    6611,    6617,  
  6618,    6623  

Corresponding lines from the LIS File.

   6483   C
   6484         INDU(IZ,IS)=(IZ-1)*NUGL+IS
   6485  

   6538         DO 41 K=1,NUGL
   6539      41 SPG(INDU(N,K))=0
   6540         DO 40 K=K1,K2 

   6542           IF(KOL(K).GT.NGL) THEN
   6543              IUGL=KOL(K)-NGL
   6544              SPG(INDU(N,IUGL))=MAT(K)
   6545           ENDIF

   6548         IF(N.GT.NGL) SPG(INDU(N,N-NGL))=DIAG(N)
   6549   C
   6550   C***** BESETZEN DER RECHTEN SEITE

   6597   C
   6598         R=-SPG(INDU(IZ,IS))/SPG(INDU(IS+NGL,IS))
   6599         VEK0(IZ)=VEK0(IZ)+VEK0(IS+NGL)*R

   6610         IF(NUGL.GT.0) THEN
   6611          X(NGL+NUGL) = VEK0(NGL+NUGL)/SPG(INDU(NGL+NUGL,NUGL))
   6612         ENDIF

   6617      50 X(IZ)=X(IZ)-SPG(INDU(IZ,IS))*X(IS+NGL)
   6618         X(IZ)=X(IZ)/SPG(INDU(IZ,IZ-NGL))
   6619      49 CONTINUE  
   6620         DO 51 IZ=NGL,1,-1
   6621         X(IZ)=VEK0(IZ)
   6622         DO 52 IS=1,NUGL
   6623         X(IZ)=X(IZ)-SPG(INDU(IZ,IS))*X(IS+NGL)

Can somebody please explain what this means? What is the function actually doing?

Edit2: OK I understand whats going on. INDU(IZ,IS) is a function, IZ and IS are parameters of that function. The function calculates '(IZ-1)*NUGL+IS' and returns the result as an integer value. 😮ops: So the line 'INDU(IZ,IS)=(IZ-1)*NUGL+IS' defines a whole function...crazy 😒hock: Thanks for your help!

Quoted from Tenobaal88

Can somebody please explain what this means? What is the function actually doing?

The statement function is:

INDU(IZ,IS)=(IZ-1)*NUGL+IS

when it is called in your code like this:

SPG(INDU(N,IUGL))=MAT(K)

The compiler replaces this by the following:

SPG((N-1)*NUGL+IUGL)=MAT(K)

So the statement function acts like a macro in C if you know C.

Statement functions are a bit limited since they are confined to one statement (one line plus any continuation lines). The modern way is to write an 'internal functon' which allows you to use multiple statements.

Quoted from mecej4 If array subscript checking is desired for the code in the subroutine, one has to put in actual limits, which can be constants or expressions involving integer dummy arguments, instead of '*'.

Doesn't Silverfrost's FTN95 do subscript checking even if * is used on dummy arguments?

Good question!

If all the sources are compiled with /check, FTN95 will certainly flag the subscript error in a subroutine.

If the subroutine containing a subscript error contains '*' as the dimension of a dummy argument, the subroutine relies on the caller to pass information on the actual array size to it. If the caller is not compiled with /check but the subroutine is, then crashes or bad results can be seen in the subroutine even when the subroutine is compiled with /check. Here is an example: the subroutine is

subroutine sub(a,n)
integer a(*),n
integer i
do i=1,n*n+5
   a(i)=2*i-3
end do
return
end

The caller contains

program blowbnds
implicit none
integer n,a(10)
n=3
call sub(a,n)
write(*,*)a(3)
end program

If I compile only the subroutine with /check and run, no error is detected. If I change 'a(*)' to 'a(n)' in the subroutine, and compile in the same way, the error is reported in the subroutine.

A related issue is about using UBOUND, etc., on a dummy array argument. With assumed size arguments, one can't apply UBOUND.

I think only Silverfrost, NAG and Lahey compilers can do this now.

Related to this, in the following there should [u:b2be17e8ba]not[/u:b2be17e8ba] be a run time error if the code is compiled with /CHECK and /OLD_ARRAYS. I have not tried it though 😃

subroutine add(a, n)
   ! Need /OLD_ARRAYS so a(1) is interpreted as (*)
   real a(1)
   print *, 'sum is ', sum(a(1:n))
end

program main
   real b(5)
   b = (/1.0, 2.0, 3.0, 4.0, 5.0/)
   call add(b, 5)
end

Generally, as all my subprograms have an explicit interface I don't much use (*) in my own work, but do come across it a lot in other peoples code.

Edit: Interesting I think there's a bug here in the compiler. Will post separately in support area.

not[/u:b2be17e8ba] be a run time error if the code is compiled with /CHECK and /OLD_ARRAYS. I have not tried it though 😃

subroutine add(a, n)
   ! Need /OLD_ARRAYS so a(1) is interpreted as (*)
   real a(1)
   print *, 'sum is ', sum(a(1:n))
end

program main
   real b(5)
   b = (/1.0, 2.0, 3.0, 4.0, 5.0/)
   call add(b, 5)
end

Generally, as all my subprograms have an explicit interface I don't much use (*) in my own work, but do come across it a lot in other peoples code.

Edit: Interesting I think there's a bug here in the compiler. Will post separately in support area.

[quote:b2be17e8ba="mecej4"] With assumed size arguments, one can't apply UBOUND.

Well yes and no, you can use UBOUND (and SIZE etc) on all ranks except for the last one which has the (*). You can't use it on rank 1 assumed size arrays at all.

subroutine report(a, n)
   real a(n,*)
   print *, ubound(a,1)
end

program main
   real b(5,2)
   b = 1.0
   call report(b, 5)
end

I do recall that a very early version of FTN95 did not support statement functions, or at least had problems with them. Since converting to F95 I have removed them from my code.

I would not use statement functions in an equation solver as shown above. The preference now would be to have the array defined as 2D.

Early fortran coding would use a 1D array for efficiency, using smarts for the array address, rather than letting the compiler do the 2D address calculation. I think replacing that with an in-line function has lost the plot.

John

John,

I think it unlikely that statement functions didn't work with FTN95, as they've been part of Fortran for more than a half century, and they certainly worked with FTN77 - unless it was one of those things that Paul calls a 'regression'.

Whether an inline function is a good idea - or any other type of subprogram, or writing the code explicitly - is a function of the compiler and the computer. Whether it is good style depends on the reader. If the origins of the code go back to 'the dawn of time' other factors may influence this, e.g. one is dissuaded from using comments if there is a limit on the number of cards in a job (e.g. Imperial College's 'instant turnround' jobs on the CDC6400 c. 1972) or how many cards could read reliably. Moreover, such limitations militate against using white space and/or indents. I have one arithmetic IF left over from programming on an IBM1130 which still works perfectly - that computer had no logical IF!

In my view, one is least likely to find a problem in any facility that has been in a compiler for a long time.

Eddie

You asked what the function actually does. It seems to me that if you were to pass a two dimensional array to the subroutine, in which it is dimensioned using the old Fortran 66 method of dimensioning it to (1). Array bounds checking did not happen in those days so just declaring it with one element was enough to locate the start of the array in memory. In the subroutine, it becomes a one dimensional array. Outside the subroutine, one of the dimensions was equal to NUGL. This function determines the position in the one dimensional array that corresponds to element (IZ,IS). For example consider the two dimensional array outside the subroutine: a1 a2 a3 a4 b1 b2 b3 b4 c1 c2 c3 c4

Within the subroutine in one dimensional arrangement, it would look like:

a1 a2 a3 a4 b1 b2 b3 b4 c1 c2 c3 c4

and NUGL would be 4. Therefore to access element b3 for example, the two dimensional subscripts would be column=3, row=2. Evaluating the function would give: (2 - 1) * 4 + 3 = 7

and looking at the seventh element of the one dimensional array, we see that it is b3, so we have located it correctly.

Thus we can assume that NUGL is one of the dimensions of the array and is passed to the subroutine. No doubt the other dimension can be located by examining the range of IS in some part of the subroutine. The dimension statements could then use these dimensions to re-create the two dimensional arrays in the subroutine and all occurances of INDU(?,?) in an array subscript could be replaced with just the ?,?.

Ian

Ian Lambley wrote:

Within the subroutine in one dimensional arrangement, it would look like:

a1 a2 a3 a4 b1 b2 b3 b4 c1 c2 c3 c4

In Fortran, the equivalent one-dimensional array would actually be

a1 b1 c1 a2 b2 c2 a3 b3 c3 a4 b4 c4

We are used to reading matrices row by row. Fortran stores matrices by columns, i.e., all elements in col-1 are stored in sequence, then col-2, etc.

C and other languages use the convention of storing by rows.

One man's row is another man's column!

Using Excel notation, I would agree with Ian, but if the array is defined as (a:c,1:4) then I would agree with mecej4.

Unfortunately the picture Ian painted was not of excel layout.

I've always thought this inconsistency is over-rated. It is just a matter of observing the correct subscript order, which is left to right in Fortran and I presume right to left in C.

After all, there is no concept of 'row' and 'column' in Fortran, which makes reference to these attributes somewhat difficult. When I first learnt matrix mathematics, I did learn about rows and columns but was not aware of Fortran. It's a bit like considering if memory addresses count up, down, left or right. Then there are the 4 bytes of an INTEGER*4, which count to the left !!

John

Actually C doesn't have 2-dimensional arrays at all. You can only mimic 2D arrays using arrays of arrays. It is then a matter of choice whether you represent these as 'arrays of rows' or 'arrays of columns'.

Conventionally the elements of a matrix A (a mathematical object) with row index r = 1 ... M and column index c = 1 ... N are references in C as A[r][c] meaning it is convenient to think of arrays of rows (read indices right-to-left).

However, if you reference elements as A[c][r] it is convenient to think of the representation of arrays of columns. The latter notation has the same layout in memory as 2D arrays in Fortran.

So you could just as well say, C uses the same memory layout for 2D arrays as Fortran but requires the indices to be given in reverse order!

Although the following comments are about C, and therefore a little bit off-topic for this forum, since C is often used with Fortran it is worth noting them.

Quoted from davidb Actually C doesn't have 2-dimensional arrays at all. You can only mimic 2D arrays using arrays of arrays. It is then a matter of choice whether you represent these as 'arrays of rows' or 'arrays of columns'.

That is not quite correct. C does have 2-dimensional arrays, although there are some limitations to their usefulness. Here is a program to illustrate the differences between a native 2D array and a 2D array implemented through an array of row pointers. Note that the 2D array, c1, has 3 X 5 = 15 elements, and occupies exactly 60 bytes -- there is no space there to hold any pointers. The array of pointers version, on the other hand, occupies a total of 72 bytes (84 on a 64-bit target).

There are other differences. (i) For example, c1[2] is a const pointer, and cannot be changed; c2[2], on the other hand, is not defined until being assigned to the pointer returned by calloc(), and could be redefined later. (ii) The array c1 must be rectangular (or square, as a special case). An array of pointers variable, such as c2, can even be a 'ragged' array, with each row containing a different number of elements than the other rows, or with a first index different from 0 (even negative if desired), and this feature is useful in working with banded matrices.

The program outputs:

(32-bit): 60 12 72 (64-bit): 60 24 84

#include <stdio.h>
#include <stdlib.h>

int main(){
int c1[3][5], *c2[3];
int i,s,*p;

for(i=0, s=0; i<3; i++){
   c2[i]=(int *)calloc(5,sizeof(int));
   s+=sizeof(c2[i])+sizeof(int)*5;
   }

printf('%d %d %d\\n',sizeof(c1),sizeof(c2),s);
}