MPI_Type_Create_Hindexed_Block 生成错误的派生数据类型范围答案

【问题标题】：MPI_Type_Create_Hindexed_Block generates wrong extent of derived datatypeMPI_Type_Create_Hindexed_Block 生成错误的派生数据类型范围
【发布时间】：2017-02-06 17:33:12
【问题描述】：

使用Fortran，我正在尝试为动态分配的结构构建派生数据类型，但它得到了新类型的错误范围，代码如下：

PROGRAM MAIN
IMPLICIT NONE
INCLUDE 'mpif.h'
INTEGER :: I
INTEGER :: MYID,NUMPROCS,IError
INTEGER :: Extent,Size,Disp(2)
INTEGER :: Status(MPI_STATUS_SIZE)
INTEGER :: New_Type, Blocks(3), Types(3), Offsets(3), POS(2)
INTEGER :: POS_(4)
INTEGER :: ElmOffset(3),Send_Type
INTEGER :: M

TYPE Struct    
    INTEGER :: N
    REAL :: A
    REAL :: B(2)
END TYPE Struct
TYPE(Struct),ALLOCATABLE :: Structs(:)

    M=9

    CALL MPI_INIT( IError )
    CALL MPI_COMM_SIZE( MPI_COMM_WORLD, NUMPROCS, IError )
    CALL MPI_COMM_RANK( MPI_COMM_WORLD, MYID,     IError )

    ALLOCATE( Structs(M) )
    DO I=1,M
        Structs(I)%N = I*1000 + MYID
        Structs(I)%A = 250.0_8 + MYID*1.0
        Structs(I)%B(1) = 10.0_8 + MYID*1.0
        Structs(I)%B(2) = 20.0_8 + MYID*1.0
    END DO

    CALL MPI_GET_ADDRESS( Structs(1)%N,    POS_(1), IError )
    CALL MPI_GET_ADDRESS( Structs(1)%A,    POS_(2), IError )
    CALL MPI_GET_ADDRESS( Structs(1)%B(1), POS_(3), IError )
    CALL MPI_GET_ADDRESS( Structs(1)%B(2), POS_(4), IError )
    POS_=POS_ - POS_(1)
    IF (MYID.EQ.0) THEN
        WRITE(*,*) MYID, POS_
    END IF

    Types(1) = MPI_INTEGER
    Types(2) = MPI_DOUBLE_PRECISION
    Types(3) = MPI_DOUBLE_PRECISION

    Offsets(1) = 0
    CALL MPI_GET_ADDRESS( Structs(1)%N, Disp(1), IError )
    CALL MPI_GET_ADDRESS( Structs(1)%A, Disp(2), IError )
    Offsets(2) = Offsets(1) + Blocks(1)*( Disp(2)-Disp(1) )
    Disp(1) = Disp(2)
    CALL MPI_GET_ADDRESS( Structs(1)%B(1), Disp(2), IError )
    Offsets(3) = Offsets(2) + Blocks(2)*( Disp(2)-Disp(1) )

    CALL MPI_TYPE_STRUCT( 3, Blocks, Offsets, Types, New_Type, IError )
    CALL MPI_TYPE_COMMIT( New_Type, IError )

    CALL MPI_TYPE_EXTENT(New_Type, Extent, IError)
    CALL MPI_TYPE_SIZE(New_Type, Size, IError)
    IF (MYID.EQ.0) THEN
        WRITE(*,*) 'New_Type extents = ', Extent
        WRITE(*,*) 'New_Type size = ', Size
    END IF

    CALL MPI_GET_ADDRESS( Structs(1)%N, ElmOffset(1), IError )
    CALL MPI_GET_ADDRESS( Structs(2)%N, ElmOffset(2), IError )
    CALL MPI_GET_ADDRESS( Structs(3)%N, ElmOffset(3), IError )
    ElmOffset=ElmOffset - ElmOffset(1)

    IF (MYID.EQ.0) THEN
        WRITE(*,*) MYID,ElmOffset
    END IF

    CALL MPI_TYPE_CREATE_HINDEXED_BLOCK( 3, 1, ElmOffset, New_Type, Send_Type, IError )
    CALL MPI_TYPE_COMMIT( Send_Type, IError )

    CALL MPI_TYPE_EXTENT( Send_Type, Extent, IError )
    CALL MPI_TYPE_SIZE( Send_Type, Size, IError )

    IF (MYID.EQ.0) THEN
        WRITE(*,*) 'Send_Type extents = ', Extent
        WRITE(*,*) 'Send_Type size = ', Size
    END IF

    CALL MPI_TYPE_FREE(Send_Type,IError)
    CALL MPI_TYPE_FREE(New_Type,IError)
    CALL MPI_FINALIZE(IError)

