[Patch, gfortran]: Enable loop unrolling for the matmul intrinsic.
Janne Blomqvist
jblomqvi@cc.hut.fi
Tue Nov 29 20:02:00 GMT 2005
Hi,
not so long ago my proposal to enable loop unrolling for the matmul
intrinsic was turned down due to fears about a negative performance
impact on register starved architectures.
However, I now believe there is sufficient benchmark data to dispel
these fears. I made a slightly improved version of my benchmark
program to also test double precision and logicals, here are the
result (well for logicals the correct unit should be gops/s not
gflops/s, and probably the ops count is wrong anyway, but it should
provide some clue about the relative speeds with and without loop
unrolling):
trunk on 1.8 GHz A64, i686-pc-linux-gnu:
Single precision matrix multiplication test
Matrix side size Matmul (Gflops/s) sgemm (Gflops/s) Loops
====================================================================
2 0.092 0.024 100000
4 0.361 0.162 100000
8 0.733 0.495 100000
16 0.687 0.768 100000
32 0.874 1.316 15500
64 0.988 1.397 1922
128 1.041 2.929 239
256 0.848 4.809 29
512 0.794 5.093 3
1024 0.800 5.173 1
2048 0.807 5.368 1
Double precision matrix multiplication test
Matrix side size Matmul (Gflops/s) dgemm (Gflops/s) Loops
====================================================================
2 0.086 0.020 100000
4 0.373 0.158 100000
8 0.762 0.522 100000
16 0.703 0.949 100000
32 0.889 1.306 15500
64 0.973 1.391 1922
128 0.777 1.643 239
256 0.698 2.259 29
512 0.521 2.673 3
1024 0.522 2.749 1
2048 0.527 2.759 1
Default kind logical matrix multiplication test
Matrix side size Matmul (Gflops/s) Loops
==============================================
2 0.109 100000
4 0.311 100000
8 0.466 100000
16 0.839 100000
32 1.493 15500
64 2.659 1922
128 5.642 239
256 10.914 29
512 22.990 3
1024 43.812 1
2048 75.675 1
matmul and matmull compiled with -funroll-loops:
Single precision matrix multiplication test
Matrix side size Matmul (Gflops/s) sgemm (Gflops/s) Loops
====================================================================
2 0.057 0.012 100000
4 0.215 0.091 100000
8 0.447 0.267 100000
16 0.683 0.750 100000
32 1.242 1.342 15500
64 1.316 1.458 1922
128 1.227 2.981 239
256 1.047 4.761 29
512 0.951 5.093 3
1024 0.947 5.173 1
2048 0.958 5.265 1
Double precision matrix multiplication test
Matrix side size Matmul (Gflops/s) dgemm (Gflops/s) Loops
====================================================================
2 0.100 0.021 100000
4 0.400 0.158 100000
8 0.781 0.513 100000
16 1.038 0.942 100000
32 1.224 1.326 15500
64 1.282 1.408 1922
128 0.957 1.653 239
256 0.738 2.285 29
512 0.547 2.664 3
1024 0.545 2.731 1
2048 0.557 2.784 1
Default kind logical matrix multiplication test
Matrix side size Matmul (Gflops/s) Loops
==============================================
2 0.120 100000
4 0.361 100000
8 0.619 100000
16 0.990 100000
32 1.852 15500
64 3.135 1922
128 6.614 239
256 12.781 29
512 26.822 3
1024 48.789 1
2048 80.649 1
As can be seen, with -funroll-loops performance is about 20-30%
better.
There are also some results posted for ppc, where loop unrolling
improved performance by 15-20%. See
http://gcc.gnu.org/ml/fortran/2005-11/msg00644.html
And once again, thanks to rth for showing the correct make
incarnation.
--
Janne Blomqvist
-------------- next part --------------
libgfortran ChangeLog:
2005-11-29 Janne Blomqvist <jb@gcc.gnu.org>
* Makefile.am: Enable loop unrolling for matmul.
* configure: Regenerated.
* Makefile.in: Regenerated.
* aclocal.m4: Regenerated.
-------------- next part --------------
Index: Makefile.am
===================================================================
--- Makefile.am (revision 107677)
+++ Makefile.am (working copy)
@@ -575,8 +575,9 @@ $(gfor_built_specific2_src) \
intrinsics/dprod_r8.f90 \
intrinsics/f2c_specifics.F90
-# Turn on vectorization for matmul.
-$(patsubst %.c,%.lo,$(notdir $(i_matmul_c))): AM_CFLAGS += -ftree-vectorize
+# Turn on vectorization and loop unrolling for matmul.
+$(patsubst %.c,%.lo,$(notdir $(i_matmul_c))): AM_CFLAGS += -ftree-vectorize -funroll-loops
+$(patsubst %.c,%.lo,$(notdir $(i_matmull_c))): AM_CFLAGS += -funroll-loops
BUILT_SOURCES=$(gfor_built_src) $(gfor_built_specific_src) \
$(gfor_built_specific2_src)
-------------- next part --------------
program matmul_bench
implicit none
integer, parameter :: sp = selected_real_kind(4), &
dp = selected_real_kind(15)
call runsbench (2500)
call rundbench (2500)
call runlbench (2500)
contains
! Run single precision matrix mult benchmark with different sized arrays.
subroutine runsbench (nmax)
integer, intent(in) :: nmax
real(sp), allocatable, dimension(:,:) :: a, b, res
integer :: n, loop
real(dp) :: time, flops, time2
print *, 'Single precision matrix multiplication test'
print *, 'Matrix side size Matmul (Gflops/s) sgemm (Gflops/s) Loops'
print *, '===================================================================='
n = 2
do
allocate (a(n,n), b(n,n), res(n,n))
call random_number (a)
call random_number (b)
res = 0.0_sp
! matmul for square matrix is (2n-1)*n**2 flops.
flops = (2.0_dp * real (n, dp) - 1.0_dp) * real (n, dp)**2
! Assuming an on average 1 gflop/s cpu, 1e9 flops takes about 1 second and
! should be enough. We also do a maximum of 1e5 loops, since
! for small arrays the overhead is large.
loop = max (min (int (1.0e9_dp / flops), 10**5), 1)
call smatmul_timing (a, b, res, loop, time)
res = 0.0_sp
call sgemm_timing (a, b, res, loop, time2)
print '(I5,15X,F6.3,15X,F6.3,15X,I7)', n, &
flops * real(loop, dp) / time / 1.0e9_dp, &
flops * real (loop, dp) / time2 / 1.0e9_dp, loop
deallocate (a, b, res)
n = n * 2
if (n > nmax) exit
end do
end subroutine runsbench
! Run double precision matrix mult benchmark with different sized arrays.
subroutine rundbench (nmax)
integer, intent(in) :: nmax
real(dp), allocatable, dimension(:,:) :: a, b, res
integer :: n, loop
real(dp) :: time, flops, time2
print *, 'Double precision matrix multiplication test'
print *, 'Matrix side size Matmul (Gflops/s) dgemm (Gflops/s) Loops'
print *, '===================================================================='
n = 2
do
allocate (a(n,n), b(n,n), res(n,n))
call random_number (a)
call random_number (b)
res = 0.0_dp
! matmul for square matrix is (2n-1)*n**2 flops.
flops = (2.0_dp * real (n, dp) - 1.0_dp) * real (n, dp)**2
! Assuming an on average 1 gflop/s cpu, 1e9 flops takes about 1 second and
! should be enough. We also do a maximum of 1e5 loops, since
! for small arrays the overhead is large.
loop = max (min (int (1.0e9_dp / flops), 10**5), 1)
call dmatmul_timing (a, b, res, loop, time)
res = 0.0_dp
call dgemm_timing (a, b, res, loop, time2)
print '(I5,15X,F6.3,15X,F6.3,15X,I7)', n, &
flops * real(loop, dp) / time / 1.0e9_dp, &
flops * real (loop, dp) / time2 / 1.0e9_dp, loop
deallocate (a, b, res)
n = n * 2
if (n > nmax) exit
end do
end subroutine rundbench
! Run logical matrix mult benchmark with different sized arrays.
subroutine runlbench (nmax)
integer, intent(in) :: nmax
logical, allocatable, dimension(:,:) :: a, b, res
real(dp), allocatable, dimension(:,:) :: rtmp
integer :: n, loop
real(dp) :: time, flops
print *, 'Default kind logical matrix multiplication test'
print *, 'Matrix side size Matmul (Gflops/s) Loops'
print *, '=============================================='
n = 2
do
allocate (a(n,n), b(n,n), res(n,n), rtmp(n,n))
call random_number (rtmp)
a = .false.
where (rtmp > 0.5)
a = .true.
end where
call random_number (rtmp)
b = .false.
where (rtmp > 0.5)
b = .true.
end where
res = .false.
! matmul for square matrix is (2n-1)*n**2 flops.
flops = (2.0_dp * real (n, dp) - 1.0_dp) * real (n, dp)**2
! Assuming an on average 1 gflop/s cpu, 1e9 flops takes about 1 second and
! should be enough. We also do a maximum of 1e5 loops, since
! for small arrays the overhead is large.
loop = max (min (int (1.0e9_dp / flops), 10**5), 1)
call lmatmul_timing (a, b, res, loop, time)
print '(I5,15X,F6.3,15X,I7)', n, &
flops * real(loop, dp) / time / 1.0e9_dp, &
loop
deallocate (a, b, res)
n = n * 2
if (n > nmax) exit
end do
end subroutine runlbench
! Actual routine, and timing.
subroutine smatmul_timing (a, b, res, loop, time)
real(sp), intent(in), dimension(:, :) :: a, b
real(sp), intent(inout) :: res(:,:)
integer, intent(in) :: loop
real(dp), intent(out) :: time
real(dp) :: t1, t2
integer :: i
call cpu_time (t1)
do i = 1, loop
res = matmul (a, b)
end do
call cpu_time (t2)
time = t2 - t1
end subroutine smatmul_timing
! Actual routine, and timing.
subroutine dmatmul_timing (a, b, res, loop, time)
real(dp), intent(in), dimension(:, :) :: a, b
real(dp), intent(inout) :: res(:,:)
integer, intent(in) :: loop
real(dp), intent(out) :: time
real(dp) :: t1, t2
integer :: i
call cpu_time (t1)
do i = 1, loop
res = matmul (a, b)
end do
call cpu_time (t2)
time = t2 - t1
end subroutine dmatmul_timing
! Actual routine, and timing.
subroutine lmatmul_timing (a, b, res, loop, time)
logical, intent(in), dimension(:, :) :: a, b
logical, intent(inout) :: res(:,:)
integer, intent(in) :: loop
real(dp), intent(out) :: time
real(dp) :: t1, t2
integer :: i
call cpu_time (t1)
do i = 1, loop
res = matmul (a, b)
end do
call cpu_time (t2)
time = t2 - t1
end subroutine lmatmul_timing
subroutine sgemm_timing (a, b, res, loop, time)
real(sp), intent(in), dimension(:, :) :: a, b
real(sp), intent(inout) :: res(:,:)
integer, intent(in) :: loop
real(dp), intent(out) :: time
real(dp) :: t1, t2
integer :: i, n
n = size (a, 1)
call cpu_time (t1)
do i = 1, loop
call sgemm('n','n',n, n, n, 1.0_sp, a, n, b, n, 0.0_sp, res, n)
end do
call cpu_time (t2)
time = t2 - t1
end subroutine sgemm_timing
subroutine dgemm_timing (a, b, res, loop, time)
real(dp), intent(in), dimension(:, :) :: a, b
real(dp), intent(inout) :: res(:,:)
integer, intent(in) :: loop
real(dp), intent(out) :: time
real(dp) :: t1, t2
integer :: i, n
n = size (a, 1)
call cpu_time (t1)
do i = 1, loop
call dgemm('n','n',n, n, n, 1.0_dp, a, n, b, n, 0.0_dp, res, n)
end do
call cpu_time (t2)
time = t2 - t1
end subroutine dgemm_timing
end program matmul_bench
-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: application/pgp-signature
Size: 185 bytes
Desc: not available
URL: <http://gcc.gnu.org/pipermail/fortran/attachments/20051129/8a443c5e/attachment.sig>
More information about the Fortran
mailing list