Matrixmultiplikation mit Blas-Routine

      PROGRAM MATMUL
C
C  MULTIPLIKATION VON MATRIZEN
C

      IMPLICIT LOGICAL(A-Z)
      
      INTEGER N
      PARAMETER (N = 1000)

      DOUBLE PRECISION   A(N, N), B(N,N), C(N,N)

      DOUBLE PRECISION ALPHA, BETA
      INTEGER NN
      CHARACTER TA, TB

C     ANFANGSWERTE SETZEN
      CALL INIT(A, B, N)

C     MATRIXMULTIPLIKATION
      CALL DGEMM('N', 'N', N, N, N, 1.0D0, A, N, B, N, 0.0D0, C, N)

C     TESTDATEN AUSGEBEN
      CALL TEST(C, N)

      END


      SUBROUTINE INIT(A, B, N)

      INTEGER            N
      DOUBLE PRECISION   A(N,N), B(N,N)
      INTEGER            I, J

      DO I=1, N
         DO J=1, N
            A(I,J) = I + J
            B(I,J) = I - J
         ENDDO
         B(I,I) = 1.0
      ENDDO

      END

      
      SUBROUTINE TEST(C, N)

      DOUBLE PRECISION   C(N,N)
      INTEGER            N

      DOUBLE PRECISION   SPUR, SPUR1, CTH
      INTEGER            I

      SPUR1 = N * (N+1)
      PRINT *, "SPUR THEORETISCH: ", SPUR1
      SPUR = 0.0
      DO I = 1, N
         SPUR = SPUR + C(I,I)
      ENDDO
      PRINT *, "SPUR            : ", SPUR

      CTH = 333331003.0D0
      PRINT *, "C(1,2) THEOR.   : ", CTH
      PRINT *, "C(1,2)          : ", C(1,2)

      END