END PROGRAM MAIN

输出如下：

            POS_ : 0  8  16  24
New_Type Extents : 32
   New_Type Size : 28

以上结果显示没有问题

      ElemOffsets :  0  32  64
Send_Type Extents : -32             <= Problem is here !!! It should be 96
   Send_Type Size :  84

我实际上想使用派生数据类型发送 3 个结构块：Send_Type

IF (MYID.EQ.0) THEN
    DO I=1,(NUMPROCS-1)
         CALL MPI_SEND( Structs(1)%N, 1, Send_Type, I, 0, MPI_COMM_WORLD, IError)       
ELSE
    CALL MPI_RECV( Structs(1)%N, 1, Send_Type, 0, 0, MPI_COMM_WORLD, Status, IError)

END IF

WRITE( (MYID+10),*) Structs(1)%N, Structs(1)%A
WRITE( (MYID+10),*) Structs(1)%B(1), Structs(1)%B(2)

WRITE( (MYID+100),*) Structs(3)%N, Structs(3)%A
WRITE( (MYID+100),*) Structs(3)%B(1), Structs(3)%B(2)

但是，出现错误：程序异常 - 访问冲突

我不知道怎么了... 但肯定是Send_Type没有正确创建

这样的问题怎么解决？

【问题讨论】：

标签： types fortran mpi

【解决方案1】：

问题是由于在 64 位操作系统上，地址的大小大于 32 位整数。因此，函数int MPI_Get_address(const void *location, MPI_Aint *address) 输出一个MPI_Aint，大到足以包含一个地址。实际上，MPI_Aint 可以大于 MPI_INT。

在 Fortran 中，MPI_Aint 写入 INTEGER (KIND=MPI_ADDRESS_KIND)。另请参阅MPI_Aint in MPI_(I)NEIGHBOR_ALLTOALLW() vs int in MPI_(I)ALLTOALLW() 和第 48 页上的MPI Standard 的第 2.5.6 节。

因此，只要涉及地址，就必须使用数据类型 INTEGER (KIND=MPI_ADDRESS_KIND)（对于 POS_、Disp、Offset、Extent 和 ElmOffset）。

一个基于你的更正示例代码，由mpif90 main.f90 -o main -Wall 编译并由mpirun -np 2 main 运行写道：

PROGRAM MAIN
IMPLICIT NONE
INCLUDE 'mpif.h'
INTEGER :: I
INTEGER :: MYID,NUMPROCS,IError
INTEGER :: Size
INTEGER :: Status(MPI_STATUS_SIZE)
INTEGER :: New_Type, Blocks(3), Types(3)
INTEGER :: Send_Type
INTEGER :: M
INTEGER (KIND=MPI_ADDRESS_KIND):: Offsets(3),POS_(4), ElmOffset(3), Disp(2),Extent

TYPE Struct    
    INTEGER :: N
    REAL*8 :: A
    REAL*8 :: B(2)
