scalapack中不一致的行分配

Question

考虑以下简单的fortran程序

program test_vec_allocation
    use mpi
    implicit none
    integer(kind=8)             :: N
    ! =========================BLACS and MPI=======================
    integer                     :: ierr, size, rank,dims(2)
    ! -------------------------------------------------------------
    integer, parameter          :: block_size = 100
    integer                     :: context, nprow, npcol, local_nprow, local_npcol
    integer                     :: numroc, indxl2g, descmat(9),descvec(9)
    integer                     :: mloc_mat ,nloc_mat ,mloc_vec ,nloc_vec

    call blacs_pinfo(rank,size)
    dims=0
    call MPI_Dims_create(size, 2, dims, ierr)
    nprow = dims(1);npcol = dims(2)
    call blacs_get(0,0,context)
    call blacs_gridinit(context, 'R', nprow, npcol)
    call blacs_gridinfo(context, nprow, npcol, local_nprow,local_npcol)

    N = 700

    mloc_vec = numroc(N,block_size,local_nprow,0, nprow)
    nloc_vec = numroc(1,block_size,local_npcol,0, npcol)
    print *,"Rank", rank, mloc_vec, nloc_vec

    call blacs_gridexit(context)
    call blacs_exit(0)

end program test_vec_allocation

当我用11个mpi级别运行它时，我得到

 Rank           0         100           1
 Rank           4         100           1
 Rank           2         100           1
 Rank           1         100           1
 Rank           3         100           1
 Rank          10           0           1
 Rank           6         100           1
 Rank           5         100           1
 Rank           9           0           1
 Rank           8           0           1
 Rank           7           0           1

这就是我期望scalapack对这个数组进行分割的方式，但是，对于我所得到的偶数级别：

 Rank           0         200           1
 Rank           8         200           0
 Rank           9         100           1
 Rank          10         100           0
 Rank           1         200           0
 Rank           6         200           1
 Rank          11         100           0
 Rank           3         200           1
 Rank           4         200           0
 Rank           2         200           0
 Rank           7         200           0
 Rank           5         200           0

这是没有意义的，为什么0将获得200个元素的块大小100和秩*块大小> N，因此，我的程序工作为mpi排名1,2,3,5,7,11，但失败的排名4,6,8,9,10,12等(我没有为什么它是失败的排名9!)。有人能解释一下我的做法有什么不对吗？

GFortran版本: 6.1.0

SCALPACK版本: 2.1.0

MacOS版本: 10.11

Answer 1

您的代码有许多错误。

第一，不要使用整数(8 )。正如弗拉迪米尔所说，请不要学到这一点。它不仅不是可移植的，因此也是非常糟糕的实践(请看这里的许多例子，例如Fortran 90 kind parameter)，这是错误的，因为numroc期望一个默认类型的整数作为它的第一个参数(参见https://software.intel.com/content/www/us/en/develop/documentation/mkl-developer-reference-fortran/top/scalapack-routines/scalapack-utility-functions-and-routines/numroc.html)

2)在调用MPI_Init之前调用MPI例程，这会导致未定义的行为。注意，https://www.netlib.org/blacs/BLACS/QRef.html#BLACS_PINFO中的描述没有提到实际调用MPI_Init。因此，我也更喜欢调用MPI_Finalise

3)你误解了MPI_Dims_create。你似乎假设你会得到一个一维分布，但实际上你要求得到一个二维分布。引用https://www.mpi-forum.org/docs/mpi-3.1/mpi31-report.pdf的标准

数组dims中的条目被设置为描述具有ndims维数和总节点数的笛卡尔网格。使用适当的可分性算法，将维度设置为尽可能接近彼此。调用方可以通过指定数组dims的元素来进一步限制此例程的操作。如果dimsi设置为正数，则例程将不会修改维度i中的节点数；只有那些dimsi =0的条目被调用修改。

您将dims设置为零，因此例程可以自由设置两个维度。因此，对于11个进程，您将得到一个1x11或11x1网格，这似乎是您所期望的。但是，对于12个进程，作为The dimensions are set to be as close to each other as possible，您将得到一个3x4或4x3网格，而不是12x1。如果每一行都是3x4，那么您期望numroc返回包含200个元素(2个块)的3个进程，以及一个带100个元素的进程。因此，由于有3行，您希望3x3=9进程返回200，3x1=3返回100。这就是你所看到的。还可以尝试15个procs -您将看到一些根据您“不工作”的进程，这是因为(高级数学警报) 15=3x5。顺便说一句，在我的机器上，9进程不返回3x3 --在我看来，这就像openmpi中的一个bug。