Community
cancel
Showing results for
Search instead for
Did you mean:
Beginner
694 Views

## Allocatable arrays of different type and EQUIVALENCE

Jump to solution

I work on a very old code which relies on two big working arrays IW and RW, which are equivalenced in the main program and then used directly or propagated to called subroutines like this:

``````      program test
implicit none
external inner
integer(4),parameter :: li = 100
integer(4) iw(li)
real(8)    rw(li/2)
equivalence (iw,rw)
rw = -10.20
print *,iw
print *,rw
!       call someroutine(iw,li,rw,lr) will also work
end
``````

When I replace the static arrays with dynamically allocated, then I of course cannot equivalence the arrays. I came up with this code:

``````      program test
implicit none
external inner,outer
integer(4),parameter :: li = 100
call outer(inner,li)
end
subroutine outer(proc,li)
implicit none
external proc
integer(4),intent(in)  :: li
integer(4),allocatable :: iw(:)
allocate(iw(li))
call proc(iw,li,iw,li/2)
deallocate(iw)
end
subroutine inner(iw,li,rw,lr)
implicit none
integer(4),intent(in)    :: li,lr
integer(4),intent(inout) :: iw(li)
real(8),intent(inout)    :: rw(lr)
rw = -10.20
print *,iw
print *,rw
!       call someroutine(iw,li,rw,lr) will also work
end
``````

So here the main program is reduced to a call to a helper subroutine outer, which allocates a single array and then passes it as two arrays to subroutine inner, which in turn contains the code of the original main program. The reason the compiler does not complain about the third parameter to inner not being a real(8) array is that in outer, inner is hidden in proc as an external subroutine with arbitrary arguments.

Is there any other, maybe more simpler or clean, way?

1 Solution
Honored Contributor I
349 Views

Also, if your interest in Fortran is limited to doing the least with this legacy codebase to get it to work with a larger problem size beyond the stack limits, you are essentially done by using the approach with the IW allocatable+target attributes, RW as pointer, and C_F_POINTER memory  mapping shown above.  No need to read any further, just "let the sleeping dogs lie" with the legacy code.

But if that ain't the case and you need to do more with Fortran, you may want to first start with a close review of a book such as this:

https://www.amazon.com/Fortran-Scientists-Engineers-Stephen-Chapman-ebook/dp/B06XCTY8KR

And pick up on modern Fortran concepts that can open up a world of options to write better new code and improve upon existing/legacy codebases for efficient and manageable applications in scientific and technical computing.

32 Replies
Beginner
187 Views

One last question, just to be sure. To solve the issue with the array argument, it is sufficient to add just a colon after the index to have a standard-compliant code? Like so:

CALL BAR(...,RW(LBFOO:),...)

Black Belt
180 Views

Isn't RW(...) argument a subsection of the blob?

As such, you would need to pass the lower bound (LBFOO) as well as the upper bound (UBFOO).

Using the RW(LBFOO:) would pass a slice starting at LBFOO and extending to the end of the blob.

Jim Dempsey

Beginner
173 Views

Remember that this is ancient code where arrays are always passed by their first element and length:

``````CALL FOO(RW(LBBAR),LBAR,IW(LBBAZ),LBAZ)
.
.
.
END
SUBROUTINE FOO(A,LA,IA,LIA)
DIMENSION A(LA),IA(LIA)
.
.
.``````

However, the moment RW becames POINTER, the compiler starts to complain about RW(LBBAR) not being an array.

Adding a colon in the argument RW(LBBAR:) makes the code compilable again. For all intents and purposes, this should have no effect at runtime.

Black Belt
164 Views

While you can use (LBBAR:) today, sometime in the future, you or your replacement may be updating the code and which at that time may introduce errors (or at least confusion). While adding ":" makes sense to you today, it might not make sense to someone else later.

So you have three choices:

Use ...,RW(LBBAR:),LBAR.... ! ":" indicates to end of blob, fix later

or use RW(LBBAR:UBBAR),LBAR,...! LBAR is passed for legacy purposes, fix later

