program zherk_test
   implicit none

   integer, parameter :: wp = kind(0d0)
   integer, parameter :: funit = 18

   real(wp), parameter :: rone = 1.0
   real(wp), parameter :: rzero = 0.0
   complex(wp), parameter :: cone = 1.0
   complex(wp), parameter :: czero = 0.0

   integer :: m, n, ldest
   complex(wp), allocatable, dimension(:) :: src, dest1, dest2

   ! Read data
   open(unit=funit, file="zherk.dat", status="old")
   read(funit, "(3i5)") m, n, ldest
   allocate(src(n**2))
   read(funit, "(3es25.17)") src(:)
   close(funit)

   ! Try zherk
   allocate(dest1(ldest**2))
   !dest1(:) = 0
   call zherk('U', 'C', m, n, -rone, src, n, rzero, dest1, ldest)
   print *, "dest1(1) = ", dest1(1)

   ! Try zgemm
   allocate(dest2(ldest**2))
   call zgemm('C', 'N', m, m, n, -cone, src, n, src, n, czero, dest2, ldest)
   print *, "dest2(1) = ", dest2(1)

end program zherk_test