!----------------------------------------------------------------------
!     Copyright (c) 2025 Shin-ichi Takehiro. All rights reserved.
!----------------------------------------------------------------------
!
! 表題  spml w_mpi_module 速度計測プログラム
!
! 履歴  2025/10/04 竹広真一  ispack3/sample/sypack-test/test.f90 を改造
!       
!       デフォルト spml ver.1.0.0 用
!       -DSPML1 で spml ver.0.8 用
!
program w_mpi_module_bench
  use w_mpi_module
  use mpi
  
  implicit none
  
  INTEGER(8),PARAMETER :: JM=2**10,NTR=10
!  INTEGER(8),PARAMETER :: JM=2**8,NTR=10
!  INTEGER(8),PARAMETER :: JM=2**6,NTR=10
  INTEGER(8),PARAMETER :: MM=JM-1,IM=JM*2
  INTEGER(8),PARAMETER :: NM=MM,NN=NM

  INTEGER(8) :: N,M,L,ICPU,ITR,MAXTD
  INTEGER :: IERR,NP,IP
  REAL(8) :: RAN,GFLOPS,RC,TIM0,TIM1  

  REAL(8),DIMENSION(:,:),ALLOCATABLE :: xv_Data
  REAL(8),DIMENSION(:),ALLOCATABLE   :: w_Data

  integer :: nwmax
  !$    INTEGER :: omp_get_max_threads

  RC=1D0*5*IM*LOG(1D0*IM)/LOG(2D0)*0.5D0*JM+1D0*(MM+1)*(MM+1)*JM
  ! 逆変換/正変換の1回あたりの演算数概算

  CALL MPI_INIT(IERR)
  CALL MPI_COMM_SIZE(MPI_COMM_WORLD,NP,IERR)
  CALL MPI_COMM_RANK(MPI_COMM_WORLD,IP,IERR)

  call w_mpi_initial(Int(MM),Int(IM),Int(JM))
  allocate(xv_Data(0:im-1,jc))

#ifdef SPML1
  nwmax = (nm+1)*(nm+1)
#else
  nwmax = nc
#endif
  allocate(w_Data(nwmax))

  IF(IP.EQ.0) THEN
     PRINT '(A,I5,A,I5,A,I5,A,I5,A,I4)','MM=',MM,', IM=',IM,' JM=',JM,', NTR=',NTR
     CALL MXGCPU(ICPU)
     IF(ICPU.EQ.0) THEN
        PRINT '(A)','SSE=fort'
     ELSE IF(ICPU.EQ.10) THEN
        PRINT '(A)','SSE=avx'
     ELSE IF(ICPU.EQ.20) THEN
        PRINT '(A)','SSE=fma'
     ELSE IF(ICPU.EQ.30) THEN
        PRINT '(A)','SSE=avx512'
     ELSE IF(ICPU.EQ.100) THEN
        PRINT '(A)','SSE=sx'
     ELSE IF(ICPU.EQ.1000) THEN
        PRINT '(A)','SSE=fx'
     END IF
     MAXTD=1
     !$       MAXTD=omp_get_max_threads()      
     PRINT '(A,I3)','number of threads =',MAXTD
     PRINT '(A,I3)','number of processes =',NP
  END IF

  DO L=1,nwmax
     call random_number(RAN)            
     w_data(L)=2*RAN-1
  END DO

  CALL MXTIME(TIM0)
  DO ITR=1,NTR
     xv_Data = xv_w(w_data)
  END DO
  CALL MXTIME(TIM1)
  GFLOPS=RC*NTR/(TIM1-TIM0)/1D9
  IF(IP.EQ.0) THEN
     PRINT '(A,ES9.2,A,F6.1,A)','S2G: ',(TIM1-TIM0)/NTR,' sec  (',  &
          &      GFLOPS,' GFlops)'
  END IF

  CALL MXTIME(TIM0)
  DO ITR=1,NTR
     w_data = w_xv(xv_Data)
  END DO
  CALL MXTIME(TIM1)
  GFLOPS=RC*NTR/(TIM1-TIM0)/1D9
  IF(IP.EQ.0) THEN
     PRINT '(A,ES9.2,A,F6.1,A)','G2S: ',(TIM1-TIM0)/NTR,' sec  (',  &
          &      GFLOPS,' GFlops)'
  END IF

  CALL MPI_FINALIZE(IERR)

END program
