2-5  VARIABLE SIZE ARRAYS 
 *************************

 FORTRAN array sizes are determined at the time of compilation, 
 the 'top declaration' of the array allocates a chunk of memory, 
 and the size of that memory chunk can't be changed during runtime. 

 Sometimes we don't know beforehand the array size, because the size 
 of the input data is unpredictable or inconstant.

 There are some methods to solve this problem in the framework of 
 FORTRAN 77, two of them standard conforming and somewhat limited, 
 the others are non-standard. 

 Fortran 90 support of automatic and allocateable arrays solved 
 this problem. 


 Standard-conforming methods
 ---------------------------
 The idea behind these methods is simple: declare in the calling 
 procedure an array larger than the maximum expected size and 
 use only part of it in the called procedure. 

 There are two options now:

    1) Pass the 'physical' and 'logical' array dimensions, 
       ('physical' means here the original dimensions, and 
        'logical' are the reduced dimensions). 
       Declare the array with the 'physical' dimensions and 
       Use the dimensional variables to explicitly limit all 
       LOOP CONTROL VARIABLEs to the 'logical' dimensions.

    2) Re-configure the array in the called subroutine to have
       the required size. Be careful, after once writing to a 
       re-configured array, all the procedures referencing it
       must re-configure it to the same dimensions (or compute
       manually all indices).

 By the way, re-configuration of the dimensions of an array passed 
 to a procedure, is valid according to the FORTRAN 77 standard. 

 From the point of view of memory utilization (minimization of 
 memory paging, at least for multi-dimensional arrays) it is 
 preferred to redefine the array dimensions. In this method we 
 use the original array as a kind of memory pool only.

 Both options use the syntax of adjustable arrays.



 Example of the physical/logical method
 ======================================

      PROGRAM PHYLOG
C     ------------------------------------------------------------------
      INTEGER     NPHYS, NLOG, I
      PARAMETER   (NPHYS = 100)
      REAL        A(NPHYS)
C     ------------------------------------------------------------------
      DO I = 1, NPHYS
        A(I) = I
      ENDDO
C     ------------------------------------------------------------------
      WRITE(*,*) 'Enter array size: '
      READ(*,*) NLOG
      IF (NLOG .GT. NPHYS) STOP 'Too large passed array '
C     ------------------------------------------------------------------
      CALL SUB(A, NPHYS, NLOG)
C     ------------------------------------------------------------------
      END


      SUBROUTINE SUB(A, NPHYS, NLOG)
C     ------------------------------------------------------------------
      INTEGER    NPHYS, NLOG, I
      REAL       A(NPHYS)
C     ------------------------------------------------------------------
      DO I = 1, NLOG
        WRITE (*,*) I, A(I)
      ENDDO
C     ------------------------------------------------------------------
      RETURN
      END




 Example of the re-configuration method
 ======================================

      PROGRAM RECNFG
C     ------------------------------------------------------------------
      INTEGER     NPHYS, NLOG, I
      PARAMETER   (NPHYS = 100)
      REAL        A(NPHYS)
C     ------------------------------------------------------------------
      DO I = 1, NPHYS
        A(I) = I
      ENDDO
C     ------------------------------------------------------------------
      WRITE(*,*) 'Enter array size: '
      READ(*,*) NLOG
      IF (NLOG .GT. NPHYS) STOP 'Too large passed array '
C     ------------------------------------------------------------------
      CALL SUB(A, NLOG)
C     ------------------------------------------------------------------
      END


      SUBROUTINE SUB(A, N)
C     ------------------------------------------------------------------
      INTEGER    N, I
      REAL       A(N)
C     ------------------------------------------------------------------
      DO I = 1, N
        WRITE (*,*) I, A(I)
      ENDDO
C     ------------------------------------------------------------------
      RETURN
      END


 A more fancy example for the reconfiguration technique: 
 =======================================================

      PROGRAM ADJARR
C     ------------------------------------------------------------------
      INTEGER
     *          I, J,
     *          M, N,
     *          A(9,9)
C     ------------------------------------------------------------------
      DO I=1,9
        DO J=1,9
          A(I,J) = 10*I + J
        ENDDO
      ENDDO