or use RW(LBBAR:LBBAR+LBAR-1),LBAR, ... ! add UBBAR later, then remove LBAR on next revision.

Jim Dempsey

Black Belt
159 Views

When you revise the code (some time later), you might consider using Defining Generic Names for Procedures (intel.com)

This way you can have your subroutine FOO become a generic name with two interfaces: old/legacy argument list, and new updated argument list.

Jim Dempsey

Honored Contributor I
146 Views

Jim,

See my comment just prior to yours and the option therein - here's a link to it.

The use of additional variables of POINTER attribute for the work arrays of interest to OP - IA1, IA2, RA1, RA2, etc. - in this legacy code is more safe and readable, maintainable in my experience.

Black Belt
138 Views

Petr_Parik,

A word of caution about being overly zealous with allocating vary large array.

At time of allocation of a huge blob, the allocation may succeed in obtaining virtual addresses, but can fail some time later at first touch when it is found that the page file has insufficient space. This failure (page file full) is not trappable.

The is a caution not to request a few terabytes of memory for your buffer.

Jim Dempsey

Honored Contributor I
167 Views
``````@Petr_Parik wrote:
.. One last question, just to be sure. To solve the issue with the array argument, it is sufficient to add just a colon after the index to have a standard-compliant code? Like so:

CALL BAR(...,RW(LBFOO:),...)``````

While you can use an array section notation to solve the issue, note you may possibly face other issues and that's why I didn't suggest it earlier.  Note it's tough to list what problems you may run into what without knowing your entire code but it can be a performance hit in the form of an array temporary in CALLs to other fault scenarios during run-time.

If you're OPEN to the notion of making changes on the CALLER side like what you show with "CALL BAR(..,RW(..:), ..", then my suggestion will be to use the POINTER facility in modern Fortran to setup SPECIFIC objects with the POINTER attribute for EACH of your "work" ARRAYs.  See an illustration of this with your SAMPLE program where you will see 4 SPECIFIC objects for your 4 arrays - IA1, IA2, RA2, and RA2.  Given your legacy code with all the existing "book-keeping" toward the indices and array sizes, you'll see it's straightforward to do this and it WILL work with Intel Fortran, gfortran, and pretty much all the standard-conforming Fortran compilers out there.  I show this Intel Fortran and gfortran below.

``````      PROGRAM TEST_PTR
USE ISO_C_BINDING
IMPLICIT INTEGER(4) (I-N), REAL(8) (A-H,O-Z)
PARAMETER (LI = 100)
ALLOCATABLE IW(:)
TARGET IW
REAL(8), POINTER :: RW(:)
C Additional pointers for array argument to WORK* subprograms
INTEGER(4), POINTER :: IA1(:)
INTEGER(4), POINTER :: IA2(:)
REAL(8), POINTER :: RA1(:)
REAL(8), POINTER :: RA2(:)
ALLOCATE(IW(LI))
CALL C_F_POINTER(C_LOC(IW),RW,(/LI/2/)) ! RW => IW
LBIA1 =  1; LIA1 = 10
LBIA2 = 21; LIA2 = 10
LBRA1 =  6; LRA1 = 5
LBRA2 = 16; LRA2 = 5
C Set up array arguments
IA1(1:LIA1) => IW(LBIA1:)
IA2(1:LIA2) => IW(LBIA2:)
RA1(1:LRA1) => RW(LBRA1:)
RA2(1:LRA2) => RW(LBRA2:)
CALL WORK1(IA1,LIA1,IA2,LIA2,
&           RA1,LRA1,RA2,LRA2)
CALL WORK2(IA1,LIA1,IA2,LIA2,
&           RA1,LRA1,RA2,LRA2)
END
SUBROUTINE WORK1(IA1,LIA1,IA2,LIA2,RA1,LRA1,RA2,LRA2)
IMPLICIT INTEGER(4) (I-N), REAL(8) (A-H,O-Z)
DIMENSION IA1(LIA1),IA2(LIA2),RA1(LRA1),RA2(LRA2)
IA1 = 1
IA2 = 2
RA1 = .123456789
RA2 = .987654321
END
SUBROUTINE WORK2(IA1,LIA1,IA2,LIA2,RA1,LRA1,RA2,LRA2)
IMPLICIT INTEGER(4) (I-N), REAL(8) (A-H,O-Z)
DIMENSION IA1(LIA1),IA2(LIA2),RA1(LRA1),RA2(LRA2)
PRINT *,"IA1 ="
PRINT *,IA1
PRINT *,"IA2 ="
PRINT *,IA2
PRINT *,"RA1 ="
PRINT *,RA1
PRINT *,"RA2 ="
PRINT *,RA2
END
``````

