scalapack 中的行分配不一致
Inconsistent rows allocation in scalapack
考虑以下简单的 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
这没有意义,为什么 rank 0 会在块大小 100 和 ranks * 块大小 > N 的情况下得到 200 个元素。
因此,我的程序适用于 mpi 等级 1、2、3、5、7、11,但等级 4、6、8、9、10、12 等失败(我不明白为什么等级 9 失败!) .任何人都可以解释我的方法有什么问题吗?
GFortran 版本:6.1.0
头皮版本:2.1.0
MacOS 版本:10.11
您的代码有很多问题
1) 首先不要使用 Integer( 8 )。正如弗拉基米尔所说,请忘记这一点。它不仅不可移植,因此也是非常糟糕的做法(请在此处查看许多示例,例如 Fortran 90 kind parameter) here it is wrong as numroc
expects an integer of default kind as its first argument (see e.g. 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
The entries in the array dims are set to describe a Cartesian grid
with ndims dimensions and a total of nnodes nodes. The dimensions are
set to be as close to each other as possible,using an appropriate
divisibility algorithm. The caller may further constrain the
operation of this routine by specifying elements of array dims. If
dims[i] is set to a positive number,the routine will not modify the
number of nodes in dimension i; only those entries where dims[i] = 0
are modified by the call.
您将 dims 设置为零,因此例程可以自由设置两个维度。因此,对于 11 个进程,您将获得一个 1x11 或 11x1 的网格,这似乎是您所期望的。但是,对于 12 个进程,如 The dimensions are set to be as close to each other as possible
,您将获得 3x4 或 4x3 网格,而不是 12x1。如果每行是 3x4,你期望 numroc
到 return 3 个进程有 200 个元素(2 个块),1 个进程有 100 个。因为有 3 行,你因此期望 3x3=9 个进程 returning 200 和 3x1=3 returning 100。这就是您所看到的。也可以尝试 15 个进程 - 根据您的说法,您会看到奇数个进程 "does not work",这是因为(高等数学警报)15=3x5。顺便说一下,在我的机器上,9 个进程不是 return 3x3 - 这对我来说看起来像是 openmpi 中的一个错误。
考虑以下简单的 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
这没有意义,为什么 rank 0 会在块大小 100 和 ranks * 块大小 > N 的情况下得到 200 个元素。 因此,我的程序适用于 mpi 等级 1、2、3、5、7、11,但等级 4、6、8、9、10、12 等失败(我不明白为什么等级 9 失败!) .任何人都可以解释我的方法有什么问题吗?
GFortran 版本:6.1.0
头皮版本:2.1.0
MacOS 版本:10.11
您的代码有很多问题
1) 首先不要使用 Integer( 8 )。正如弗拉基米尔所说,请忘记这一点。它不仅不可移植,因此也是非常糟糕的做法(请在此处查看许多示例,例如 Fortran 90 kind parameter) here it is wrong as numroc
expects an integer of default kind as its first argument (see e.g. 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
The entries in the array dims are set to describe a Cartesian grid with ndims dimensions and a total of nnodes nodes. The dimensions are set to be as close to each other as possible,using an appropriate divisibility algorithm. The caller may further constrain the operation of this routine by specifying elements of array dims. If dims[i] is set to a positive number,the routine will not modify the number of nodes in dimension i; only those entries where dims[i] = 0 are modified by the call.
您将 dims 设置为零,因此例程可以自由设置两个维度。因此,对于 11 个进程,您将获得一个 1x11 或 11x1 的网格,这似乎是您所期望的。但是,对于 12 个进程,如 The dimensions are set to be as close to each other as possible
,您将获得 3x4 或 4x3 网格,而不是 12x1。如果每行是 3x4,你期望 numroc
到 return 3 个进程有 200 个元素(2 个块),1 个进程有 100 个。因为有 3 行,你因此期望 3x3=9 个进程 returning 200 和 3x1=3 returning 100。这就是您所看到的。也可以尝试 15 个进程 - 根据您的说法,您会看到奇数个进程 "does not work",这是因为(高等数学警报)15=3x5。顺便说一下,在我的机器上,9 个进程不是 return 3x3 - 这对我来说看起来像是 openmpi 中的一个错误。