C     ------------------------------------------------------------------
      WRITE(*,*) ' '
      WRITE(*,*) '  We have an integer array of size 9X9 '
      WRITE(*,*) '  Array(i,j) = (10 * i) + j            '
      WRITE(*,*) ' '
      WRITE(*,*) '  We will pass it to a subroutine with '
      WRITE(*,*) '  the adjustable array mechanism.      '
      WRITE(*,*) '  The array dimensions will be         '
      WRITE(*,*) '  RE-CONFIGURED to  M x N              '
      WRITE(*,*) ' '
      WRITE(*,*) '  Please supply M,N the new            '
      WRITE(*,*) '  dimensions in the range [1-9]        '
C     ------------------------------------------------------------------
      WRITE(*,*) '  '
      WRITE(*,*) 'Enter M: '
      READ(*,*) M
      WRITE(*,*) 'Enter N: '
      READ(*,*) N
C     ------------------------------------------------------------------
      CALL SUB(M,N,A)
      END


      SUBROUTINE SUB(M,N,A)
C     ------------------------------------------------------------------
      INTEGER    M, N, A(M,N)
C     ------------------------------------------------------------------
      WRITE(*,*) 
      WRITE(*,*) '  The re-configured array is: '
      WRITE(*,*)
      WRITE(*,'(1X,I4)') A 
C     ------------------------------------------------------------------
      RETURN
      END

 The variable format used in the last WRITE statement is non-standard.



 Non-standard methods
 --------------------
 All these solutions are based on the ability of most operating 
 systems to allocate a contiguous area of memory to a running 
 program upon request. Operating systems supply non-standard 
 (from FORTRAN point of view) 'system calls' that a program 
 can utilize to request a memory chunk while running.

 A better solution is to use the memory allocating functions 
 of C, of course they are implemented with those system calls, 
 but the C interface is standard (C standard).

 The following examples will work when certain requirements
 are met, we will try to make them explicit as possible.


 The classical 'non-contiguous array' trick 
 ------------------------------------------
 The main-program allocates the memory for the 'client procedure', 
 where a few simple manipulation are performed on the newly created
 array, just to prove the method works. Passing the array element 
 A(OFFSET) to a subroutine, make it easier to access the 'array extension'. 

 We declare an array A with one element, the 'mem' routine allocates
 memory for 100 elements the size of each is 4 bytes. 'mem' returns 
 an 'array index' OFFSET, computed with 'pointer arithmetic', that 
 integer is the 'start index' of the 'array extension'. The routine 
 'unmem', of course, deallocates the new memory area.

 Here we probably assume a call-by-reference argument passing
 mechanism, and our INTEGERs are 4 bytes long. 


      PROGRAM DYNARR
C     ------------------------------------------------------------------
      INTEGER*4
     *          A(1), OFFSET
C     ------------------------------------------------------------------
      CALL mem(A, 100, 4, OFFSET)
      CALL SUB(A(OFFSET))
      CALL unmem(A, 100, 4, OFFSET)
C     ------------------------------------------------------------------
      END


      SUBROUTINE SUB(A)
C     ------------------------------------------------------------------
      INTEGER
     *              A(100),
     *              I
C     ------------------------------------------------------------------
      DO I = 1, 100
       A(I) = I
      ENDDO
C     ------------------------------------------------------------------
      WRITE(*,*)
      WRITE(*,*) A
C     ------------------------------------------------------------------
      RETURN
      END


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

void mem(int *a, int *m, int *n, int *b)
    {
    int     *new;

    new = (int *)calloc(*m, *n);

    if (new == NULL)
        {
        printf("\n not enough memory \n");
        exit(1);
        }
    else
        *b = (int)(new - a) + 1;
    }


