      PROGRAM MPI_TEST

      include 'mpif.h'
      parameter (inpes=2, jnpes=2, npes_mod=inpes*jnpes)
!      parameter (iquilt_group=1)
      parameter (tag1=111, tag2=112, tag3 = 113, tag4=114, tag5 = 115)
      parameter (tagu=123, tagv=1234, tagw=12345, im=247, jm=231)
      integer, PARAMETER:: maxpes=6
      integer :: ierr, status(MPI_STATUS_SIZE)
      integer :: etaserverandwam,ieven,jeven
      integer :: npes, mype,yourrank,source,comdup,mpi_comm_comp
      integer :: npesinter, mypeinter, groupinter
      integer :: parentsize, parentgroup
      integer :: j,iworld_minus,iwind,iwave
      integer :: iworldgroup,newcom,mpi_comm_test,comm2
      integer :: info, ihrst,ihr
      integer :: root=0
      integer :: INTERCOMM
      integer :: errCodes(1:maxpes), inumq(100)
      integer :: mpi_comm_inter_array(100), idat(3)
      real :: etantsd, waventsd
      REAL, ALLOCATABLE, DIMENSION(:,:) :: U10, V10, wavepar_global
      integer, allocatable :: irank ( : )
      logical yes

! -------------------------------------
!  Initialize MPI
! -------------------------------------
      CALL MPI_INIT(ierr)
      CALL MPI_COMM_SIZE(MPI_COMM_WORLD,npes,ierr)
      CALL MPI_COMM_RANK(MPI_COMM_WORLD,mype,ierr)
      write(6,*)'MPI_TEST npes,mype ',npes, mype
      iwave=500
C
      ieven=-1
      jeven=-1
      etaserverandwam=-1
      if (mod(inpes,2).eq.0) ieven=inpes  !reserve even threads for I/O servers
      if (mod(jnpes,2).eq.0) jeven=jnpes
      etaserverandwam = max(ieven,jeven)
      if (etaserverandwam.eq.-1) then
         print*,'SETUP_SERVERS ERROR ',etaserverandwam
         print*,'RECOMPILE WITH EVEN INPES OR JNPES'
         CALL MPI_ABORT(MPI_COMM_WORLD,1,IERR)
      end if
      iqserver = etaserverandwam
      iquilt_group = 1
cpk      iqserver = NPES - NPES_MOD - 1
cpk      iqserver = NPES - NPES_MOD
      print*,'MPI_TEST iqserver ',iqserver
      if ( iqserver .eq. 0 ) then
         if ( mype .eq. 0 ) then
           print *, ' *** you specified 0 I/O servers '
           print *, ' CHKOUT will write a file'
         end if
         iquilt_group = 0
      end if
      if ( iquilt_group .gt. iqserver )  then
          iquilt_group = iqserver
          print *, ' ***** NOT ENOUGH SERVERS'
          print *, ' ***** WE NEED TO REDUCE THE NUMB OF SERVER GROUPS'
          print *, ' ***** NUMB OF SERVER GROUPS IS ', iquilt_group
      end if
      do i = 0, iquilt_group - 1
!         call para_range(1,iqserver,iquilt_group,i,istaq,iendq)
         iwork1 = ( iqserver - 1 + 1 ) / iquilt_group
         iwork2 = mod ( iqserver - 1 + 1, iquilt_group )
         istaq = i * iwork1 + 1 + min ( i, iwork2 )
         iendq = istaq + iwork1 - 1
         if ( iwork2 .gt. i ) iendq = iendq + 1
         inumq(i+1) = iendq-istaq+1
         if ( mype .eq. 0 ) print *, ' i, inumq = ',i+1,inumq(i+1)
      end do
C
      if ( mype .lt. NPES_MOD ) then
         icolor = 0
      else
         istaxx = NPES_MOD
         do i = 1, iquilt_group
            iendxx = istaxx + inumq(i) - 1
            if ( mype .ge. istaxx .and. mype .le. iendxx ) then
               icolor = i
            end if
            istaxx = iendxx + 1
         end do
      end if
