(j3.2006) If SIZE(B) > SIZE(A), what part of B does RESHAPE(B, SHAPE(A)) contain, if both A, B are rank N arrays ?

Toon Moene toon
Fri Dec 9 13:59:23 EST 2011


It might sound academic, but this question results from a real 
(standards' violation) problem I solved yesterday.

Somewhere in the 3.5 million lines of code of our weather forecasting 
code, someone had written:

A(:,:,:,:) = B(:,:,:,:)

Looked perfectly OK, has worked for years on various platforms, until I 
tried it at home using gfortran 4.6.2.

BWAM - segmentation fault.

Due to the way the code was set up, I knew both arrays were 
allocatables, so my first thought was: one of them is not allocated. 
The second thought was: B is larger than A.  So I added some debug 
prints that showed that both were allocated, and B was larger than A 
(it's shape differed in one position).

So why did this "work" with the other compilers ?  Fortunately, I had 
seen *that* one before; the answer is very simple: The other compilers 
use the shape of the LHS of the assignment to construct their 
scalarization loops, and gfortran uses the shape of the RHS; after all, 
if you have to scalarize expressions, in case of

PRINT*, A+B+C

why wouldn't you re-use that code for assignments ?

So now my task was to change this non-standard code (that was working 
correctly under the assumption that the shape of the LHS determined the 
copying) into something that would work independent of this assumption.

I came up with:

   ISHAPE = SHAPE(A)  ! ISHAPE is a local 4 element integer array
   A(:,:,:,:)=B(LBOUND(B,1):ISHAPE(1)-LBOUND(B,1)+1,&
   &            LBOUND(B,2):ISHAPE(2)-LBOUND(B,2)+1,&
   &            LBOUND(B,3):ISHAPE(3)-LBOUND(B,3)+1,&
   &            LBOUND(B,4):ISHAPE(4)-LBOUND(B,4)+1)

However, today I thought: I'm just doing RESHAPE's job here - by hand.

Certainly, one must be able to replace the above with the following:

   A = RESHAPE(B, SHAPE(A))

Unfortunately, perusing paragraph 13.7.140 of the 2008 Standard 
(actually the 09-007r2 pdf), I can't convince myself that this is the 
case.  It says remarkably little about RESHAPEing a larger array to a 
smaller: "The elements of the result ... are those of SOURCE in normal 
array element order ...".

This doesn't sound like the rank N analogon of the "upper-left hand 
corner" of a rank 2 array, which the above hand-coded assignment 
certainly is.

Can RESHAPE replace the above convolution ?

Thanks in advance,

-- 
Toon Moene - e-mail: toon at moene.org - phone: +31 346 214290
Saturnushof 14, 3738 XG  Maartensdijk, The Netherlands
At home: http://moene.org/~toon/; weather: http://moene.org/~hirlam/
Progress of GNU Fortran: http://gcc.gnu.org/wiki/GFortran#news



More information about the J3 mailing list