Build and execute with Intel Fortran and gfortran:

``````C:\Temp>ifort test_ptr.for
Intel(R) Fortran Intel(R) 64 Compiler Classic for applications running on Intel(R) 64, Version 2021.1.2 Build 20201208_000000
Copyright (C) 1985-2020 Intel Corporation.  All rights reserved.

Microsoft (R) Incremental Linker Version 14.26.28806.0
Copyright (C) Microsoft Corporation.  All rights reserved.

-out:test_ptr.exe
-subsystem:console
test_ptr.obj

C:\Temp>test_ptr.exe
IA1 =
1           1           1           1           1           1
1           1           1           1
IA2 =
2           2           2           2           2           2
2           2           2           2
RA1 =
0.123456791043282       0.123456791043282       0.123456791043282
0.123456791043282       0.123456791043282
RA2 =
0.987654328346252       0.987654328346252       0.987654328346252
0.987654328346252       0.987654328346252

C:\Temp>gfortran test_ptr.for -o gcc-test_ptr.exe

C:\Temp>gcc-test_ptr.exe
IA1 =
1           1           1           1           1           1           1           1           1           1
IA2 =
2           2           2           2           2           2           2           2           2           2
RA1 =
0.12345679104328156       0.12345679104328156       0.12345679104328156       0.12345679104328156       0.12345679104328156
RA2 =
0.98765432834625244       0.98765432834625244       0.98765432834625244       0.98765432834625244       0.98765432834625244

C:\Temp>``````

Beginner
113 Views

I know sub-array pointers would work. Unfortunately, that is even more error-prone work than adding ":" to calls, as I have found. So I have settled with a slightly bloated but simple "loader" approach:

``````      PROGRAM LOADER
C     allocation/equivalence stuff ...
CALL FOO(IW,LI,RW,LI)
END
SUBROUTINE FOO(IW,LI,RW,LR)
DIMENSION IW(LI),RW(LR)
C     original PROGRAM FOO code with almost no changes
C     'RW' behaves like a normal array, no need to modify calls
C     ...
END``````

Adding ":" or upper bound would be too much error-prone work, see above.

However, interfaces for incremental refactoring are an excellent idea, thank you.

Regarding the virtual memory; allocating a large blob then failing later at runtime due OOM is no different to allocating arrays as needed and then still failing later at runtime when virtual memory is depleted. But fear not, LI is INTEGER(4) so we are talking 8 GiB blob max. Some new programs have 64-bit LI, but even in that case the maximum needed so far has been less than the physical memory (128 GiB). Generally, if I need to do computations, I need to make sure nobody else is using the computer.

Black Belt
98 Views

>> allocating a large blob then failing later at runtime due OOM is no different to allocating arrays as needed and then still failing later at runtime when virtual memory is depleted