C
      call mpi_comm_dup(MPI_COMM_WORLD,comdup,ierr)
      call mpi_comm_split(comdup,icolor,mype,mpi_comm_comp,ierr)
      print*,'MPI_TEST comdup,icolor,mype ',comdup,icolor,mype
C
      allocate ( irank ( iqserver ) )
      ixx = NPES_MOD
      do i = 1, iquilt_group
         yes = .true.
         if ( mype .lt. NPES_MOD ) then
            irlr = ixx
         else
            irlr = 0
         end if
      icc = 0
      iss = NPES_MOD
C     THIS IS THE FIRST POSSIBLE TASK ID THAT COULD BE EXCLUDED
      do jj = 1, iquilt_group
         if ( jj .ne. i ) then
            issl = iss
            do kk = 1, inumq(jj)
               icc = icc + 1
               irank(icc)= issl
               if ( mype .eq. issl ) yes = .false.
               issl = issl + 1
            end do
         end if
         iss = iss + inumq(jj)
      end do
C
      iworld = MPI_COMM_WORLD
      call mpi_comm_group(iworld,igroup,ierr)
      print*,'MPI_TEST igroup ',igroup
      call mpi_group_excl(igroup,icc,irank,igroup_x,ierr)
      print*,'MPI_TEST icc,irank,igroup_x ',icc,irank,igroup_x
      call mpi_comm_create(iworld,igroup_x,iworld_minus,ierr)
      print*,'MPI_TEST iworld_minus ',iworld_minus
      call mpi_group_free(igroup,ierr)
      call mpi_group_free(igroup_x,ierr)
C
!      if ( yes ) then
!      print*,'MPI_TEST i,irlr ',i,irlr
!      call mpi_intercomm_create(mpi_comm_comp,0,iworld_minus,irlr,0,
!     *   mpi_comm_inter_array(i),ierr)
!      mpi_comm_inter = mpi_comm_inter_array(i)
!      print*,'MPI_TEST i,irlr,mpi_comm_inter',i,irlr,mpi_comm_inter
!      end if
C
      call mpi_barrier(MPI_COMM_WORLD,ierr)
C
      end do
!       print*,"MPI_TEST iwave ",iwave
!       call MPI_BCAST(iwave,1,MPI_INTEGER,0,MPI_COMM_WORLD,ierr)
C
      ALLOCATE(U10(im,jm))
      ALLOCATE(V10(im,jm))
      ALLOCATE(wavepar_global(im,jm))
C  RECEIVE FROM QUILT
      do ntsd=1,2400,20
         print*,"MPI_TEST recv ntsd ",ntsd
         call mpi_recv
     &     (idat,3,MPI_INTEGER,4,tag1,MPI_COMM_WORLD,status,ierr)
         call mpi_recv
     &     (ihrst,1,MPI_INTEGER,4,tag2,MPI_COMM_WORLD,status,ierr)
         call mpi_recv
     &     (ihr,1,MPI_INTEGER,4,tag3,MPI_COMM_WORLD,status,ierr)
         call mpi_recv
     &     (etantsd,1,MPI_REAL,4,tag4,MPI_COMM_WORLD,status,ierr)
         call mpi_recv
     &     (waventsd,1,MPI_REAL,4,tag5,MPI_COMM_WORLD,status,ierr)
         call mpi_recv
     &     (U10,im*jm,MPI_REAL,4,tagu,MPI_COMM_WORLD,status,ierr)
         call mpi_recv
     &     (V10,im*jm,MPI_REAL,4,tagv,MPI_COMM_WORLD,status,ierr)
         print*,"MPI_TEST recv idat= ",idat
         print*,"MPI_TEST recv ihrst= ",ihrst
         print*,"MPI_TEST recv ihr= ",ihr
         print*,"MPI_TEST recv etantsd= ",etantsd
         print*,"MPI_TEST recv waventsd= ",waventsd
         print*,"MPI_TEST recv U10 ",U10(im/2,jm/2), ntsd
         print*,"MPI_TEST recv V10 ",V10(im/2,jm/2), ntsd
