passing allocatable array section into subroutine and return error 1

[anex]

Hi, I code fortran in AIX with xlf90 compiler, find a error,
the source code is as follows:

MODULE modaa
CONTAINS
!----------------------------
SUBROUTINE SUBAA
integer,intent(in) :: n
integer,parameter :: m=10
real,dimension(1:m,1:n) :: ary
call subcc(ary,5,n)
print *, ary(3,1:n)
END SUBROUTINE SUBAA
!------------------------------
SUBROUTINE SUBCC(ary,m,n)
integer,intent(in) :: m,n
real,dimension(1:m,1:n),INTENT(OUT) :: ary
integer :: i,j
do i= 1, m
do j = 1, n
ary(i,j) = 5
end do
end do
END SUBROUTINE SUBCC
!---------------------------------
END MODULE modaa
!=======================================
PROGRAM MAIN
USE modaa
call subaa(10)
END PROGRAM MAIN

I run this program, it will lead some error, output print number is error,such as:
5.00000000000000000 5.00000000000000000 5.00000000000000000 5.00000000000000000 5.00000000000000000 0.381534781347033038E-269 0.128822975391891636E-230 0.000000000000000000E+00 0.988131291682493088E-323 -0.310503970101894453E+232
I hope experts can help me? thanks

 

[FJacq]

I don't think that you really know what your are doing with your array :

- in subaa, you create an automatic array, ary(10,10)

- in subcc, you reshape the same array in dimensioning it a(5,10). Notice the verb reshape ! You don't pass really an array section : you pass just the whole array but you decide to change its size and to use only the 50 first elements as if they belong to a matrix(5,10).

This kind of programming is correct as fas as you know what you are doing.

To understand what happened, you have to consider the matrix as a simple vector of 100 values.

- the routine subcc initializes an array of 50 values which simply correspond to the 50 first values of ary

- the routine tries to write the third row of the initial array. But the third row is not composed of contiguous elements ! FORTRAN organizes a matrix column by column. So the elements you are printing are located at the positions 3,13,23,33,43,53,63,73,83,93 is you consider the matrix as as a vector.

Unfortunately, you initialized only the 50 first elemenst of the matrix. This is way the 5 first elements of the third row are OK and not the 5 last ones.

To go a little bit further, I think that you expected to program as follows :

MODULE modaa
  CONTAINS
  !----------------------------
   SUBROUTINE SUBAA(n)
      integer,intent(in)  ::  n
      integer,parameter   ::  m=10
      real,dimension(1:m,1:n)  ::  ary
      call subcc(ary(1:5,1:n))
      print *, ary(3,1:n)
   END SUBROUTINE SUBAA
   !------------------------------
   SUBROUTINE SUBCC(ary)
      real,dimension(:,:),INTENT(OUT)  ::  ary
      integer  ::  i,j
      do i= 1, size(ary,1)
       do j = 1, size(ary,2)
          ary(i,j) = 5
       end do
      end do
   END SUBROUTINE SUBCC
   !---------------------------------
END MODULE modaa
!=======================================
PROGRAM MAIN
    USE modaa
   call subaa(10)
END PROGRAM MAIN

With the result :

[lcoul@localhost test]$ ifort t20.f90
[lcoul@localhost test]$ ./a.out
   5.000000       5.000000       5.000000       5.000000       5.000000    
   5.000000       5.000000       5.000000       5.000000       5.000000

I have modified few things :

- first of all, I really pass an array section to the routine subcc (ary(1:5,1:n))

- after that, I changed the dimensioning in the routine subbc in omitting the true dimensions : dimension,. This is the way to inherit the dimensions from the calling routine;

- there is no need to pass the dimensions themselves because the intrinsic function SIZE allows you to retrieve them. But passing them as argument is correct as well.

Again, I could go deeper in my explanations. But as exercise, try yourself to understand the differences, in the routine subcc, of the following declarations :

   SUBROUTINE SUBCC(ary,n,m)
      INTEGER,INTENT(i) :: m,n
      real,dimension(:,:),INTENT(OUT)  ::  ary
      real,dimension(1:,1:),INTENT(OUT)  ::  ary
      real,dimension(m,*),INTENT(OUT)  ::  ary
      real,dimension(m,n),INTENT(OUT)  ::  ary
      real,dimension(1:m,1:n),INTENT(OUT)  ::  ary
   END SUBROUTINE

 

