Passing subsections of arrays in f77

Passing subsections of arrays in f77

Post by rjmagya » Thu, 21 Sep 2006 03:20:55


Hi,

How do I pass a subsection of an array to a subroutine? Thanks, Rudy

For example, suppose I have:

Real*8 a(10,10)

call subrout(a(,10))

. . .

subroutine subrout(a)

real*8 a(10)
 
 
 

Passing subsections of arrays in f77

Post by beliavsk » Thu, 21 Sep 2006 03:44:58


Fortran arrays are stored in column-major order, so you can just write

call subrout(a(1,10))

Here is an example in Fortran 95 using an integer matrix:

subroutine sub(ivec)
integer :: ivec(2)
print*,"in sub, ivec =",ivec
end subroutine sub

program xsub
integer, parameter :: n1 = 2, n2 = 3
integer :: imat(n1,n2),i1,i2
do i1=1,n1
do i2=1,n2
imat(i1,i2) = i1 + 10*i2
end do
end do
print*,"imat=",imat
call sub(imat(1,n2))
end program xsub

output:

imat= 11 12 21 22 31 32
in sub, ivec = 31 32

With few changes this code could be F77. I don't think there is a
simple way to pass a row subsection in F77. Fortran 95 offers more
powerful array handling, including array sections, than F77, and has
many other advantages. Why don't you use it?

 
 
 

Passing subsections of arrays in f77

Post by glen herrm » Thu, 21 Sep 2006 03:45:34


In standard Fortran 77 you pass A(1,10) and realize that the
called routine gets the ten elements in memory starting at A(1,10).

Fortunately for you, those are A(1,10), A(2,10), A(3,10), ...

(The called routine is allowed to read any element from the specified
element through the end. If you pass A(1,5) 60 elements are
available to the called routine. You can also pass the number
of elements, and use that in the DIMENSION statement in the
called routine.)

If you instead wanted to pass A(10,) you would have to make a copy
yourself first.

This is still true for assumed size arrays in versions since 1977.

-- glen
 
 
 

Passing subsections of arrays in f77

Post by nospa » Thu, 21 Sep 2006 03:47:43

m. Ooops. I just finished typing this whole thing and then noticed that
the subject line said f77. Almost all of my reply was about features new
to f90. I guess I'll go ahead and post it anyway, along with a short f77
answer and the suggestion that the OP consider moving to f90/f95.

To answer the original question in a pure f77 context, f77 basically has
no such thing as an array subsection, so you can't readily pass one.
There is a hack you can do for some cases, but it is really quite quirky
and only works at all for very limited cases.

The hack is that you pass an array element as the actual argument. Note
that the OP's a(,10) is not a proper syntax for an array element. THat
looks like an attempt at inventing a syntax for a section, but again,
f77 has no such thing. The syntax for the appropriate array element here
would be a(1,10).

Yes, the syntax is highly confusing and irregular in that there is no
way to tell by looking at it whether you are passing just that array
element or using this hack to pass what amounts to a section.

Anyway, if you pass an array element as an actual argument for an array
dummy argument, the dummy argument "gets" consecutive memory locations
starting with the specified element. You hav eto know about Fortran
array layout in memory in order to figure out what that means in any
particular case. Generally, it ends up meaning that this hack can be
used to pass an array column, but not to pass an array row.

So just changing your a(,10) to a(1,10) will probably do what you want.
But I recommend looking into f90/f95 instead. See below for that.

rjmagyar < XXXX@XXXXX.COM > wrote:


Well, there are 2 basic parts to it.


The first part is how to write a subsection of an array in any context.
This has nothing in particular to do with subroutines. I'm not 100% sure
what you mean by a(,10) above, as that is not a standard Fortran syntax
and you didn't explain what you intended, but I am reasonably confident
that what you meant was the 10th column of the array. That can be
written in any of several ways, the simplest of which is a(:,10). The
":" is the important part here; it indicates that there is a slice along
that dimension. You could more verbosely write a(1:10,10), explicitly
indicating the range of the slice, but the simple colon uses the
declared upper and lower bounds as defaults.


While this will work, it can have several undesired consequences.

1. Of course, if the actual argument ever has a size other than 10,
there will be problems. That might or might not be an isssue for your
code.

2. This might trigger copy-in/copy-out in some cases with some
compilers. Basically, this f77-style (I'm ignoring the non-standard
real*8 bit) array declaration implies that the array is contiguous in
memory (like all f77 arrays are in practice). With array slices, you can
make non-contiguous actual arguments. The compiler is likely to make a
temporary contiguous copy to work around this.

I generally recommend declaring dummy arguments with assumed shape in
f90. That is, use the declaration form (again, I don't recommend the
real*8 part, but I'll copy your usage of that as it is not the main
point here)

real*8 a(:)

This form says to get the array size from whatever is passed as the
actual argument. It also directly accomodates non-contiguous slices; the
information needed to deal with them is also passed from th eactual
argument.

Ho