!
         aminu=100000.
         amaxu=-100000.
         aminv=100000.
         amaxv=-100000.
         do i=1,im
         do j=1,jm
            if (U10(i,j).lt.aminu) aminu=U10(i,j)
            if (U10(i,j).gt.amaxu) amaxu=U10(i,j)
            if (V10(i,j).lt.aminv) aminv=V10(i,j)
            if (V10(i,j).gt.amaxv) amaxv=V10(i,j)
         enddo
         enddo
         print*,"MPI_TEST recv min-max U10 ",aminu,amaxu, ntsd
         print*,"MPI_TEST recv min-max V10 ",aminv,amaxv, ntsd
!         call mpi_barrier(MPI_COMM_WORLD,ierr)
C  SEND TO TURBL
         aminw=100000.
         amaxw=-100000.
         do i=1,im
         do j=1,jm
            wavepar_global(i,j)=ntsd*(i+j)/2.
            if (wavepar_global(i,j).lt.aminw) aminw=wavepar_global(i,j)
            if (wavepar_global(i,j).gt.amaxw) amaxw=wavepar_global(i,j)
         enddo
         enddo
         print*,"MPI_TEST send ntsd ",ntsd
         print*,"MPI_TEST send min-max wavepar ",aminw,amaxw, ntsd
         call mpi_send
     &     (wavepar_global,im*jm,MPI_REAL,2,tagw,MPI_COMM_WORLD,ierr)
         print*,"MPI_TEST send wavepar_global",wavepar_global(im/2,jm/2)
     &           ,ntsd
!         call mpi_barrier(MPI_COMM_WORLD,ierr)
      enddo