END TYPE Struct
TYPE(Struct),ALLOCATABLE :: Structs(:)
    WRITE(*,*) 'Size of Integer = ',SIZEOF(M)
    WRITE(*,*) 'Size of Integer (KIND=MPI_ADDRESS_KIND)= ',SIZEOF(Extent)
    M=9

    CALL MPI_INIT( IError )
    CALL MPI_COMM_SIZE( MPI_COMM_WORLD, NUMPROCS, IError )
    CALL MPI_COMM_RANK( MPI_COMM_WORLD, MYID,     IError )

    ALLOCATE( Structs(M) )
    DO I=1,M
        Structs(I)%N = I*1000 + MYID
        Structs(I)%A = 250.0_8 + MYID*1.0
        Structs(I)%B(1) = 10.0_8 + MYID*1.0
        Structs(I)%B(2) = 20.0_8 + MYID*1.0
    END DO

    Blocks(1)=1
    Blocks(2)=1
    Blocks(3)=2

    CALL MPI_GET_ADDRESS( Structs(1)%N,    POS_(1), IError )
    CALL MPI_GET_ADDRESS( Structs(1)%A,    POS_(2), IError )
    CALL MPI_GET_ADDRESS( Structs(1)%B(1), POS_(3), IError )
    CALL MPI_GET_ADDRESS( Structs(1)%B(2), POS_(4), IError )
    POS_=POS_ - POS_(1)
    IF (MYID.EQ.0) THEN
        WRITE(*,*) MYID, POS_
    END IF

    Types(1) = MPI_INTEGER
    Types(2) = MPI_DOUBLE_PRECISION
    Types(3) = MPI_DOUBLE_PRECISION

    Offsets(1) = 0
    CALL MPI_GET_ADDRESS( Structs(1)%N, Disp(1), IError )
    CALL MPI_GET_ADDRESS( Structs(1)%A, Disp(2), IError )
    !Offsets(2) = Offsets(1) + Blocks(1)*( Disp(2)-Disp(1) )
    Offsets(2) = Offsets(1) + ( Disp(2)-Disp(1) )
    Disp(1) = Disp(2)
    CALL MPI_GET_ADDRESS( Structs(1)%B(1), Disp(2), IError )
    !Offsets(3) = Offsets(2) + Blocks(2)*( Disp(2)-Disp(1) )
    Offsets(3) = Offsets(2) + ( Disp(2)-Disp(1) )

    CALL MPI_TYPE_CREATE_STRUCT( 3, Blocks, Offsets, Types, New_Type, IError )
    CALL MPI_TYPE_COMMIT( New_Type, IError )

    CALL MPI_TYPE_GET_EXTENT(New_Type, Extent, IError)
    CALL MPI_TYPE_SIZE(New_Type, Size, IError)
    IF (MYID.EQ.0) THEN
        WRITE(*,*) 'New_Type extents = ', Extent
        WRITE(*,*) 'New_Type size = ', Size
    END IF

    CALL MPI_GET_ADDRESS( Structs(1)%N, ElmOffset(1), IError )
    CALL MPI_GET_ADDRESS( Structs(2)%N, ElmOffset(2), IError )
    CALL MPI_GET_ADDRESS( Structs(3)%N, ElmOffset(3), IError )
    ElmOffset=ElmOffset - ElmOffset(1)

    IF (MYID.EQ.0) THEN
        WRITE(*,*) MYID,ElmOffset
    END IF

    CALL MPI_TYPE_CREATE_HINDEXED_BLOCK( 3, 1, ElmOffset, New_Type, Send_Type, IError )
    CALL MPI_TYPE_COMMIT( Send_Type, IError )

    CALL MPI_TYPE_GET_EXTENT( Send_Type, Extent, IError )
    CALL MPI_TYPE_SIZE( Send_Type, Size, IError )

    IF (MYID.EQ.0) THEN
        WRITE(*,*) 'Send_Type extents = ', Extent
        WRITE(*,*) 'Send_Type size = ', Size
    END IF


    IF (MYID.EQ.0) THEN
        DO I=1,(NUMPROCS-1)
            CALL MPI_SEND( Structs(1)%N, 1, Send_Type, I, 0, MPI_COMM_WORLD, IError)
        END DO       
    ELSE
        CALL MPI_RECV( Structs(1)%N, 1, Send_Type, 0, 0, MPI_COMM_WORLD, Status, IError)

    END IF

    WRITE( (MYID+10),*) Structs(1)%N, Structs(1)%A
    WRITE( (MYID+10),*) Structs(1)%B(1), Structs(1)%B(2)

    WRITE( (MYID+100),*) Structs(3)%N, Structs(3)%A
    WRITE( (MYID+100),*) Structs(3)%B(1), Structs(3)%B(2)

    CALL MPI_TYPE_FREE(Send_Type,IError)
    CALL MPI_TYPE_FREE(New_Type,IError)
    CALL MPI_FINALIZE(IError)