void unmem(int *a, int *m, int *n, int *b)
    {
    free(a + *b - 1);
    }



 Another example, a 5x5 two-dimensional INTEGER*2 array is allocated 
 and initialized in a C routine:

      program test
      integer		m, n, offset
      Integer*2		matrix(1)
      call mem(matrix, offset)
      m = 5
      n = 5
      call sub(m, n, matrix(offset))
      end


 A FORTRAN routine (just for testing):

      subroutine sub(m, n, matrix)
      integer		m,n, i
      integer*2		matrix(m,n)
      write (*,*)
      write (*,*) ' In Fortran: '
      write (*,*)
      write (*,'(1X,5I11)') ((matrix(i,j), j=1,5), i=1,5)
      return
      end


 The C routine can't modify the value of the pointer to "matrix", 
 so we use the non-contiguous array technique, and compute the 
 "offset" between a one-element array declared in the Fortran 
 code and the new allocated array.

 We use the allocated array in sub-procedures of the procedure 
 that called the C routine, by passing to them "matrix(offset)" 
 which is a "pointer" to the new array computed by f77.

 The whole method can be easily made flexible, and you can have 
 different size for the array on each run (with the adjustable 
 array syntax).


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

mem(short **matrix, long *offset)
{
	long	m, i, j;
	short	*ptr, **aux;

	ptr = (short *) malloc(25 * sizeof(short));
	*offset = (long)(ptr - (short *)matrix + 1);

	aux = (short **) malloc(5 * sizeof(short *));
	for(m = 0 ; m < 5 ; m++)
    		*(aux + m) = (short *)(ptr + (5 * m));

	printf("\n In the C routine: \n");

        for(i = 0 ; i < 5 ; i++)
		{
                printf("\n");
                for(j = 0 ; j < 5 ; j++)
                        {
                        aux[i][j] = (i+1) * 10 + (j+1);
                        printf("%11d", aux[i][j]);
                        }
		}

	printf("\n");

        return;
}


 Note that the C routine allocate the storage actually used in the 
 Fortran code in the first "malloc", it computes the offset between 
 the one-element Fortran array and the new array with a bit of 
 pointer arithmetic and passes it back to the Fortran code.

 The output looks like this:

 In the C routine: 

         11         12         13         14         15
         21         22         23         24         25
         31         32         33         34         35
         41         42         43         44         45
         51         52         53         54         55

 In Fortran: 

         11         21         31         41         51
         12         22         32         42         52
         13         23         33         43         53
         14         24         34         44         54
         15         25         35         45         55

 Note that printing "row by row" produces
 transposed results!


 Some remarks on the C routines used above
 -----------------------------------------
 On most UNIX machines, you will have to add a trailing underscore to 
 the routine names (mem_, unmem_). See the chapter on FORTRAN/C 
 interface, or just add the underscore if the compiler/linker complains.

 Old C compilers (non-ANSI), declare the dummy (formal) arguments types
 not inside the argument list, but in a separate line just below it.


 Explicitly controlling the parameter-passing-mechanisms 
 -------------------------------------------------------
 Here we probably assume a call-by-reference argument passing
 mechanism, and: 

   1) The existence of language extensions that control 
      the mechanism of parameter-passing (In this example 
      DEC fortran extensions are assumed).

   2) The pointer returned by 'malloc' can be stored
      in an INTEGER variable and then passed intact
      to a FORTRAN procedure (valid in most machines?).


      PROGRAM DYNARR
C     ------------------------------------------------------------------
      INTEGER
     *			m, n,
     *			ptr
C     ------------------------------------------------------------------
      WRITE(*,*) ' Enter two integers: '
      READ(*,*) m, n
C     ------------------------------------------------------------------
      ptr = malloc(%VAL(m*n))
      IF (ptr .EQ. 0) STOP 'Not enough memory '
C     ------------------------------------------------------------------
      CALL SUB(%VAL(ptr), m, n)
C     ------------------------------------------------------------------
      END


      SUBROUTINE SUB(arr, m, n)
C     ------------------------------------------------------------------
      INTEGER
     *			m, n,
     *			arr(m,n),
     *			i, j
C     ------------------------------------------------------------------
      DO j = 1, n
        DO i = 1, m
          arr(i,j) = 10*i + j
        ENDDO
      ENDDO
C     ------------------------------------------------------------------
      WRITE(*, '(1X,I5)') arr
C     ------------------------------------------------------------------
      RETURN
      END


 Many compilers support Cray pointers, these beasts may be used
 for memory allocation, but be careful, using Cray pointers for
 in computation may cause erroneous results (when automatic
 compiler optimization is turned on - the default).


Return to contents page