C  SEND TO TURBL
!      do ntsd=2,1000,4
!         amin=100000.
!         amax=-100000.
!         do i=1,im
!         do j=1,jm
!            wavepar_global(i,j)=ntsd*i*j/2.
!            if (wavepar_global(i,j).lt.amin) amin=wavepar_global(i,j)
!            if (wavepar_global(i,j).gt.amax) amax=wavepar_global(i,j)
!         enddo
!         enddo
!         print*,"MPI_TEST send min-max wavepar ",amin,amax
!         print*,"MPI_TEST ntsd ",ntsd
!!         call mpi_recv
!!     &     (U10,im*jm,MPI_REAL,2,NTSD,MPI_COMM_WORLD,status,ierr)
!!         print*,"MPI_TEST recv U10 ",U10(im/2,jm/2)
!         call mpi_send
!     &     (wavepar_global,im*jm,MPI_REAL,2,ntsd,MPI_COMM_WORLD,ierr)
!       print*,"MPI_TEST send wavepar_global ",wavepar_global(im/2,jm/2)
!      enddo
C
!      itest=555
!      call mpi_send(itest,1,MPI_INTEGER,2,tag,MPI_COMM_WORLD,ierr)
!       call MPI_BCAST(iwave,1,MPI_INTEGER,0,MPI_COMM_COMP,ierr)
!      call mpi_comm_group(MPI_COMM_WORLD,iworldgroup,ierr)
!      write(6,*)'MPI_TEST iworldgroup ',iworldgroup
!      call mpi_comm_create(MPI_COMM_WORLD,iworldgroup,newcom,ierr)
!      print*,'MPI_TEST newcom ',newcom
!      call MPI_COMM_SPAWN(cmd, MPI_ARGV_NULL, maxpes, info, root, 
!     &                    MPI_COMM_WORLD, INTERCOMM, errCodes, ierr)
!      write(6,*)'TURBL_RECV MPI_COMM_WORLD ',MPI_COMM_WORLD
!      write(6,*)'TURBL_RECV INTERCOMM ',INTERCOMM
!      CALL MPI_COMM_SIZE(INTERCOMM,npesinter,ierr)
!      CALL MPI_COMM_RANK(INTERCOMM,mypeinter,ierr)
!      CALL MPI_COMM_GROUP(INTERCOMM,groupinter,ierr)
!      write(6,*)'TURBL_RECV npesinter,mypeinter,groupinter ',
!     &           npesinter, mypeinter, groupinter
!      call MPI_COMM_REMOTE_SIZE(INTERCOMM, parentsize, ierr)
!      call MPI_COMM_REMOTE_GROUP(INTERCOMM, parentgroup, ierr)
!      write(6,*)'TURBL_RECV parentsize, parentgroup ',
!     &           parentsize, parentgroup
cpk      call mpi_comm_dup(MPI_COMM_WORLD,comdup,ierr)
cpk      print*,'MPI_TEST comdup ',comdup
cpk      call mpi_comm_split(comdup,0,0,mpi_comm_comp,ierr)
cpk      print*,'MPI_TEST mpi_comm_comp ',mpi_comm_comp
!      call mpi_send(i,1,MPI_INTEGER,0,tag,MPI_COMM_WORLD,ierr)
!      iworld = mpi_comm_comp2
cpk      iworld = MPI_COMM_WORLD
cpk      write(6,*)'MPI_TEST iworld',iworld
cpk      call mpi_comm_group(iworld,igroup,ierr)
cpk      write(6,*)'MPI_TEST igroup',igroup
!      call mpi_group_excl(igroup,0,0,igroup_x,ierr)
!!      call mpi_group_excl(igroup,icc,irank,igroup_x,ierr)
!      write(6,*)'TURBL_RECV igroup,icc,irank,igroup_x ',
!     &           igroup,icc,irank,igroup_x
cpk      call mpi_comm_create(iworld,igroup,iworld_minus,ierr)
!      call mpi_comm_create(iworld,0,iwamcomm,ierr)
!      call mpi_comm_create(mpi_comm_comp2,0,iwamcomm,ierr)
!      write(6,*)'TURBL_RECV iwamcomm ',iwamcomm
cpk      write(6,*)'MPI_TEST iworld_minus ',iworld_minus
!      call mpi_group_free(igroup,ierr)
!      write(6,*)'TURBL_RECV igroup',igroup
!      call mpi_group_free(igroup_x,ierr)
cpk      call mpi_intercomm_create(mpi_comm_comp,0,iworld_minus,0,0,
cpk     *   mpi_comm_inter_array(1),ierr)
cpk      write(6,*)'MPI_TEST mpi_comm_inter_array ',mpi_comm_inter_array
cpk      inter=mpi_comm_inter_array(1)
!      call mpi_send(i,1,MPI_INTEGER,0,tag,inter,ierr)
!      call mpi_send(i,1,MPI_INTEGER,0,tag,MPI_COMM_WORLD,ierr)
!      call mpi_intercomm_create(mpi_comm_comp,0,iworld_minus,irlr,0,
!     *   mpi_comm_inter_array(i),ierr)
!      
!
!      if (myrank .EQ. 0) then
!            yourrank=0
!        elseif (myrank .EQ. 1) then
!            yourrank=1
!        endif
!       i=500
!       do i=0,6
!       if (mype .EQ. 0) then
!         source=0
    
!         write(6,*)'TURBL_RECV source ',source
!
!       print*,"MPI_TEST send2TURBL",i
!       call mpi_send(i,1,MPI_INTEGER,1,tag,MPI_COMM_COMP,ierr)
!       call mpi_send(i,1,MPI_INTEGER,0,tag,INTERCOMM,ierr)
!       print*,"MPI_TEST send2TURBL",i
!       call MPI_BCAST(i,1,MPI_INTEGER,0,MPI_COMM_WORLD,ierr)
!       call mpi_recv
!     &     (windfrometa,nx*ny,MPI_REAL,source,tag,INTERCOMM,
!     &      status,ierr)
!       print*, "TURBL_RECV windfrometa=",windfrometa
!       endif
!       enddo
!      call mpi_barrier(MPI_COMM_WORLD,ierr)
! -------------------------------------
!  Finalize MPI
! -------------------------------------
      CALL MPI_FINALIZE(ierr)

      STOP  
      END PROGRAM MPI_TEST
