1! -*- Mode: Fortran; -*-
2!
3!  (C) 2003 by Argonne National Laboratory.
4!      See COPYRIGHT in top-level directory.
5!
6! This test contributed by Kim McMahon, Cray
7!
8      program main
9      use mpi
10      implicit none
11      integer ierr, i, j, type, count,errs
12      parameter (count = 4)
13      integer rank, size, xfersize
14      integer status(MPI_STATUS_SIZE)
15      integer blocklens(count), displs(count)
16      double precision,dimension(:,:),allocatable :: sndbuf, rcvbuf
17      logical verbose
18
19      verbose = .false.
20      call mtest_init ( ierr )
21      call mpi_comm_size( MPI_COMM_WORLD, size, ierr )
22      call mpi_comm_rank( MPI_COMM_WORLD, rank, ierr )
23      if (size .lt. 2) then
24         print *, "Must have at least 2 processes"
25         call MPI_Abort( 1, MPI_COMM_WORLD, ierr )
26      endif
27
28      errs = 0
29      allocate(sndbuf(7,100))
30      allocate(rcvbuf(7,100))
31
32      do j=1,100
33        do i=1,7
34           sndbuf(i,j) = (i+j) * 1.0
35         enddo
36      enddo
37
38      do i=1,count
39         blocklens(i) = 7
40      enddo
41
42! bug occurs when first two displacements are 0
43      displs(1) = 0
44      displs(2) = 0
45      displs(3) = 10
46      displs(4) = 10
47
48      call mpi_type_indexed( count, blocklens, displs*blocklens(1),  &
49      &                         MPI_DOUBLE_PRECISION, type, ierr )
50
51      call mpi_type_commit( type, ierr )
52
53! send using this new type
54
55      if (rank .eq. 0) then
56
57          call mpi_send( sndbuf(1,1), 1, type, 1, 0, MPI_COMM_WORLD,ierr )
58
59      else if (rank .eq. 1) then
60
61          xfersize=count * blocklens(1)
62          call mpi_recv( rcvbuf(1,1), xfersize, MPI_DOUBLE_PRECISION, 0, 0, &
63           &   MPI_COMM_WORLD,status, ierr )
64
65
66! Values that should be sent
67
68        if (verbose) then
69!       displacement = 0
70            j=1
71            do i=1, 7
72               print*,'sndbuf(',i,j,') = ',sndbuf(i,j)
73            enddo
74
75!       displacement = 10
76            j=11
77            do i=1,7
78               print*,'sndbuf(',i,j,') = ',sndbuf(i,j)
79            enddo
80            print*,' '
81
82! Values received
83            do j=1,count
84                do i=1,7
85                    print*,'rcvbuf(',i,j,') = ',rcvbuf(i,j)
86                enddo
87            enddo
88        endif
89
90! Error checking
91        do j=1,2
92           do i=1,7
93             if (rcvbuf(i,j) .ne. sndbuf(i,1)) then
94                print*,'ERROR in rcvbuf(',i,j,')'
95                print*,'Received ', rcvbuf(i,j),' expected ',sndbuf(i,11)
96                errs = errs+1
97             endif
98           enddo
99        enddo
100
101        do j=3,4
102           do i=1,7
103              if (rcvbuf(i,j) .ne. sndbuf(i,11)) then
104                print*,'ERROR in rcvbuf(',i,j,')'
105                print*,'Received ', rcvbuf(i,j),' expected ',sndbuf(i,11)
106                errs = errs+1
107              endif
108           enddo
109        enddo
110      endif
111!
112      call mpi_type_free( type, ierr )
113      call mtest_finalize( errs )
114      call mpi_finalize( ierr )
115
116      end
117