Re: Is it a bug: data aligment of matrix returned by gemm

[Tim Prince <TimothyPrince at sbcglobal dot net>]

Maciej Zwierzycki wrote:
> On Thu, Jun 25, 2009 at 06:04:24AM -0700, Tim Prince wrote:
> [...] 
>> The compiler can't diagnose your syntax error, where you declare a with
>> intent(in), although I suppose it could in principle if it could scan the
>> source of dgemm.  It should diagnose a syntax error when you don't set
>> rank, and I can't see how you could get results without rank.  
> 
> Why not? Assumed shape arrays are part of the standard as far as I know.
> also what is wrong with specifying intent(in) for a matrices not to be
> changed by the function?
> 
You declared the result array intent(in).  The compiler would be justified
in discarding your results, as you forbade it to return the results to the
caller.
>> I don't see
>> your complaint about how dgemm fills a; it seems, contrary to what you
> 
> a is _not_ being filled. Is it a typo or we are talking about
> different things?
> 
>> show, that you may have declared a such that size(a,1) is 2, and the
>> result then is what you requested.
>> I don't see how any of the stuff you mention could be deemed "correct" if
>> you don't provide correct source code.
> 
> Can you be more specific? What exactly is wrong with the code?
> The function itself, or the way it is called?
I don't see how we can guess what you actually did from posted pseudo-fortran.
> 
> The reason I do not specify ranks in the definition of the function
> this is meant to be a part of generic function - a drop-in replacement
> for intrinsic matmul. Which works as intended with two other
> compilers - three actualy if I remember the bahaviour of
> Lahey-Fujitsu from few years back correctly.
gfortran does what it can to substitute ?gemm when implementing matmul, if
you set -fexternal-blas.  That does generate a BLAS call which uses the
arrays consistent with the way they are declared by the caller.  It
doesn't take advantage in general of the way the BLAS call is designed to
avoid allocation of an additional result matrix.
You appear to be asking both the compiler and the readers on this forum to
guess what you want, and apparently what you want is different from what
you wrote.
> 
> If the function returns 2x2 matrix why should I not be allowed to 
> use it to fill a 2x2 sub-block of a multi-dimensional matrix?
> BTW. the problem is also present for 2x2 blocks of larger 2D matrix.
> 
> I understand what is happening. The function returns a
> vector of 4 values. These are then being used to fill the bigger matrix
> staring from the pointer defined in the 
> 
> c(1,:,:)=matmul_func(....)
> 
> The filling however is done column-wise for the whole c matrix, not the
> c(1,:,:) sub-matrix. BTW. specifying something like c(1,1:2,1:2) does
> not help.
> 
You can't get a returned result in an array section.  It may run OK, but
then your syntax implicitly discards the result.  But you already told the
compiler it could discard the result by intent(in).