------- Additional Comments From paulthomas2 at wanadoo dot fr 2005-04-29
09:51 -------
Tobi,
The reason that it works is that the value of *data that is passed to matmul
points to the first element of the array or subarray. matmul then calculates
the number of elements from {stride, lbound, ubound} and C-indexing is used for
athe aritmetic. *base is therefore never needed. This is illustrated by the
following example:
Program pr18857a
integer, parameter :: N = 5
integer, parameter :: T = 4
real(kind=T), dimension(:,:), allocatable, Target :: a, b, c
real(kind=T), dimension(:,:), POINTER :: d, e
real(kind=T), dimension(N,N) :: x, y, z
allocate (a(2*N, N), b(N, N), c(2*N, N))
a = 1.0_T
a(1:N,:) = 2.0_T
b = 4.0_T
x = 1.0_T
y = 2.0_T
z = 0.0_T
z = matmul (x, y)
if (sum (z) /= 250.0_T) call abort ()
c = 0.0_T
c = matmul (a(:,3:N), b(3:N,:))
if (sum (c) /= 900.0_T) call abort ()
deallocate (a, b, c)
end program pr18857a
I will not pollute Bugzilla with the entire tree dump but limit myself to
reproducing the part just before and including the second call to matmul (lines
384-386 of pr18857a.f90.t02.original)
parm.21.data = (real4[0:] *) (real4[0:] *) &(*b.data)[(3 - b.dim[0].lbound)
* D.595 + (D.588 - b.dim[1].lbound) * D.596];
parm.21.offset = 0;
_gfortran_matmul_r4 (&c, &parm.20, &parm.21);
The (3 - b.dim[0].lbound) provides the offset to the data being multiplied.
Earlier on in the tree dump, D.588 = b.dim[1].lbound;, so the second part of
the offset is zero. By some miracle, the front end does this correctly for
every case that I can find.
--
http://gcc.gnu.org/bugzilla/show_bug.cgi?id=18857
