Hauptprogramm

program matmulf90

!  Multiplikation von Matrizen

  use kinds
  implicit none
      
  integer, parameter ::  n = 1000 
  
  real(kind=REAL8), dimension(n, n) :: a, b, c
  real(kind=REAL8)                  :: spur, spur1, cth
  integer                       :: i, j, k

  
  ! Anfangswerte setzen
  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
  

  ! Matrixmultiplikation
  c = matmul(a, b)

  ! Testdaten ausgeben
  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 program matmulf90