END PROGRAM MAIN

我将REAL :: A 更改为REAL*8 :: A 以删除Structs(I)%A = 250.0_8 + MYID*1.0 线上关于双精度到浮点转换的警告。正如 Hristo Iliev 所注意到的，它与使用 MPI_DOUBLE_PRECISION 的新数据类型一致。

【讨论】：

感谢您的友好解释！数据已成功发送！但输出让我感到困惑：New_Type extent = 160; New_Type 大小 = 208；为什么字体大小大于其范围？
是因为我设置Block的时候出错了...应该是Blocks(1)=1而不是Blocks(1)=SIZEOF(Structs(1)%N)，是块的大小，单位是字节，不是项数在每个块中。对不起！今晚我会纠正这个。
我终于成功了！函数 MPI_TYPE_STRUCT 和 MPI_TYPE_EXTENT 都已弃用！相反，我应该使用 MPI_TYPE_CREATE_STRUCT 和 MPI_TYPE_GET_EXTENT。
将REAL 更改为REAL*8 不仅消除了警告，而且还减少了程序错误，因为类型构造函数中的数据类型被指定为MPI_DOUBLE_PRECISION。此外，MPI_TYPE_STRUCT 和 MPI_TYPE_SIZE 在 MPI-3.0 中被删除
其实默认的REAL KIND设置为8，@HristoIliev

【解决方案2】：

实现你想要的正确方法如下。

1) 创建表示一条记录的结构化数据类型。

CALL MPI_GET_ADDRESS(Structs(1)%N,    POS_(1), IError)
CALL MPI_GET_ADDRESS(Structs(1)%A,    POS_(2), IError)
CALL MPI_GET_ADDRESS(Structs(1)%B(1), POS_(3), IError)
Offsets = POS_ - POS_(1)

Types(1) = MPI_INTEGER
Types(2) = MPI_REAL
Types(3) = MPI_REAL

Blocks(1) = 1
Blocks(2) = 1
Blocks(3) = 2

CALL MPI_TYPE_CREATE_STRUCT(3, Blocks, Offsets, Types, Elem_Type, IError)

此数据类型现在可用于发送该结构的一个记录：

CALL MPI_TYPE_COMMIT(Elem_Type, IError)
CALL MPI_SEND(Structs(1), 1, Elem_Type, ...)

2) 要发送多个记录，首先调整新数据类型的大小（强制其范围为特定大小）以匹配结构的真实大小。这样做是为了说明编译器可能在记录末尾插入的任何填充。

CALL MPI_TYPE_GET_EXTENT(Elem_Type, Lb, Extent, IError)
CALL MPI_GET_ADDRESS(Structs(1)%N, POS_(1), IError)
CALL MPI_GET_ADDRESS(Structs(2)%N, POS_(2), IError)
Extent = POS_(2) - POS_(1)
CALL MPI_TYPE_CREATE_RESIZED(Elem_Type, Lb, Extent, ElemSized_Type, IError)

3) 您现在可以使用新的数据类型来发送结构的多条记录：

CALL MPI_TYPE_COMMIT(ElemSized_Type, IError)
CALL MPI_SEND(Structs(1), 3, ElemSized_Type, ...)

或者，您可以创建一个同时包含三个元素的连续数据类型：

CALL MPI_TYPE_CONTIGUOUS(3, ElemSized_Type, BunchOfElements_Type, IError)
CALL MPI_TYPE_COMMMIT(BunchOfElements_Type, IError)
CALL MPI_SEND(Structs(1), 1, BunchOfElements_Type, ...)

注意：不必提交未用于通信或 I/O 操作的数据类型。

【讨论】：

感谢您的建议！但是还有一个问题...Struct with dynamically allocated memeber