With the exception that one can test for failure at ALLOCATE (assuming the heap isn't larger than remaining space in page file).

>> But fear not, LI is INTEGER(4) so we are talking 8 GiB blob max

From your description (which I may have misunderstood) is your code allocates once and anything becoming unused is lost. Presumably the allocations reach a stasis with a not unreasonable amount of dead space in the blob.

Jim Dempsey

Beginner
91 Views

>> With the exception that one can test for failure at ALLOCATE (assuming the heap isn't larger than remaining space in page file).

You are right, but if we distinguish only between a successful computation and a failed computation, it does not matter if the program fails gracefully.

>> From your description (which I may have misunderstood) is your code allocates once and anything becoming unused is lost. Presumably the allocations reach a stasis with a not unreasonable amount of dead space in the blob.

As you pointed out earlier, modern operating systems allocate only the memory accessed by the program. But yes, there can be several gigs of unused space in the blob.

Honored Contributor I
374 Views
``````@Petr_Parik wrote:

.. If I am not mistaken this problem is nonexistent in C where you can simply recast any variable to the desired type while passing it (by reference) to a function (subroutine).``````

To the extent you're using modern Fortran compilers such as Intel Fortran (IFORT) as part of Intel oneAPI or recent incarnations of GCC/gfortran, NAG compiler, etc. - all of which work with a companion C processor, you can take advantage of modern Fortran standard - Fortran 2018 - and its facilities to start introducing INCREMENTAL REFACTORING such as dynamic memory and associated management and build toward transforming such legacy codebase into a modern app.

Shown below is a simple prelude you can review - note the EQUIVALENCE statement has been replaced with a C processor based trick that maps the ALLOCATED memory to an array of floating-point type (REAL with kind of R8) with a POINTER attribute.  Then with suitable use of indices in this legacy codebase - assuming they are set up to keep the integer and real "data" of this program separate, you should be able to work this with the rest of your code.  Note the "editorial" liberty I've taken is to replace the non-standard "INTEGER*4" and "REAL*8" types and kinds with a tab more palatable options from the Fortran standard.

``````      PROGRAM SAMPLE
USE, INTRINSIC :: ISO_FORTRAN_ENV, ONLY : I4 => INT32,
&                                          R8 => REAL64,
&                                          STDIN => INPUT_UNIT,
&                                          STDOUT => OUTPUT_UNIT
USE, INTRINSIC :: ISO_C_BINDING, ONLY : C_LOC, C_F_POINTER
IMPLICIT INTEGER(KIND=I4)(I-N), REAL(KIND=R8)(A-H,O-Z)
INTEGER(I4), ALLOCATABLE, TARGET :: IW(:)
REAL(R8), POINTER :: RW(:)
WRITE(STDOUT,*) "Enter program size parameter value: LI = "
WRITE(STDOUT,*) "(for demo, enter an integer value >= 6)"
READ(STDIN,*) LI
ALLOCATE( IW(LI), STAT=ISTAT )
IF ( ISTAT /= 0 ) THEN
WRITE(STDOUT,*) "Error allocating memory: ISTAT=", ISTAT
STOP
END IF
IW = 0
CALL C_F_POINTER( CPTR=C_LOC(IW), FPTR=RW, SHAPE=[ LI/2 ] )
C Assuming program logic with indexes ensures integer and real data
C do not overlap
LBEGIW = 1 ; LENDIW = 2
LBEGRW = 2 ; LENDRW = 3
IW(LBEGIW) = 42 ; IW(LBEGIW+1) = 43
RW(LBEGRW) = 10.0_r8 ; RW(LBEGRW+1) = 20.0_r8
WRITE(STDOUT,*) "Integer data: ", IW(LBEGIW:LENDIW)
WRITE(STDOUT,*) "Real data: ", RW(LBEGRW:LENDRW)
RW => null()
END
``````

Upon execution with Intel Fortran,

``````C:\temp>ifort /standard-semantics p.for
Intel(R) Fortran Intel(R) 64 Compiler Classic for applications running on Intel(R) 64, Version 2021.1.2 Build 20201208_000000
Copyright (C) 1985-2020 Intel Corporation.  All rights reserved.

Microsoft (R) Incremental Linker Version 14.27.29112.0
Copyright (C) Microsoft Corporation.  All rights reserved.

-out:p.exe
-subsystem:console
p.obj

C:\temp>p.exe
Enter program size parameter value: LI =
(for demo, enter an integer value >= 6)
6
Integer data:  42 43
Real data:  10.0000000000000 20.0000000000000

C:\temp>``````