[anex]

thanks, I seem to understand your explanation. If I want to get my value, I must specify the dimension of array ary when calling subroutine subcc, such as CALL subcc(ary(1:5,1:n),5,n) or CALL subcc(ary,1:6),10,6). If I pass array ary into subroutine subcc with 'CALL subcc(ary,5,n) or subcc(ary,10,6)', it just pass 5*n or 10*6 first element of array ary into subcc, so after return value from subcc, it just change the value of matching address. is that so?

 

[FJacq]

Partly right, but not totally. You cannot consider only the "call" statement. The definition of the arguments in the called subroutine plays an important role too.

Example 1 :

CALL subcc(ary(1:5,1:n),5,n)

SUBROUTINE subcc(ary,m,n)
  INTEGER,INTENT(in) :: m,n
  real,INTENT(inout) :: ary(:,:)
END SUBROUTINE

This is the solution I proposed. It is even not necessary to pass the dimensions because you can get them with the SIZE function (5 == SIZE(ary,1) and n==SIZE(ary,2)).

Example 2 :

CALL subcc(ary(1:5,1:n),5,n)

SUBROUTINE subcc(ary,m,n)
  INTEGER,INTENT(in) :: m,n
  real,INTENT(inout) :: ary(m,n)
END SUBROUTINE

Again you pass the array section. But the definition ary(m,n) in the routine means that "subcc" wants an array with contiguous elements. Unfortunately, the array section is not composed of contiguous elements. Despite that, the compiler will do the correct job in allocating a temporary array with contiguous elements, copying inside the section, calling the subroutine passing to it the temporary array, and finally, copying back the temporary array into the array section. One calls that the "copy in / copy out" technique. As you see, this solution 2 works but can be much more expensive than the first one.

Example 3 :

CALL subcc(ary,5,n)

SUBROUTINE subcc(ary,m,n)
  INTEGER,INTENT(in) :: m,n
  real,INTENT(inout) :: ary(n,m)
END SUBROUTINE

This is your initial programming : you pass the whole array and you reshape it; So you use simply the 5*n first elements... and you get a strange result when printing the third row of the matrix which original shape is (10,10).

Example 4 :

CALL subcc(ary,5,n)

SUBROUTINE subcc(ary,m,n)
  INTEGER,INTENT(in) :: m,n
  real,INTENT(inout) :: ary(:,:)
  integer :: i,j
  DO i=1,n
    DO j=1,m
      a(i,j)=5
    ENDDO
  ENDDO
END SUBROUTINE

This is another correct possibility : you pass the whole array and you inherit automatically of the right dimensions. After that, you use only the 5 first rows. In that case SIZE(ary,1) returns 10. This will give you also the right answer and without the need of the expensive "copy in/copy out".

But this last solution does not work if you want to pass the section (3:8,2:4) for instance... except if you decide to pass additional information to describe more precisely the array section. So only the solutions 1 and 2 are OK with 3 arguments, with a preference for the solution 1 avoiding "copy in/copy out" and limiting the number of needed arguments (only 1 is enough !).

 

[anex]

thank you very much! you are most kind!

I think your explanation is very good for me, thanks again!

if I use example 2, such as CALL subcc(ary(1:5,1:n),5,n),
when running this program , is it need additional memory to store temporary array? that is, example 2 comparing with example 1, it need more memory to run. is that so?

 

[FJacq]

if I use example 2, such as CALL subcc(ary(1:5,1:n),5,n),
when running this program , is it need additional memory to store temporary array? that is, example 2 comparing with example 1, it need more memory to run. is that so?

Yes. More memory and more time consuming!

But I repeat again that each example I gave was the association of a specific calling sequence with a particular subroutine argument definition. For instance, the examples 1 and 2 share the same "call" statement. Considering this one alone is therefore not sufficient to understand what could happen.