MPI_REDUCE returns 1000x1000x6 REAL 数组的错误答案

Question

我在gfortran大气传输模型代码中用MPI_REDUCE调用了几天，输入参数很好，但在master的recvbuf中返回了非常不合理的结果。我已经能够在一个简单的示例中重现该问题，如下所示：

PROGRAM TEST

    USE mpi

    IMPLICIT NONE

    INTEGER my_rank, size, ierror
    INTEGER, PARAMETER :: nx=1000, ny=1000, nz=6
    INTEGER :: buffsize

    REAL, DIMENSION(nx, ny, nz) :: u, v

    call MPI_INIT(ierror)
    call MPI_COMM_SIZE(MPI_COMM_WORLD, size, ierror)
    call MPI_COMM_RANK(MPI_COMM_WORLD, my_rank, ierror)

    PRINT *, 'my_rank, size: ', my_rank, size

    buffsize = nx*ny*nz

    u = my_rank + 1

    PRINT *, 'PE: ', my_rank, ', Before reduce, SHAPE(u): ', SHAPE(u)
    PRINT *, 'PE: ', my_rank, ', Before reduce, SUM(u): ', SUM(u)

    CALL MPI_REDUCE(u, v, buffsize, MPI_REAL, &
&                   MPI_SUM, 0, MPI_COMM_WORLD, ierror)

    CALL MPI_BARRIER(MPI_COMM_WORLD, ierror)

    PRINT *, 'PE: ', my_rank, ', after reduce, ierror: ', ierror
    PRINT *, 'PE: ', my_rank, ', after reduce, SUM(u): ', SUM(u)
    PRINT *, 'PE: ', my_rank, ', after reduce, SUM(v): ', SUM(v)

    CALL MPI_FINALIZE(ierror)

END PROGRAM test

它returns:

mpirun -np 2 ./test3
 my_rank, size:            0           2
 my_rank, size:            1           2
 PE:            1 , Before reduce, SHAPE(u):         1000        1000           6
 PE:            0 , Before reduce, SHAPE(u):         1000        1000           6
 PE:            0 , Before reduce, SUM(u):    6000000.00    
 PE:            1 , Before reduce, SUM(u):    12000000.0    
 PE:            0 , after reduce, ierror:            0
 PE:            1 , after reduce, ierror:            0
 PE:            1 , after reduce, SUM(u):    12000000.0    
 PE:            0 , after reduce, SUM(u):    6000000.00    
 PE:            1 , after reduce, SUM(v):    0.00000000    
 PE:            0 , after reduce, SUM(v):    18407592.0

PE0 "should" 在最后一行显示 18000000.0 作为 SUM(v)。

如果我将代码中的 nz 参数设置为 6 到 5，则运行会产生正确的结果。真正令人困惑的是，它以这种方式运行，在 a) 带有 gfortran 5.3 和 openmpi 的 AWS EC2 实例，b) 我笔记本电脑的带有 mpich 的 gfortran 5.4，以及 c) 带有 openmpi 的工作站的 gfortran 4.4 上返回相同的减少值总和。

如果我将数组类型更改为 DOUBLE PRECISION（并在 MPI_REDUCE 调用中指定）它工作正常，即使对于更大的数组也是如此。如果我使用 REAL4 而不是 REAL，它会产生同样糟糕的结果。

我知道这一定很简单，我在这里真是个白痴，但我就是不明白这一点。我读过一些建议，我的缓冲区大小需要是一个小于 2^31-1 的整数值，但这里确实是这样。

Answer 1

这与MPI无关，只是求和精度问题：

PROGRAM TEST
    IMPLICIT NONE
    INTEGER, PARAMETER :: nx=1000, ny=1000, nz=6
    REAL, DIMENSION(nx, ny, nz) :: u
    u = 3
    PRINT *, SUM(u)
END PROGRAM test

Returns同样的结果。如果将一个大数加到一个小数上，可能会出现舍入问题，在许多小数的总和中，这种影响可能会累积成一个重大错误。有一些求和算法可以防止这种影响，比如 Kahan summation，显然 Fortran 的 SUM 不是这样实现的。

MPI_REDUCE returns 1000x1000x6 REAL 数组的错误答案

MPI_REDUCE returns wrong answer for 1000x1000x6 REAL array

precision

fortran

sum

mpi