diff --git a/docs/make.jl b/docs/make.jl index 600cd2a4f..7047b36be 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -51,6 +51,7 @@ makedocs( "library.md", "environment.md", "comm.md", + "buffers.md", "pointtopoint.md", "collective.md", "onesided.md", diff --git a/docs/src/advanced.md b/docs/src/advanced.md index bda44f533..c7e8d066e 100644 --- a/docs/src/advanced.md +++ b/docs/src/advanced.md @@ -6,14 +6,6 @@ MPI.free ``` -## Buffers - -```@docs -MPI.Buffer -MPI.Buffer_send -MPI.MPIPtr -``` - ## Datatype objects ```@docs diff --git a/docs/src/buffers.md b/docs/src/buffers.md new file mode 100644 index 000000000..f6aab9c8c --- /dev/null +++ b/docs/src/buffers.md @@ -0,0 +1,13 @@ +# Buffers + +Buffers are used for sending and receiving data. MPI.jl provides the following buffer types: + +```@docs +MPI.IN_PLACE +MPI.Buffer +MPI.Buffer_send +MPI.UBuffer +MPI.VBuffer +MPI.RBuffer +MPI.MPIPtr +``` diff --git a/docs/src/collective.md b/docs/src/collective.md index afb72e87f..0ef08e836 100644 --- a/docs/src/collective.md +++ b/docs/src/collective.md @@ -10,6 +10,7 @@ MPI.Barrier ```@docs MPI.Bcast! +MPI.bcast ``` ## Gather/Scatter @@ -31,9 +32,7 @@ MPI.Gatherv ```@docs MPI.Scatter! -MPI.Scatter MPI.Scatterv! -MPI.Scatterv ``` ### All-to-all diff --git a/src/MPI.jl b/src/MPI.jl index 0e940cf15..76052ffec 100644 --- a/src/MPI.jl +++ b/src/MPI.jl @@ -4,7 +4,7 @@ using Libdl, Serialization using Requires using DocStringExtensions -export mpiexec +export mpiexec, UBuffer, VBuffer function serialize(x) s = IOBuffer() diff --git a/src/buffers.jl b/src/buffers.jl index c824c460d..78bee906f 100644 --- a/src/buffers.jl +++ b/src/buffers.jl @@ -57,11 +57,42 @@ Additionally, certain sentinel values can be used, e.g. `MPI_IN_PLACE` or `MPI_B """ MPIPtr +# MPI_IN_PLACE + +struct InPlace +end +Base.cconvert(::Type{MPIPtr}, ::InPlace) = MPI_IN_PLACE + + +""" + MPI.IN_PLACE + +A sentinel value that can be passed as a buffer argument for certain collective operations +to use the same buffer for send and receive operations. + +- [`Scatter!`](@ref) and [`Scatterv!`](@ref): can be used as the `recvbuf` argument on the + root process. + +- [`Gather!`](@ref) and [`Gatherv!`](@ref): can be used as the `sendbuf` argument on the + root process. + +- [`Allgather!`](@ref), [`Allgatherv!`](@ref), [`Alltoall!`](@ref) and + [`Alltoallv!`](@ref): can be used as the `sendbuf` argument on all processes. + +- [`Reduce!`](@ref) (root only), [`Allreduce!`](@ref), [`Scan!`](@ref) and + [`Exscan!`](@ref): can be used as `sendbuf` argument. + +""" +const IN_PLACE = InPlace() + +# TODO: MPI_BOTTOM + """ MPI.Buffer -An MPI buffer for communication operations. +An MPI buffer for communication with a single rank. It is used for point-to-point +communication and some collective operations. # Fields $(DocStringExtensions.FIELDS) @@ -83,6 +114,10 @@ and `datatype`. Methods are provided for - `SubArray`s of an `Array` or `CUDA.CuArray` where the layout is contiguous, sequential or blocked. +# See also + +- [`Buffer_send`](@ref) + """ struct Buffer{A} """a Julia object referencing a region of memory to be used for communication. It is @@ -125,6 +160,9 @@ function Buffer(sub::SubArray{T,N,P,I,false}) where {T,N,P,I<:Tuple{Vararg{Union Buffer(parent(sub), Cint(1), datatype) end +Buffer(::InPlace) = Buffer(IN_PLACE, 0, DATATYPE_NULL) +Buffer(::Nothing) = Buffer(nothing, 0, DATATYPE_NULL) + """ Buffer_send(data) @@ -133,7 +171,194 @@ Construct a [`Buffer`](@ref) object for a send operation from `data`, allowing c """ Buffer_send(data) = isbits(data) ? Buffer(Ref(data)) : Buffer(data) Buffer_send(str::String) = Buffer(str, sizeof(str), MPI.CHAR) +Buffer_send(::InPlace) = Buffer(InPlace()) +Buffer_send(::Nothing) = Buffer(nothing) + + + + + + +""" + MPI.UBuffer + +An MPI buffer for chunked collective communication, where all chunks are of uniform size. + +# Fields +$(DocStringExtensions.FIELDS) + +# Usage + + UBuffer(data, count::Integer, nchunks::Union{Nothing, Integer}, datatype::Datatype) + +Generic constructor. + + UBuffer(data, count::Integer) + +Construct a `UBuffer` backed by `data`, where `count` is the number of elements in each chunk. + +# See also + +- [`VBuffer`](@ref): similar, but supports chunks of non-uniform sizes. +""" +struct UBuffer{A} + """A Julia object referencing a region of memory to be used for communication. It is + required that the object can be `cconvert`ed to an [`MPIPtr`](@ref).""" + data::A + + """The number of elements of `datatype` in each chunk.""" + count::Cint + + """The maximum number of chunks stored in the buffer. This is used only for + validation, and can be set to `nothing` to disable checks.""" + nchunks::Union{Nothing,Cint} + + """The [`MPI.Datatype`](@ref) stored in the buffer.""" + datatype::Datatype +end +UBuffer(data, count::Integer, nchunks::Union{Integer, Nothing}, datatype::Datatype) = + UBuffer(data, Cint(count), nchunks isa Integer ? Cint(nchunks) : nothing, datatype) + +function UBuffer(arr::AbstractArray, count::Integer) + @assert stride(arr, 1) == 1 + UBuffer(arr, count, div(length(arr), count), Datatype(eltype(arr))) +end +Base.similar(buf::UBuffer) = + UBuffer(similar(buf.data), buf.count, buf.nchunks, buf.datatype) + +UBuffer(::Nothing) = UBuffer(nothing, 0, nothing, DATATYPE_NULL) +UBuffer(::InPlace) = UBuffer(IN_PLACE, 0, nothing, DATATYPE_NULL) + + + +""" + MPI.VBuffer + +An MPI buffer for chunked collective communication, where chunks can be of different sizes and at different offsets. + + +# Fields +$(DocStringExtensions.FIELDS) + +# Usage + + VBuffer(data, counts[, displs[, datatype]]) + +Construct a `VBuffer` backed by `data`, where `counts[j]` is the number of elements in the +`j`th chunk, and `displs[j]` is the 0-based displacement. In other words, the `j`th chunk +occurs in indices `displs[j]+1:displs[j]+counts[j]`. + +The default value for `displs[j] = sum(counts[1:j-1])`. + +# See also + +- [`UBuffer`](@ref) when chunks are all of the same size. +""" +struct VBuffer{A} + """A Julia object referencing a region of memory to be used for communication. It is + required that the object can be `cconvert`ed to an [`MPIPtr`](@ref).""" + data::A + + """An array containing the length of each chunk.""" + counts::Vector{Cint} + + """An array containing the (0-based) displacements of each chunk.""" + displs::Vector{Cint} + + """The [`MPI.Datatype`](@ref) stored in the buffer.""" + datatype::Datatype +end +VBuffer(data, counts, displs, datatype::Datatype) = + VBuffer(data, convert(Vector{Cint}, counts), + convert(Vector{Cint}, displs), datatype) +VBuffer(data, counts, displs) = + VBuffer(data, counts, displs, Datatype(eltype(data))) + +function VBuffer(arr::AbstractArray, counts) + @assert stride(arr,1) == 1 + counts = convert(Vector{Cint}, counts) + displs = similar(counts) + d = zero(Cint) + for i in eachindex(displs) + displs[i] = d + d += counts[i] + end + @assert length(arr) >= d + VBuffer(arr, counts, displs, Datatype(eltype(arr))) +end + +VBuffer(::Nothing) = VBuffer(nothing, Cint[], Cint[], DATATYPE_NULL) +VBuffer(::InPlace) = VBuffer(IN_PLACE, Cint[], Cint[], DATATYPE_NULL) + + +""" + MPI.RBuffer + +An MPI buffer for reduction operations ([`MPI.Reduce!`](@ref), [`MPI.Allreduce!`](@ref), [`MPI.Scan!`](@ref), [`MPI.Exscan!`](@ref)). + +# Fields +$(DocStringExtensions.FIELDS) +# Usage + + RBuffer(senddata, recvdata[, count, datatype]) + +Generic constructor. + + RBuffer(senddata, recvdata) + +Construct a `Buffer` backed by `senddata` and `recvdata`, automatically determining the +appropriate `count` and `datatype`. + +- `senddata` can be [`MPI.IN_PLACE`](@ref) +- `recvdata` can be `nothing` on a non-root node with [`MPI.Reduce!`](@ref) +""" +struct RBuffer{S,R} + """A Julia object referencing a region of memory to be used for the send buffer. It is + required that the object can be `cconvert`ed to an [`MPIPtr`](@ref).""" + senddata::S + + """A Julia object referencing a region of memory to be used for the receive buffer. It is + required that the object can be `cconvert`ed to an [`MPIPtr`](@ref).""" + recvdata::R + + """the number of elements of `datatype` in the buffer. Note that this may not + correspond to the number of elements in the array if derived types are used.""" + count::Cint + + """the [`MPI.Datatype`](@ref) stored in the buffer.""" + datatype::Datatype +end + +RBuffer(senddata, recvdata, count::Integer, datatype::Datatype) = + RBuffer(senddata, recvdata, Cint(count), datatype) + +function RBuffer(senddata::AbstractArray{T}, recvdata::AbstractArray{T}) where {T} + @assert (count = length(senddata)) == length(recvdata) + @assert stride(senddata,1) == stride(recvdata,1) == 1 + RBuffer(senddata, recvdata, count, Datatype(T)) +end +function RBuffer(::InPlace, recvdata::AbstractArray{T}) where {T} + count = length(recvdata) + @assert stride(recvdata,1) == 1 + RBuffer(IN_PLACE, recvdata, count, Datatype(T)) +end +function RBuffer(senddata::AbstractArray{T}, recvdata::Nothing) where {T} + count = length(senddata) + @assert stride(senddata,1) == 1 + RBuffer(senddata, nothing, count, Datatype(T)) +end + +function RBuffer(senddata::Ref{T}, recvdata::Ref{T}) where {T} + RBuffer(senddata, recvdata, 1, Datatype(T)) +end +function RBuffer(senddata::InPlace, recvdata::Ref{T}) where {T} + RBuffer(IN_PLACE, recvdata, 1, Datatype(T)) +end +function RBuffer(senddata::Ref{T}, recvdata::Nothing) where {T} + RBuffer(senddata, nothing, 1, Datatype(T)) +end -const BUFFER_NULL = Buffer(C_NULL, 0, DATATYPE_NULL) +Base.eltype(rbuf::RBuffer) = eltype(rbuf.senddata) +Base.eltype(rbuf::RBuffer{InPlace}) = eltype(rbuf.recvdata) diff --git a/src/collective.jl b/src/collective.jl index ca4f717b6..3080a898f 100644 --- a/src/collective.jl +++ b/src/collective.jl @@ -1,5 +1,3 @@ -const IN_PLACE = MPI_IN_PLACE - """ Barrier(comm::Comm) @@ -19,28 +17,37 @@ function Barrier(comm::Comm) end """ - Bcast!(buf[, count=length(buf)], root::Integer, comm::Comm) + Bcast!(buf, root::Integer, comm::Comm) + +Broadcast the buffer `buf` from `root` to all processes in `comm`. -Broadcast the first `count` elements of the buffer `buf` from `root` to all processes. +# See also +- [`bcast`](@ref) # External links $(_doc_external("MPI_Bcast")) """ -function Bcast!(buffer, count::Integer, - root::Integer, comm::Comm) - +function Bcast!(buf::Buffer, root::Integer, comm::Comm) # int MPI_Bcast(void* buffer, int count, MPI_Datatype datatype, int root, # MPI_Comm comm) @mpichk ccall((:MPI_Bcast, libmpi), Cint, (MPIPtr, Cint, MPI_Datatype, Cint, MPI_Comm), - buffer, count, Datatype(eltype(buffer)), root, comm) - buffer + buf.data, buf.count, buf.datatype, root, comm) + return buf.data end - -function Bcast!(buffer::Union{Ref, AbstractArray}, root::Integer, comm::Comm) - Bcast!(buffer, length(buffer), root, comm) +function Bcast!(data, root::Integer, comm::Comm) + Bcast!(Buffer(data), root, comm) end +""" + bcast(obj, root::Integer, comm::Comm) + +Broadcast the object `obj` from rank `root` to all processes on `comm`. This is able to handle arbitrary data. + +# See also + +- [`Bcast!`](@ref) +""" function bcast(obj, root::Integer, comm::Comm) isroot = Comm_rank(comm) == root count = Ref{Cint}() @@ -56,168 +63,130 @@ function bcast(obj, root::Integer, comm::Comm) if !isroot obj = MPI.deserialize(buf) end - obj + return obj end """ - Scatter!(sendbuf, recvbuf[, count=length(recvbuf)], root::Integer, comm::Comm) + Scatter!(sendbuf::Union{UBuffer,Nothing}, recvbuf, root::Integer, comm::Comm) -Splits the buffer `sendbuf` in the `root` process into `Comm_size(comm)` chunks of length -`count`, and sends the `j`-th chunk to the process of rank `j` into the `recvbuf` buffer. +Splits the buffer `sendbuf` in the `root` process into `Comm_size(comm)` chunks, +sending the `j`-th chunk to the process of rank `j-1` into the `recvbuf` buffer. -`sendbuf` on the root process should be a buffer of length `count*Comm_size(comm)`, and -on non-root processes it is ignored and can be `nothing`. +`sendbuf` on the root process should be a [`UBuffer`](@ref) (an `Array` can also be +passed directly if the sizes can be determined from `recvbuf`). On non-root processes it +is ignored, and `nothing` can be passed instead. -`recvbuf` can be `nothing` on the root process, in which case it is unmodified (this -corresponds the behaviour of `MPI_IN_PLACE` in `MPI_Scatter`). For example +`recvbuf` is a [`Buffer`](@ref) object, or any object for which `Buffer(recvbuf)` is +defined. On the root process, it can also be [`MPI.IN_PLACE`](@ref), in which case it is +unmodified. For example: ``` if root == MPI.Comm_rank(comm) - Scatter!(buf, nothing, count, root, comm) + MPI.Scatter!(UBuffer(buf, count), MPI.IN_PLACE, root, comm) else - Scatter!(nothing, buf, count, root, comm) + MPI.Scatter!(nothing, buf, root, comm) end ``` -`count` should be the same for all processes. - # See also -- [`Scatter`](@ref) to allocate the output buffer. - [`Scatterv!`](@ref) if the number of elements varies between processes. # External links $(_doc_external("MPI_Scatter")) """ -function Scatter!(sendbuf, recvbuf, count::Integer, root::Integer, comm::Comm) - isroot = Comm_rank(comm) == root - if isroot - @assert sendbuf !== nothing - @assert_minlength sendbuf count*Comm_size(comm) - if recvbuf === nothing - recvbuf = MPI_IN_PLACE - end +function Scatter!(sendbuf::UBuffer, recvbuf::Buffer, root::Integer, comm::Comm) + if sendbuf.nchunks !== nothing && Comm_rank(comm) == root + @assert sendbuf.nchunks >= Comm_size(comm) end - @assert_minlength recvbuf count - T = recvbuf isa SentinelPtr ? eltype(sendbuf) : eltype(recvbuf) # int MPI_Scatter(const void* sendbuf, int sendcount, MPI_Datatype sendtype, # void* recvbuf, int recvcount, MPI_Datatype recvtype, int root, # MPI_Comm comm) @mpichk ccall((:MPI_Scatter, libmpi), Cint, - (MPIPtr, Cint, MPI_Datatype, MPIPtr, Cint, MPI_Datatype, Cint, MPI_Comm), - sendbuf, count, Datatype(T), recvbuf, count, Datatype(T), root, comm) - recvbuf -end - -function Scatter!(sendbuf::AbstractArray, recvbuf::Union{Ref,AbstractArray}, root::Integer, comm::Comm) - Scatter!(sendbuf, recvbuf, length(recvbuf), root, comm) + (MPIPtr, Cint, MPI_Datatype, + MPIPtr, Cint, MPI_Datatype, Cint, MPI_Comm), + sendbuf.data, sendbuf.count, sendbuf.datatype, + recvbuf.data, recvbuf.count, recvbuf.datatype, root, comm) + return recvbuf.data end +Scatter!(sendbuf::UBuffer, recvbuf, root::Integer, comm::Comm) = + Scatter!(sendbuf, Buffer(recvbuf), root, comm) +Scatter!(sendbuf::Nothing, recvbuf, root::Integer, comm::Comm) = + Scatter!(UBuffer(nothing), recvbuf, root, comm) -""" - Scatter(sendbuf, count, root, comm) - -Splits the buffer `sendbuf` in the `root` process into `Comm_size(comm)` chunks -and sends the `j`-th chunk to the process of rank `j`, allocating the output buffer. - -# See also -- [`Scatter!`](@ref) for the mutating operation. -- [`Scatterv!`](@ref) if the number of elements varies between processes. +# determine UBuffer count from recvbuf +Scatter!(sendbuf::AbstractArray, recvbuf::Union{Ref,AbstractArray}, root::Integer, comm::Comm) = + Scatter!(UBuffer(sendbuf,length(recvbuf)), recvbuf, root, comm) -# External links -$(_doc_external("MPI_Scatter")) """ -function Scatter(sendbuf, count::Integer, root::Integer, comm::Comm) - Scatter!(sendbuf, similar(sendbuf, count), count, root, comm) -end + Scatterv!(sendbuf::Union{VBuffer,Nothing}, recvbuf, root, comm) +Splits the buffer `sendbuf` in the `root` process into `Comm_size(comm)` chunks and sends +the `j`th chunk to the process of rank `j-1` into the `recvbuf` buffer. -""" - Scatterv!(sendbuf, recvbuf, counts, root, comm) - -Splits the buffer `sendbuf` in the `root` process into `Comm_size(comm)` chunks -of length `counts[j]` and sends the j-th chunk to the process of rank j into the -`recvbuf` buffer, which must be of length at least `count`. +`sendbuf` on the root process should be a [`VBuffer`](@ref). On non-root processes it is +ignored, and `nothing` can be passed instead. -`recvbuf` can be `nothing` on the root process, in which case it is unmodified (this -corresponds the behaviour of `MPI_IN_PLACE` in `MPI_Scatterv`). For example +`recvbuf` is a [`Buffer`](@ref) object, or any object for which `Buffer(recvbuf)` is +defined. On the root process, it can also be [`MPI.IN_PLACE`](@ref), in which case it is +unmodified. For example: ``` if root == MPI.Comm_rank(comm) - Scatterv!(buf, nothing, counts, root, comm) + MPI.Scatterv!(VBuffer(buf, counts), MPI.IN_PLACE, root, comm) else - Scatterv!(nothing, buf, counts, root, comm) + MPI.Scatterv!(nothing, buf, root, comm) end ``` # See also -- [`Scatterv`](@ref) for the allocating operation -- [`Scatter!`](@ref)/[`Scatter`](@ref) if the counts are the same for all processes +- [`Scatter!`](@ref) if the number of elements are the same for all processes # External links $(_doc_external("MPI_Scatterv")) """ -function Scatterv!(sendbuf, recvbuf, counts::Vector, root::Integer, comm::Comm) - rank = Comm_rank(comm) - isroot = rank == root - if isroot - @assert sendbuf !== nothing - @assert_minlength sendbuf sum(counts) - end - if recvbuf === nothing - recvbuf = MPI_IN_PLACE +function Scatterv!(sendbuf::VBuffer, recvbuf::Buffer, root::Integer, comm::Comm) + if Comm_rank(comm) == root + @assert length(sendbuf.counts) >= Comm_size(comm) end - @assert_minlength recvbuf counts[rank+1] - T = recvbuf isa SentinelPtr ? eltype(sendbuf) : eltype(recvbuf) - recvcnt = counts[Comm_rank(comm) + 1] - disps = accumulate(+, counts) - counts # int MPI_Scatterv(const void* sendbuf, const int sendcounts[], # const int displs[], MPI_Datatype sendtype, void* recvbuf, # int recvcount, MPI_Datatype recvtype, int root, MPI_Comm comm) @mpichk ccall((:MPI_Scatterv, libmpi), Cint, - (MPIPtr, Ptr{Cint}, Ptr{Cint}, MPI_Datatype, MPIPtr, Cint, MPI_Datatype, Cint, MPI_Comm), - sendbuf, counts, disps, Datatype(T), recvbuf, recvcnt, Datatype(T), root, comm) - recvbuf + (MPIPtr, Ptr{Cint}, Ptr{Cint}, MPI_Datatype, + MPIPtr, Cint, MPI_Datatype, Cint, MPI_Comm), + sendbuf.data, sendbuf.counts, sendbuf.displs, sendbuf.datatype, + recvbuf.data, recvbuf.count, recvbuf.datatype, + root, comm) + return recvbuf.data end +Scatterv!(sendbuf::VBuffer, recvbuf, root::Integer, comm::Comm) = + Scatterv!(sendbuf, Buffer(recvbuf), root, comm) +Scatterv!(sendbuf::Nothing, recvbuf, root::Integer, comm::Comm) = + Scatterv!(VBuffer(nothing), recvbuf, root, comm) -""" - Scatterv(sendbuf, counts, root, comm) - -Splits the buffer `sendbuf` in the `root` process into `Comm_size(comm)` chunks -of length `counts[j]` and sends the j-th chunk to the process of rank j, which -allocates the output buffer - -# See also -- [`Scatterv!`](@ref) for the mutating operation -- [`Scatter!`](@ref)/[`Scatter`](@ref) if the counts are the same for all processes -# External links -$(_doc_external("MPI_Scatterv")) """ -function Scatterv(sendbuf, counts::Vector, root::Integer, comm::Comm) - recvbuf = similar(sendbuf, counts[Comm_rank(comm) + 1]) - Scatterv!(sendbuf, recvbuf, counts, root, comm) -end - + Gather!(sendbuf, recvbuf::Union{UBuffer,Nothing}, root::Integer, comm::Comm) -""" - Gather!(sendbuf, recvbuf[, count::Integer=length(sendbuf)], root::Integer, comm::Comm) +Each process sends the contents of the buffer `sendbuf` to the `root` process. The `root` +process stores elements in rank order in the buffer buffer `recvbuf`. -Each process sends the first `count` elements of the buffer `sendbuf` to the -`root` process. The `root` process stores elements in rank order in the buffer -buffer `recvbuf`. +`sendbuf` should be a [`Buffer`](@ref) object, or any object for which +[`Buffer_send`](@ref) is defined, with the same length on all processes, and should be the +same length on all processes. -`sendbuf` can be `nothing` on the root process, in which case the corresponding entries in -`recvbuf` are assumed to be already in place (this corresponds the behaviour of -`MPI_IN_PLACE` in `MPI_Gather`). For example +On the root process, `sendbuf` can be [`MPI.IN_PLACE`](@ref) on the root process, in which +case the corresponding entries in `recvbuf` are assumed to be already in place (this +corresponds the behaviour of `MPI_IN_PLACE` in `MPI_Gather`). For example: ``` if root == MPI.Comm_rank(comm) - Gather!(nothing, buf, count, root, comm) + MPI.Gather!(MPI.IN_PLACE, UBuffer(buf, count), root, comm) else - Gather!(buf, nothing, count, root, comm) + MPI.Gather!(buf, nothing, root, comm) end ``` -`recvbuf` on the root process should be a buffer of length `count*Comm_size(comm)`, and -on non-root processes it is ignored and can be `nothing`. - -`count` should be the same for all processes. +`recvbuf` on the root process should be a [`UBuffer`](@ref), or can be an `AbstractArray` +if the length can be determined from `sendbuf`. On non-root processes it is ignored and +can be `nothing`. # See also - [`Gather`](@ref) for the allocating operation. @@ -227,36 +196,38 @@ on non-root processes it is ignored and can be `nothing`. # External links $(_doc_external("MPI_Gather")) """ -function Gather!(sendbuf, recvbuf, count::Integer, root::Integer, comm::Comm) - isroot = Comm_rank(comm) == root - if isroot - @assert recvbuf !== nothing - @assert_minlength recvbuf count*Comm_size(comm) - if sendbuf === nothing - sendbuf = MPI_IN_PLACE - end +function Gather!(sendbuf::Buffer, recvbuf::UBuffer, root::Integer, comm::Comm) + if recvbuf.nchunks !== nothing && Comm_rank(comm) == root + @assert recvbuf.nchunks >= Comm_size(comm) end - @assert_minlength sendbuf count - T = sendbuf isa SentinelPtr ? eltype(recvbuf) : eltype(sendbuf) # int MPI_Gather(const void* sendbuf, int sendcount, MPI_Datatype sendtype, # void* recvbuf, int recvcount, MPI_Datatype recvtype, int root, # MPI_Comm comm) @mpichk ccall((:MPI_Gather, libmpi), Cint, (MPIPtr, Cint, MPI_Datatype, MPIPtr, Cint, MPI_Datatype, Cint, MPI_Comm), - sendbuf, count, Datatype(T), recvbuf, count, Datatype(T), root, comm) - isroot ? recvbuf : nothing -end -function Gather!(sendbuf, recvbuf, root::Integer, comm::Comm) - Gather!(sendbuf, recvbuf, length(sendbuf), root, comm) + sendbuf.data, sendbuf.count, sendbuf.datatype, + recvbuf.data, recvbuf.count, recvbuf.datatype, root, comm) + return recvbuf.data end +Gather!(sendbuf, recvbuf::UBuffer, root::Integer, comm::Comm) = + Gather!(Buffer_send(sendbuf), recvbuf, root, comm) +Gather!(sendbuf::Union{Ref,AbstractArray}, recvbuf::AbstractArray, root::Integer, comm::Comm) = + Gather!(sendbuf, UBuffer(recvbuf, length(sendbuf)), root, comm) +Gather!(sendbuf, recvbuf::Nothing, root::Integer, comm::Comm) = + Gather!(sendbuf, UBuffer(nothing), root, comm) +Gather!(sendbuf::Nothing, recvbuf, root::Integer, comm::Comm) = + Gather!(IN_PLACE, recvbuf, root, comm) + + """ - Gather(sendbuf[, count=length(sendbuf)], root, comm) + Gather(sendbuf, root, comm::Comm) -Each process sends the first `count` elements of the buffer `sendbuf` to the `root` -process. The `root` allocates the output buffer and stores elements in rank order. +Each process sends the contents of the buffer `sendbuf` to the `root` process. The `root` +allocates the output buffer and stores elements in rank order. -`count` should be the same for all processes. +`sendbuf` can be an `AbstractArray` or a scalar, and should be the same length on all +processes. # See also - [`Gather!`](@ref) for the mutating operation. @@ -266,27 +237,77 @@ process. The `root` allocates the output buffer and stores elements in rank orde # External links $(_doc_external("MPI_Gather")) """ -function Gather(sendbuf, count::Integer, root::Integer, comm::Comm) - Gather!(sendbuf, Comm_rank(comm) == root ? similar(sendbuf, Comm_size(comm) * count) : nothing, count, root, comm) -end -function Gather(sendbuf::AbstractArray, root::Integer, comm::Comm) - Gather(sendbuf, length(sendbuf), root, comm) +Gather(sendbuf::AbstractArray, root::Integer, comm::Comm) = + Gather!(sendbuf, Comm_rank(comm) == root ? similar(sendbuf, Comm_size(comm) * length(sendbuf)) : nothing, root, comm) +Gather(object::T, root::Integer, comm::Comm) where {T} = + Gather!(Ref(object), Comm_rank(comm) == root ? Array{T}(undef, Comm_size(comm)) : nothing, root, comm) + +""" + Gatherv!(sendbuf, recvbuf::Union{VBuffer,Nothing}, root, comm) + +Each process sends the contents of the buffer `sendbuf` to the `root` process. The `root` +stores elements in rank order in the buffer `recvbuf`. + +`sendbuf` should be a [`Buffer`](@ref) object, or any object for which +[`Buffer_send`](@ref) is defined, with the same length on all processes. + +On the root process, `sendbuf` can be [`MPI.IN_PLACE`](@ref), in which case the +corresponding entries in `recvbuf` are assumed to be already in place. For example +``` +if root == MPI.Comm_rank(comm) + Gatherv!(MPI.IN_PLACE, VBuffer(buf, counts), root, comm) +else + Gatherv!(buf, nothing, root, comm) end -function Gather(object::T, root::Integer, comm::Comm) where {T} - Gather!(Ref(object), Comm_rank(comm) == root ? Vector{T}(undef, Comm_size(comm)) : nothing, 1, root, comm) +``` + +`recvbuf` on the root process should be a [`VBuffer`](@ref), or can be an `AbstractArray` +if the length can be determined from `sendbuf`. On non-root processes it is ignored and +can be `nothing`. + +# See also +- [`Gatherv`](@ref) for the allocating operation +- [`Gather!`](@ref) +- [`Allgatherv!`](@ref)/[`Allgatherv`](@ref) to send the result to all processes + +# External links +$(_doc_external("MPI_Gatherv")) +""" +function Gatherv!(sendbuf::Buffer, recvbuf::VBuffer, root::Integer, comm::Comm) + if Comm_rank(comm) == root + @assert length(recvbuf.counts) >= Comm_size(comm) + end + # int MPI_Gatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, + # void* recvbuf, const int recvcounts[], const int displs[], + # MPI_Datatype recvtype, int root, MPI_Comm comm) + @mpichk ccall((:MPI_Gatherv, libmpi), Cint, + (MPIPtr, Cint, MPI_Datatype, MPIPtr, Ptr{Cint}, Ptr{Cint}, MPI_Datatype, Cint, MPI_Comm), + sendbuf.data, sendbuf.count, sendbuf.datatype, + recvbuf.data, recvbuf.counts, recvbuf.displs, recvbuf.datatype, root, comm) + return recvbuf.data end +Gatherv!(sendbuf, recvbuf::VBuffer, root::Integer, comm::Comm) = + Gatherv!(Buffer_send(sendbuf), recvbuf, root, comm) +Gatherv!(sendbuf, recvbuf::Nothing, root::Integer, comm::Comm) = + Gatherv!(sendbuf, VBuffer(nothing), root, comm) + + """ - Allgather!(sendbuf, recvbuf[, count::Integer=length(sendbuf)], comm::Comm) - Allgather!(sendrecvbuf, count::Integer, comm::Comm) + Allgather!(sendbuf, recvbuf::UBuffer, comm::Comm) + Allgather!(sendrecvbuf::UBuffer, comm::Comm) + +Each process sends the contents of `sendbuf` to the other processes, the result of which +is stored in rank order into `recvbuf`. -Each process sends the first `count` elements of `sendbuf` to the other processes, who -store the results in rank order into `recvbuf`. +`sendbuf` can be a [`Buffer`](@ref) object, or any object for which [`Buffer_send`](@ref) +is defined, and should be the same length on all processes. -`count` should be the same for all processes. +`recvbuf` can be a [`UBuffer`](@ref), or can be an `AbstractArray` if the length can be +determined from `sendbuf`. -If only one buffer `sendrecvbuf` is provided, then each process send data is assumed to be -in the area where it would receive it's own contribution. +If only one buffer `sendrecvbuf` is provided, then on each process the data to send is +assumed to be in the area where it would receive its own contribution. # See also - [`Allgather`](@ref) for the allocating operation @@ -296,31 +317,39 @@ in the area where it would receive it's own contribution. # External links $(_doc_external("MPI_Allgather")) """ -function Allgather!(sendbuf, recvbuf, count::Integer, comm::Comm) - @assert recvbuf !== nothing - @assert_minlength recvbuf count*Comm_size(comm) - @assert_minlength sendbuf count - T = eltype(recvbuf) +function Allgather!(sendbuf::Buffer, recvbuf::UBuffer, comm::Comm) + if recvbuf.nchunks !== nothing + @assert recvbuf.nchunks >= Comm_size(comm) + end # int MPI_Allgather(const void* sendbuf, int sendcount, # MPI_Datatype sendtype, void* recvbuf, int recvcount, # MPI_Datatype recvtype, MPI_Comm comm) @mpichk ccall((:MPI_Allgather, libmpi), Cint, (MPIPtr, Cint, MPI_Datatype, MPIPtr, Cint, MPI_Datatype, MPI_Comm), - sendbuf, count, Datatype(T), recvbuf, count, Datatype(T), comm) - recvbuf + sendbuf.data, sendbuf.count, sendbuf.datatype, + recvbuf.data, recvbuf.count, recvbuf.datatype, comm) + return recvbuf.data end -function Allgather!(sendrecvbuf, count::Integer, comm::Comm) - Allgather!(MPI.IN_PLACE, sendrecvbuf, count, comm) +Allgather!(sendbuf, recvbuf::UBuffer, comm::Comm) = + Allgather!(Buffer_send(sendbuf), recvbuf, comm) + +Allgather!(sendbuf::Union{Ref,AbstractArray}, recvbuf::AbstractArray, comm::Comm) = + Allgather!(sendbuf, UBuffer(recvbuf, length(sendbuf)), comm) + + +function Allgather!(sendrecvbuf::UBuffer, comm::Comm) + Allgather!(IN_PLACE, sendrecvbuf, comm) end """ - Allgather(sendbuf[, count=length(sendbuf)], comm) + Allgather(sendbuf, comm) -Each process sends the first `count` elements of `sendbuf` to the other processes, who -store the results in rank order allocating the output buffer. +Each process sends the contents of `sendbuf` to the other processes, who store the results +in rank order allocating the output buffer. -`count` should be the same for all processes. +`sendbuf` can be an `AbstractArray` or a scalar, and should be the same size on all +processes. # See also - [`Allgather!`](@ref) for the mutating operation @@ -330,95 +359,26 @@ store the results in rank order allocating the output buffer. # External links $(_doc_external("MPI_Allgather")) """ -function Allgather(sendbuf, count::Integer, comm::Comm) - Allgather!(sendbuf, similar(sendbuf, Comm_size(comm) * count), count, comm) -end function Allgather(sendbuf::AbstractArray, comm::Comm) - Allgather(sendbuf, length(sendbuf), comm) + Allgather!(sendbuf, similar(sendbuf, Comm_size(comm) * length(sendbuf)), comm) end function Allgather(object::T, comm::Comm) where T - Allgather!(Ref(object), Vector{T}(undef, Comm_size(comm)), 1, comm) + Allgather!(Ref(object), Array{T}(undef, Comm_size(comm)), comm) end -""" - Gatherv!(sendbuf, recvbuf, counts, root, comm) - -Each process sends the first `counts[rank]` elements of the buffer `sendbuf` to -the `root` process. The `root` stores elements in rank order in the buffer -`recvbuf`. - -`sendbuf` can be `nothing` on the root process, in which case the corresponding entries in -`recvbuf` are assumed to be already in place (this corresponds the behaviour of -`MPI_IN_PLACE` in `MPI_Gatherv`). For example -``` -if root == MPI.Comm_rank(comm) - Gatherv!(nothing, buf, counts, root, comm) -else - Gatherv!(buf, nothing, counts, root, comm) -end -``` - -# See also -- [`Gatherv`](@ref) for the allocating operation -- [`Gather!`](@ref) -- [`Allgatherv!`](@ref)/[`Allgatherv`](@ref) to send the result to all processes - -# External links -$(_doc_external("MPI_Gatherv")) -""" -function Gatherv!(sendbuf, recvbuf, counts::Vector{Cint}, root::Integer, comm::Comm) - isroot = Comm_rank(comm) == root - displs = accumulate(+, counts) - counts - sendcnt = counts[Comm_rank(comm) + 1] - if isroot - @assert recvbuf !== nothing - @assert_minlength recvbuf sum(counts) - if sendbuf === nothing - sendbuf = MPI_IN_PLACE - end - end - @assert_minlength sendbuf sendcnt - T = sendbuf isa SentinelPtr ? eltype(recvbuf) : eltype(sendbuf) - - # int MPI_Gatherv(const void* sendbuf, int sendcount, MPI_Datatype sendtype, - # void* recvbuf, const int recvcounts[], const int displs[], - # MPI_Datatype recvtype, int root, MPI_Comm comm) - @mpichk ccall((:MPI_Gatherv, libmpi), Cint, - (MPIPtr, Cint, MPI_Datatype, MPIPtr, Ptr{Cint}, Ptr{Cint}, MPI_Datatype, Cint, MPI_Comm), - sendbuf, sendcnt, Datatype(T), recvbuf, counts, displs, Datatype(T), root, comm) - isroot ? recvbuf : nothing -end """ - Gatherv(sendbuf, counts, root, comm) + Allgatherv!(sendbuf, recvbuf::VBuffer, comm::Comm) + Allgatherv!(sendrecvbuf::VBuffer, comm::Comm) -Each process sends the first `counts[rank]` elements of the buffer `sendbuf` to -the `root` process. The `root` allocates the output buffer and stores elements -in rank order. +Each process sends the contents of `sendbuf` to all other process. Each process stores the +received in the [`VBuffer`](@ref) `recvbuf`. -# See also -- [`Gatherv!`](@ref) for the mutating operation -- [`Gather!`](@ref)/[`Gather`](@ref) -- [`Allgatherv!`](@ref)/[`Allgatherv`](@ref) to send the result to all processes - -# External links -$(_doc_external("MPI_Gatherv")) -""" -function Gatherv(sendbuf, counts::Vector{Cint}, root::Integer, comm::Comm) - Gatherv!(sendbuf, Comm_rank(comm) == root ? similar(sendbuf, sum(counts)) : nothing, counts, root, comm) -end - - -""" - Allgatherv!(sendbuf, recvbuf, counts, comm) - Allgatherv!(sendrecvbuf, counts, comm) - -Each process sends the first `counts[rank]` elements of the buffer `sendbuf` to -all other process. Each process stores the received data in rank order in the -buffer `recvbuf`. +`sendbuf` can be a [`Buffer`](@ref) object, or any object for which [`Buffer_send`](@ref) +is defined. If only one buffer `sendrecvbuf` is provided, then for each process, the data to be sent -is taken from the interval of `recvbuf` where it would store it's own data. +is taken from the interval of `recvbuf` where it would store its own data. # See also - [`Allgatherv`](@ref) for the allocating operation @@ -427,54 +387,35 @@ is taken from the interval of `recvbuf` where it would store it's own data. # External links $(_doc_external("MPI_Allgatherv")) """ -function Allgatherv!(sendbuf, recvbuf, counts::Vector{Cint}, comm::Comm) - displs = accumulate(+, counts) - counts - sendcnt = counts[Comm_rank(comm) + 1] - @assert recvbuf !== nothing - @assert_minlength recvbuf sum(counts) - @assert_minlength sendbuf sendcnt - T = eltype(recvbuf) +function Allgatherv!(sendbuf::Buffer, recvbuf::VBuffer, comm::Comm) + @assert length(recvbuf.counts) >= Comm_size(comm) # int MPI_Allgatherv(const void* sendbuf, int sendcount, # MPI_Datatype sendtype, void* recvbuf, const int recvcounts[], # const int displs[], MPI_Datatype recvtype, MPI_Comm comm) @mpichk ccall((:MPI_Allgatherv, libmpi), Cint, - (MPIPtr, Cint, MPI_Datatype, MPIPtr, Ptr{Cint}, Ptr{Cint}, MPI_Datatype, MPI_Comm), - sendbuf, sendcnt, Datatype(T), recvbuf, counts, displs, Datatype(T), comm) - recvbuf + (MPIPtr, Cint, MPI_Datatype, MPIPtr, Ptr{Cint}, Ptr{Cint}, MPI_Datatype, MPI_Comm), + sendbuf.data, sendbuf.count, sendbuf.datatype, + recvbuf.data, recvbuf.counts, recvbuf.displs, recvbuf.datatype, + comm) + return recvbuf.data end -function Allgatherv!(sendrecvbuf, counts::Vector{Cint}, comm::Comm) - Allgatherv!(MPI_IN_PLACE, sendrecvbuf, counts, comm) +function Allgatherv!(sendbuf, recvbuf::VBuffer, comm::Comm) + Allgatherv!(Buffer_send(sendbuf), recvbuf, comm) end - -""" - Allgatherv(sendbuf, counts, comm) - -Each process sends the first `counts[rank]` elements of the buffer `sendbuf` to -all other process. Each process allocates an output buffer and stores the -received data in rank order. - -# See also -- [`Allgatherv!`](@ref) for the mutating operation. -- [`Gatherv!`](@ref)/[`Gatherv`](@ref) to send the result to a single process. - -# External links -$(_doc_external("MPI_Allgatherv")) -""" -function Allgatherv(sendbuf, counts::Vector{Cint}, comm::Comm) - recvbuf = similar(sendbuf, sum(counts)) - Allgatherv!(sendbuf, recvbuf, counts, comm) +function Allgatherv!(sendrecvbuf::VBuffer, comm::Comm) + Allgatherv!(IN_PLACE, sendrecvbuf, comm) end + """ - Alltoall!(sendbuf, recvbuf, count::Integer, comm::Comm) - Alltoall!(sendrecvbuf, count::Integer, comm::Comm) + Alltoall!(sendbuf::UBuffer, recvbuf::UBuffer, comm::Comm) + Alltoall!(sendrecvbuf::UBuffer, comm::Comm) -Every process divides the buffer `sendbuf` into `Comm_size(comm)` chunks of -length `count`, sending the `j`-th chunk to the `j`-th process. -Every process stores the data received from the `j`-th process in the `j`-th -chunk of the buffer `recvbuf`. +Every process divides the [`UBuffer`](@ref) `sendbuf` into `Comm_size(comm)` chunks of +equal size, sending the `j`-th chunk to the process of rank `j-1`. Every process stores +the data received from rank `j-1` process in the `j`-th chunk of the buffer `recvbuf`. ``` rank send buf recv buf @@ -484,7 +425,7 @@ rank send buf recv buf 2 α,β,γ,ψ,η,ν e,f,E,F,η,ν ``` -If only one buffer `sendrecvbuf` then data is overwritten. +If only one buffer `sendrecvbuf` is used, then data is overwritten. # See also - [`Alltoall`](@ref) for the allocating operation @@ -492,31 +433,35 @@ If only one buffer `sendrecvbuf` then data is overwritten. # External links $(_doc_external("MPI_Alltoall")) """ -function Alltoall!(sendbuf, recvbuf, count::Integer, comm::Comm) - buflength = count * Comm_size(comm) - @assert_minlength recvbuf buflength - @assert_minlength sendbuf buflength - sendbuf isa SentinelPtr || @assert eltype(sendbuf) == eltype(recvbuf) - T = eltype(recvbuf) +function Alltoall!(sendbuf::UBuffer, recvbuf::UBuffer, comm::Comm) + if sendbuf.data !== MPI_IN_PLACE && sendbuf.nchunks !== nothing + @assert sendbuf.nchunks >= Comm_size(comm) + end + if recvbuf.nchunks !== nothing + @assert recvbuf.nchunks >= Comm_size(comm) + end # int MPI_Alltoall(const void* sendbuf, int sendcount, MPI_Datatype sendtype, # void* recvbuf, int recvcount, MPI_Datatype recvtype, # MPI_Comm comm) @mpichk ccall((:MPI_Alltoall, libmpi), Cint, (MPIPtr, Cint, MPI_Datatype, MPIPtr, Cint, MPI_Datatype, MPI_Comm), - sendbuf, count, Datatype(T), recvbuf, count, Datatype(T), comm) - recvbuf -end -function Alltoall!(sendrecvbuf, count::Integer, comm::Comm) - Alltoall!(MPI_IN_PLACE, recvbuf, count, comm) + sendbuf.data, sendbuf.count, sendbuf.datatype, + recvbuf.data, recvbuf.count, recvbuf.datatype, + comm) + return recvbuf.data end +Alltoall!(sendbuf::InPlace, recvbuf::UBuffer, comm::Comm) = + Alltoall!(UBuffer(IN_PLACE), recvbuf, comm) +Alltoall!(sendrecvbuf::UBuffer, comm::Comm) = + Alltoall!(IN_PLACE, recvbuf, comm) """ - Alltoall(sendbuf, count::Integer, comm::Comm) + Alltoall(sendbuf::UBuffer, comm::Comm) -Every process divides the buffer `sendbuf` into `Comm_size(comm)` chunks of -length `count`, sending the `j`-th chunk to the `j`-th process. -Every process allocates the output buffer and stores the data received from the -`j`-th process in the `j`-th chunk. +Every process divides the [`UBuffer`](@ref) `sendbuf` into `Comm_size(comm)` chunks of +equal size, sending the `j`-th chunk to the process of rank `j-1`. Every process allocates +the output buffer and stores the data received from the process on rank `j-1` in the +`j`-th chunk. ``` rank send buf recv buf @@ -532,55 +477,39 @@ rank send buf recv buf # External links $(_doc_external("MPI_Alltoall")) """ -function Alltoall(sendbuf, count::Integer, comm::Comm) - recvbuf = similar(sendbuf, Comm_size(comm)*count) - Alltoall!(sendbuf, recvbuf, count, comm) -end +Alltoall(sendbuf::UBuffer, comm::Comm) = + Alltoall!(sendbuf, similar(sendbuf), comm) + """ - Alltoallv!(sendbuf, recvbuf, scounts::Vector, rcounts::Vector, comm::Comm) + Alltoallv!(sendbuf::VBuffer, recvbuf::VBuffer, comm::Comm) -Similar to [`Alltoall!`](@ref), except with different size chunks per process. +Similar to [`Alltoall!`](@ref), except with different size chunks per process. # See also -- [`Alltoallv`](@ref) for the allocating operation +- [`VBuffer`](@ref) # External links $(_doc_external("MPI_Alltoallv")) """ -function Alltoallv!(sendbuf, recvbuf, scounts::Vector{Cint}, rcounts::Vector{Cint}, comm::Comm) - @assert_minlength sendbuf sum(scounts) - @assert_minlength recvbuf sum(rcounts) - @assert eltype(sendbuf) == eltype(recvbuf) - T = eltype(sendbuf) - - sdispls = accumulate(+, scounts) - scounts - rdispls = accumulate(+, rcounts) - rcounts +function Alltoallv!(sendbuf::VBuffer, recvbuf::VBuffer, comm::Comm) + if sendbuf.data !== MPI_IN_PLACE + @assert length(sendbuf.counts) >= Comm_size(comm) + end + @assert length(recvbuf.counts) >= Comm_size(comm) # int MPI_Alltoallv(const void* sendbuf, const int sendcounts[], # const int sdispls[], MPI_Datatype sendtype, void* recvbuf, # const int recvcounts[], const int rdispls[], # MPI_Datatype recvtype, MPI_Comm comm) @mpichk ccall((:MPI_Alltoallv, libmpi), Cint, - (MPIPtr, Ptr{Cint}, Ptr{Cint}, MPI_Datatype, MPIPtr, Ptr{Cint}, Ptr{Cint}, MPI_Datatype, MPI_Comm), - sendbuf, scounts, sdispls, Datatype(T), recvbuf, rcounts, rdispls, Datatype(T), comm) - recvbuf -end - -""" - Alltoallv(sendbuf, recvbuf, scounts::Vector, rcounts::Vector, comm::Comm) - -Similar to [`Alltoall`](@ref), except with different size chunks per process. - -# See also -- [`Alltoallv!`](@ref) for the mutating operation - -# External links -$(_doc_external("MPI_Alltoallv")) -""" -function Alltoallv(sendbuf, scounts::Vector{Cint}, - rcounts::Vector{Cint}, comm::Comm) - recvbuf = similar(sendbuf, sum(rcounts)) - Alltoallv!(sendbuf, recvbuf, scounts, rcounts, comm) + (MPIPtr, Ptr{Cint}, Ptr{Cint}, MPI_Datatype, + MPIPtr, Ptr{Cint}, Ptr{Cint}, MPI_Datatype, + MPI_Comm), + sendbuf.data, sendbuf.counts, sendbuf.displs, sendbuf.datatype, + recvbuf.data, recvbuf.counts, recvbuf.displs, recvbuf.datatype, + comm) + + return recvbuf.data end @@ -590,11 +519,11 @@ end # mutating """ - Reduce!(sendbuf, recvbuf[, count::Integer=length(sendbuf)], op, root::Integer, comm::Comm) + Reduce!(sendbuf, recvbuf, op, root::Integer, comm::Comm) Reduce!(sendrecvbuf, op, root::Integer, comm::Comm) -Performs elementwise reduction using the operator `op` on the first `count` elements of -the buffer `sendbuf` and stores the result in `recvbuf` on the process of rank `root`. +Performs elementwise reduction using the operator `op` on the buffer `sendbuf` and stores +the result in `recvbuf` on the process of rank `root`. On non-root processes `recvbuf` is ignored, and can be `nothing`. @@ -608,38 +537,27 @@ To perform the reduction in place, provide a single buffer `sendrecvbuf`. # External links $(_doc_external("MPI_Reduce")) """ -function Reduce!(sendbuf, recvbuf, count::Integer, op::Union{Op,MPI_Op}, root::Integer, comm::Comm) - isroot = Comm_rank(comm) == root - @assert_minlength sendbuf count - if isroot - @assert recvbuf !== nothing - @assert_minlength recvbuf count - end - T = sendbuf isa SentinelPtr ? eltype(recvbuf) : eltype(sendbuf) +function Reduce!(rbuf::RBuffer, op::Union{Op,MPI_Op}, root::Integer, comm::Comm) # int MPI_Reduce(const void* sendbuf, void* recvbuf, int count, # MPI_Datatype datatype, MPI_Op op, int root, MPI_Comm comm) @mpichk ccall((:MPI_Reduce, libmpi), Cint, (MPIPtr, MPIPtr, Cint, MPI_Datatype, MPI_Op, Cint, MPI_Comm), - sendbuf, recvbuf, count, Datatype(T), op, root, comm) - recvbuf + rbuf.senddata, rbuf.recvdata, rbuf.count, rbuf.datatype, op, root, comm) + return rbuf.recvdata end # Convert user-provided functions to MPI.Op -function Reduce!(sendbuf, recvbuf, count::Integer, opfunc, root::Integer, comm::Comm) - T = sendbuf isa SentinelPtr ? eltype(recvbuf) : eltype(sendbuf) - Reduce!(sendbuf, recvbuf, count, Op(opfunc, T), root, comm) -end -function Reduce!(sendbuf::AbstractArray, recvbuf::Union{Nothing,AbstractArray}, - op, root::Integer, comm::Comm) - Reduce!(sendbuf, recvbuf, length(sendbuf), op, root, comm) -end +Reduce!(rbuf::RBuffer, op, root::Integer, comm::Comm) = + Reduce!(rbuf, Op(op, eltype(rbuf)), root, comm) +Reduce!(sendbuf, recvbuf, op, root::Integer, comm::Comm) = + Reduce!(RBuffer(sendbuf, recvbuf), op, root, comm) # inplace function Reduce!(buf, op, root::Integer, comm::Comm) if Comm_rank(comm) == root - Reduce!(MPI_IN_PLACE, buf, length(buf), op, root, comm) + Reduce!(IN_PLACE, buf, op, root, comm) else - Reduce!(buf, nothing, length(buf), op, root, comm) + Reduce!(buf, nothing, op, root, comm) end end @@ -661,13 +579,17 @@ the result `recvbuf` on the process of rank `root`, and `nothing` on non-root pr $(_doc_external("MPI_Reduce")) """ function Reduce(sendbuf::AbstractArray, op, root::Integer, comm::Comm) - Reduce!(sendbuf, Comm_rank(comm) == root ? similar(sendbuf) : nothing, length(sendbuf), op, root, comm) + if Comm_rank(comm) == root + Reduce!(sendbuf, similar(sendbuf), op, root, comm) + else + Reduce!(sendbuf, nothing, op, root, comm) + end end function Reduce(object::T, op, root::Integer, comm::Comm) where {T} if Comm_rank(comm) == root - Reduce!(Ref(object), Ref{T}(), 1, op, root, comm)[] + Reduce!(Ref(object), Ref{T}(), op, root, comm)[] else - Reduce!(Ref(object), nothing, 1, op, root, comm) + Reduce!(Ref(object), nothing, op, root, comm) end end @@ -675,11 +597,11 @@ end # mutating """ - Allreduce!(sendbuf, recvbuf[, count=length(sendbuf)], op, comm) - Allreduce!(sendrecvbuf, op, comm) + Allreduce!(sendbuf, recvbuf, op, comm::Comm) + Allreduce!(sendrecvbuf, op, comm::Comm) -Performs elementwise reduction using the operator `op` on the first `count` elements of -the buffer `sendbuf`, storing the result in the `recvbuf` of all processes in the group. +Performs elementwise reduction using the operator `op` on the buffer `sendbuf`, storing +the result in the `recvbuf` of all processes in the group. `Allreduce!` is equivalent to a [`Reduce!`](@ref) operation followed by a [`Bcast!`](@ref), but can lead to better performance. @@ -694,30 +616,22 @@ If only one `sendrecvbuf` buffer is provided, then the operation is performed in # External links $(_doc_external("MPI_Allreduce")) """ -function Allreduce!(sendbuf, recvbuf, count::Integer, op::Union{Op,MPI_Op}, comm::Comm) - @assert_minlength sendbuf count - @assert_minlength recvbuf count - sendbuf isa SentinelPtr || @assert eltype(sendbuf) == eltype(recvbuf) - T = eltype(recvbuf) +function Allreduce!(rbuf::RBuffer, op::Union{Op,MPI_Op}, comm::Comm) # int MPI_Allreduce(const void* sendbuf, void* recvbuf, int count, # MPI_Datatype datatype, MPI_Op op, MPI_Comm comm) @mpichk ccall((:MPI_Allreduce, libmpi), Cint, (MPIPtr, MPIPtr, Cint, MPI_Datatype, MPI_Op, MPI_Comm), - sendbuf, recvbuf, count, Datatype(T), op, comm) - recvbuf -end -function Allreduce!(sendbuf, recvbuf, count::Integer, opfunc, comm::Comm) - Allreduce!(sendbuf, recvbuf, count, Op(opfunc, eltype(recvbuf)), comm) -end - -function Allreduce!(sendbuf::AbstractArray, recvbuf::AbstractArray, op, comm::Comm) - Allreduce!(sendbuf, recvbuf, length(recvbuf), op, comm) + rbuf.senddata, rbuf.recvdata, rbuf.count, rbuf.datatype, op, comm) + rbuf.recvdata end +Allreduce!(rbuf::RBuffer, op, comm::Comm) = + Allreduce!(rbuf, Op(op, eltype(rbuf)), comm) +Allreduce!(sendbuf, recvbuf, op, comm::Comm) = + Allreduce!(RBuffer(sendbuf, recvbuf), op, comm) # inplace -function Allreduce!(buf, op, comm::Comm) - Allreduce!(MPI_IN_PLACE, buf, length(buf), op, comm) -end +Allreduce!(buf, op, comm::Comm) = + Allreduce!(IN_PLACE, buf, op, comm) # allocating """ @@ -736,24 +650,23 @@ the result in the `recvbuf` of all processes in the group. # External links $(_doc_external("MPI_Allreduce")) """ -function Allreduce(sendbuf::AbstractArray, op, comm::Comm) - Allreduce!(sendbuf, similar(sendbuf), length(sendbuf), op, comm) -end -function Allreduce(obj::T, op, comm::Comm) where {T} - Allreduce!(Ref(obj), Ref{T}(), 1, op, comm)[] -end +Allreduce(sendbuf::AbstractArray, op, comm::Comm) = + Allreduce!(sendbuf, similar(sendbuf), op, comm) +Allreduce(obj::T, op, comm::Comm) where {T} = + Allreduce!(Ref(obj), Ref{T}(), op, comm)[] ## Scan # mutating """ - Scan!(sendbuf, recvbuf[, count::Integer], op, comm::Comm) - Scan!(buf[, count::Integer], op, comm::Comm) + Scan!(sendbuf, recvbuf, op, comm::Comm) + Scan!(sendrecvbuf, op, comm::Comm) Inclusive prefix reduction (analagous to `accumulate` in Julia): `recvbuf` on rank `i` will contain the the result of reducing `sendbuf` by `op` from ranks `0:i`. -If only a single buffer is provided, then operations will be performed in-place in `buf`. +If only a single buffer `sendrecvbuf` is provided, then operations will be performed +in-place. # See also - [`Scan`](@ref) to handle allocation of the output buffer @@ -763,30 +676,22 @@ If only a single buffer is provided, then operations will be performed in-place # External links $(_doc_external("MPI_Scan")) """ -function Scan!(sendbuf, recvbuf, count::Integer, - op::Union{Op,MPI_Op}, comm::Comm) - T = eltype(recvbuf) +function Scan!(rbuf::RBuffer, op::Union{Op,MPI_Op}, comm::Comm) # int MPI_Scan(const void* sendbuf, void* recvbuf, int count, # MPI_Datatype datatype, MPI_Op op, MPI_Comm comm) @mpichk ccall((:MPI_Scan, libmpi), Cint, (MPIPtr, MPIPtr, Cint, MPI_Datatype, MPI_Op, MPI_Comm), - sendbuf, recvbuf, count, Datatype(T), op, comm) - recvbuf -end -function Scan!(sendbuf, recvbuf, count::Integer, opfunc, comm::Comm) - Scan!(sendbuf, recvbuf, count, Op(opfunc, eltype(recvbuf)), comm) -end -function Scan!(sendbuf::AbstractArray, recvbuf, op, comm::Comm) - Scan!(sendbuf, recvbuf, length(sendbuf), op, comm) + rbuf.senddata, rbuf.recvdata, rbuf.count, rbuf.datatype, op, comm) + rbuf.recvdata end +Scan!(rbuf::RBuffer, op, comm::Comm) = + Scan!(rbuf, Op(op, eltype(rbuf)), comm) +Scan!(sendbuf, recvbuf, op, comm::Comm) = + Scan!(RBuffer(sendbuf, recvbuf), op, comm) # inplace -function Scan!(buf, count::Integer, opfunc, comm::Comm) - Scan!(MPI_IN_PLACE, buf, count, Op(opfunc, eltype(sendbuf)), comm) -end -function Scan!(buf, opfunc, comm::Comm) - Scan!(buf, length(buf), Op(opfunc, eltype(sendbuf)), comm) -end +Scan!(buf, op, comm::Comm) = + Scan!(IN_PLACE, buf, op, comm) # allocating """ @@ -806,28 +711,24 @@ type. # External links $(_doc_external("MPI_Scan")) """ -function Scan(sendbuf::AbstractArray, op, comm::Comm) +Scan(sendbuf::AbstractArray, op, comm::Comm) = Scan!(sendbuf, similar(sendbuf), op, comm) -end -function Scan(object::T, op, comm::Comm) where {T} - Scan!(Ref(object), Ref{T}(), 1, op, comm)[] -end +Scan(object::T, op, comm::Comm) where {T} = + Scan!(Ref(object), Ref{T}(), op, comm)[] ## Exscan - # mutating - """ - Exscan!(sendbuf, recvbuf[, count::Integer], op, comm::Comm) - Exscan!(buf[, count::Integer], op, comm::Comm) + Exscan!(sendbuf, recvbuf, op, comm::Comm) + Exscan!(sendrecvbuf, op, comm::Comm) Exclusive prefix reduction (analagous to `accumulate` in Julia): `recvbuf` on rank `i` will contain the the result of reducing `sendbuf` by `op` from ranks `0:i-1`. The `recvbuf` on rank `0` is ignored, and the `recvbuf` on rank `1` will contain the contents of `sendbuf` on rank `0`. -If only a single `buf` is provided, then operations are performed in-place, and `buf` on -rank 0 will remain unchanged. +If only a single `sendrecvbuf` is provided, then operations are performed in-place, and +`buf` on rank 0 will remain unchanged. # See also - [`Exscan`](@ref) to handle allocation of the output buffer @@ -837,30 +738,22 @@ rank 0 will remain unchanged. # External links $(_doc_external("MPI_Exscan")) """ -function Exscan!(sendbuf, recvbuf, count::Integer, - op::Union{Op,MPI_Op}, comm::Comm) - T = eltype(recvbuf) +function Exscan!(rbuf::RBuffer, op::Union{Op,MPI_Op}, comm::Comm) # int MPI_Exscan(const void* sendbuf, void* recvbuf, int count, # MPI_Datatype datatype, MPI_Op op, MPI_Comm comm) @mpichk ccall((:MPI_Exscan, libmpi), Cint, (MPIPtr, MPIPtr, Cint, MPI_Datatype, MPI_Op, MPI_Comm), - sendbuf, recvbuf, count, Datatype(T), op, comm) - recvbuf -end -function Exscan!(sendbuf, recvbuf, count::Integer, opfunc, comm::Comm) - Exscan!(sendbuf, recvbuf, count, Op(opfunc, eltype(recvbuf)), comm) -end -function Exscan!(sendbuf::AbstractArray, recvbuf::AbstractArray, op, comm::Comm) - Exscan!(sendbuf, recvbuf, length(sendbuf), op, comm) + rbuf.senddata, rbuf.recvdata, rbuf.count, rbuf.datatype, op, comm) + rbuf.recvdata end +Exscan!(rbuf::RBuffer, op, comm::Comm) = + Exscan!(rbuf, Op(op, eltype(rbuf)), comm) +Exscan!(sendbuf, recvbuf, op, comm::Comm) = + Exscan!(RBuffer(sendbuf, recvbuf), op, comm) # inplace -function Exscan!(buf, count::Integer, opfunc, comm::Comm) - Exscan!(MPI_IN_PLACE, buf, count, Op(opfunc, eltype(buf)), comm) -end -function Exscan!(buf, op, comm::Comm) - Exscan!(buf, length(buf), op, comm) -end +Exscan!(buf, op, comm::Comm) = + Exscan!(IN_PLACE, buf, op, comm) # allocating @@ -880,9 +773,7 @@ of `sendbuf` on rank `0`. # External links $(_doc_external("MPI_Scan")) """ -function Exscan(sendbuf::AbstractArray, op, comm::Comm) +Exscan(sendbuf::AbstractArray, op, comm::Comm) = Exscan!(sendbuf, similar(sendbuf), op, comm) -end -function Exscan(object::T, op, comm::Comm) where {T} - Exscan!(Ref(object), Ref{T}(), 1, op, comm)[] -end +Exscan(object::T, op, comm::Comm) where {T} = + Exscan!(Ref(object), Ref{T}(), op, comm)[] diff --git a/src/deprecated.jl b/src/deprecated.jl index e99a83852..a2a3bb8a2 100644 --- a/src/deprecated.jl +++ b/src/deprecated.jl @@ -1,12 +1,167 @@ import Base: @deprecate -## deprecated in v0.13 -## deprecated in v0.12 +## deprecated in v0.16 -@deprecate(Cart_get(comm, maxdims), Cart_get(comm), false) +@deprecate(Bcast!(buffer, count::Integer, root::Integer, comm::Comm), + Bcast!(view(buffer, 1:count), root, comm), false) -@deprecate(Cart_coords!(comm, rank, maxdims::Integer, coords), - Cart_coords!(comm, rank, coords), false) +@deprecate(Scatter!(sendbuf, recvbuf, count::Integer, root::Integer, comm::Comm), + Scatter!(UBuffer(sendbuf, count), recvbuf, root, comm), false) +@deprecate(Scatter!(sendbuf, recvbuf::AbstractArray, count::Integer, root::Integer, comm::Comm), + Scatter!(UBuffer(sendbuf, count), view(recvbuf, 1:count), root, comm), false) +@deprecate(Scatter!(sendbuf, recvbuf::Nothing, count::Integer, root::Integer, comm::Comm), + Scatter!(UBuffer(sendbuf, count), MPI.IN_PLACE, root, comm), false) +@deprecate(Scatter!(sendbuf::Nothing, recvbuf, count::Integer, root::Integer, comm::Comm), + Scatter!(nothing, recvbuf, root, comm), false) +@deprecate(Scatter!(sendbuf::Nothing, recvbuf::AbstractArray, count::Integer, root::Integer, comm::Comm), + Scatter!(nothing, view(recvbuf, 1:count), root, comm), false) -@deprecate(Cart_coords(comm, maxdims::Integer), Cart_coords(comm), false) +@deprecate(Scatter(sendbuf, count::Integer, root::Integer, comm::Comm), + MPI.Comm_rank(comm) == root ? + Scatter!(UBuffer(sendbuf, count), similar(sendbuf, count), root, comm) : + Scatter!(nothing, similar(sendbuf, count), root, comm), false) + +@deprecate(Scatterv!(sendbuf, recvbuf, counts::Vector, root::Integer, comm::Comm), + MPI.Comm_rank(comm) == root ? + Scatterv!(VBuffer(sendbuf, counts), recvbuf, root, comm) : + Scatterv!(nothing, recvbuf, root, comm), false) +@deprecate(Scatterv!(sendbuf, recvbuf::AbstractArray, counts::Vector, root::Integer, comm::Comm), + MPI.Comm_rank(comm) == root ? + Scatterv!(VBuffer(sendbuf, counts), view(recvbuf,1:counts[MPI.Comm_rank(comm)+1]), root, comm) : + Scatterv!(nothing, view(recvbuf,1:counts[MPI.Comm_rank(comm)+1]), root, comm), false) +@deprecate(Scatterv!(sendbuf, recvbuf::Nothing, counts::Vector, root::Integer, comm::Comm), + MPI.Comm_rank(comm) == root ? + Scatterv!(VBuffer(sendbuf, counts), MPI.IN_PLACE, root, comm) : + Scatterv!(nothing, MPI.IN_PLACE, root, comm), false) +@deprecate(Scatterv!(sendbuf::Nothing, recvbuf, counts::Vector, root::Integer, comm::Comm), + Scatterv!(nothing, recvbuf, root, comm), false) +@deprecate(Scatterv!(sendbuf::Nothing, recvbuf::AbstractArray, counts::Vector, root::Integer, comm::Comm), + Scatterv!(nothing, view(recvbuf,1:counts[MPI.Comm_rank(comm)+1]), root, comm), false) + +@deprecate(Scatterv(sendbuf, counts::Vector, root::Integer, comm::Comm), + MPI.Comm_rank(comm) == root ? + Scatterv!(VBuffer(sendbuf, counts), similar(sendbuf, counts[Comm_rank(comm) + 1]), root, comm) : + Scatterv!(nothing, similar(sendbuf, counts[Comm_rank(comm) + 1]), root, comm), false) + +@deprecate(Gather!(sendbuf, recvbuf, count::Integer, root::Integer, comm::Comm), + Gather!(sendbuf, UBuffer(recvbuf, count), root, comm), false) +@deprecate(Gather!(sendbuf::AbstractArray, recvbuf, count::Integer, root::Integer, comm::Comm), + Gather!(view(sendbuf, 1:count), UBuffer(recvbuf, count), root, comm), false) +@deprecate(Gather!(sendbuf::Nothing, recvbuf, count::Integer, root::Integer, comm::Comm), + Gather!(MPI.IN_PLACE, UBuffer(recvbuf, count), root, comm), false) +@deprecate(Gather!(sendbuf, recvbuf::Nothing, count::Integer, root::Integer, comm::Comm), + Gather!(sendbuf, nothing, root, comm), false) +@deprecate(Gather!(sendbuf::AbstractArray, recvbuf::Nothing, count::Integer, root::Integer, comm::Comm), + Gather!(view(sendbuf, 1:count), nothing, root, comm), false) + +@deprecate(Gather(sendbuf, count::Integer, root::Integer, comm::Comm), + Gather(sendbuf, root, comm), false) +@deprecate(Gather(sendbuf::AbstractArray, count::Integer, root::Integer, comm::Comm), + Gather(view(sendbuf, 1:count), root, comm), false) + +@deprecate(Gatherv!(sendbuf, recvbuf, counts::Vector{Cint}, root::Integer, comm::Comm), + Gatherv!(sendbuf, VBuffer(recvbuf, counts), root, comm), false) +@deprecate(Gatherv!(sendbuf::AbstractArray, recvbuf, counts::Vector{Cint}, root::Integer, comm::Comm), + Gatherv!(view(sendbuf,1:counts[MPI.Comm_rank(comm)+1]), VBuffer(recvbuf, counts), root, comm), false) +@deprecate(Gatherv!(sendbuf::Nothing, recvbuf, counts::Vector{Cint}, root::Integer, comm::Comm), + Gatherv!(MPI.IN_PLACE, VBuffer(recvbuf, counts), root, comm), false) +@deprecate(Gatherv!(sendbuf, recvbuf::Nothing, counts::Vector{Cint}, root::Integer, comm::Comm), + Gatherv!(sendbuf, nothing, root, comm), false) +@deprecate(Gatherv!(sendbuf::AbstractArray, recvbuf::Nothing, counts::Vector{Cint}, root::Integer, comm::Comm), + Gatherv!(view(sendbuf,1:counts[MPI.Comm_rank(comm)+1]), nothing, root, comm), false) + +@deprecate(Gatherv(sendbuf, counts::Vector{Cint}, root::Integer, comm::Comm), + Gatherv!(sendbuf, Comm_rank(comm) == root ? VBuffer(similar(sendbuf, sum(counts)), counts) : nothing, root, comm), false) +@deprecate(Gatherv(sendbuf::AbstractArray, counts::Vector{Cint}, root::Integer, comm::Comm), + Gatherv!(view(sendbuf,1:counts[MPI.Comm_rank(comm)+1]), Comm_rank(comm) == root ? VBuffer(similar(sendbuf, sum(counts)), counts) : nothing, root, comm), false) + +@deprecate(Allgather!(sendbuf, recvbuf, count::Integer, comm::Comm), + Allgather!(sendbuf, UBuffer(recvbuf, count), comm::Comm), false) +@deprecate(Allgather!(sendbuf::AbstractArray, recvbuf, count::Integer, comm::Comm), + Allgather!(view(sendbuf,1:count), UBuffer(recvbuf, count), comm::Comm), false) +@deprecate(Allgather!(sendbuf::Nothing, recvbuf, count::Integer, comm::Comm), + Allgather!(MPI.IN_PLACE, UBuffer(recvbuf, count), comm::Comm), false) +@deprecate(Allgather!(sendrecvbuf, count::Integer, comm::Comm), + Allgather!(UBuffer(sendrecvbuf, count), comm), false) + +@deprecate(Allgather(sendbuf, count::Integer, comm::Comm), + Allgather(view(sendbuf, 1:count), comm), false) + +@deprecate(Allgatherv!(sendbuf, recvbuf, counts::Vector{Cint}, comm::Comm), + Allgatherv!(sendbuf, VBuffer(recvbuf, counts), comm), false) +@deprecate(Allgatherv!(sendbuf::AbstractArray, recvbuf, counts::Vector{Cint}, comm::Comm), + Allgatherv!(view(sendbuf,1:counts[MPI.Comm_rank(comm)+1]), VBuffer(recvbuf, counts), comm), false) +@deprecate(Allgatherv!(sendrecvbuf, counts::Vector{Cint}, comm::Comm), + Allgatherv!(VBuffer(sendrecvbuf, counts), comm), false) + +@deprecate(Allgatherv(sendbuf, counts::Vector{Cint}, comm::Comm), + Allgatherv!(sendbuf, VBuffer(similar(sendbuf, sum(counts)), counts), comm), false) + +@deprecate(Alltoall!(sendbuf, recvbuf, count::Integer, comm::Comm), + Alltoall!(UBuffer(sendbuf, count), UBuffer(recvbuf, count), comm), false) +@deprecate(Alltoall!(sendbuf::InPlace, recvbuf, count::Integer, comm::Comm), + Alltoall!(sendbuf, UBuffer(recvbuf, count), comm), false) +@deprecate(Alltoall!(sendrecvbuf, count::Integer, comm::Comm), + Alltoall!(UBuffer(sendrecvbuf, count), comm), false) + +@deprecate(Alltoall(sendbuf, count::Integer, comm::Comm), + Alltoall(UBuffer(sendbuf, count), comm), false) + +@deprecate(Alltoallv!(sendbuf, recvbuf, sendcounts, recvcounts, comm::Comm), + Alltoallv!(VBuffer(sendbuf, sendcounts), VBuffer(recvbuf, recvcounts), comm), false) + +@deprecate(Alltoallv(buf, sendcounts, recvcounts, comm::Comm), + Alltoallv!(VBuffer(buf, sendcounts), VBuffer(similar(buf,sum(recvcounts)), recvcounts), comm), false) + +@deprecate(Reduce!(send::AbstractArray, recv::AbstractArray, count::Integer, op, root::Integer, comm::Comm), + Reduce!(view(send, 1:count), view(recv, 1:count), op, root, comm), false) +@deprecate(Reduce!(send, recv::AbstractArray, count::Integer, op, root::Integer, comm::Comm), + Reduce!(send, view(recv, 1:count), op, root, comm), false) +@deprecate(Reduce!(send::AbstractArray, recv, count::Integer, op, root::Integer, comm::Comm), + Reduce!(view(send, 1:count), recv, op, root, comm), false) +@deprecate(Reduce!(send, recv, count::Integer, op, root::Integer, comm::Comm), + Reduce!(send, recv, op, root, comm), false) +@deprecate(Reduce!(sendrecvbuf::AbstractArray, count::Integer, op, root::Integer, comm::Comm), + Reduce!(view(sendrecvbuf, 1:count), op, root, comm), false) +@deprecate(Reduce!(sendrecvbuf, count::Integer, op, root::Integer, comm::Comm), + Reduce!(sendrecvbuf, op, root, comm), false) + +@deprecate(Allreduce!(send::AbstractArray, recv::AbstractArray, count::Integer, op, comm::Comm), + Allreduce!(view(send, 1:count), view(recv, 1:count), op, comm), false) +@deprecate(Allreduce!(send, recv::AbstractArray, count::Integer, op, comm::Comm), + Allreduce!(send, view(recv, 1:count), op, comm), false) +@deprecate(Allreduce!(send::AbstractArray, recv, count::Integer, op, comm::Comm), + Allreduce!(view(send, 1:count), recv, op, comm), false) +@deprecate(Allreduce!(send, recv, count::Integer, op, comm::Comm), + Allreduce!(send, recv, op, comm), false) +@deprecate(Allreduce!(sendrecvbuf::AbstractArray, count::Integer, op, comm::Comm), + Allreduce!(view(sendrecvbuf, 1:count), op, comm), false) +@deprecate(Allreduce!(sendrecvbuf, count::Integer, op, comm::Comm), + Allreduce!(sendrecvbuf, op, comm), false) + +@deprecate(Scan!(sendbuf::AbstractArray, recvbuf::AbstractArray, count::Integer, op, comm::Comm), + Scan!(view(sendbuf, 1:count), view(recvbuf, 1:count), op, comm), false) +@deprecate(Scan!(sendbuf, recvbuf::AbstractArray, count::Integer, op, comm::Comm), + Scan!(sendbuf, view(recvbuf, 1:count), op, comm), false) +@deprecate(Scan!(sendbuf::AbstractArray, recvbuf, count::Integer, op, comm::Comm), + Scan!(view(sendbuf, 1:count), recvbuf, op, comm), false) +@deprecate(Scan!(sendbuf, recvbuf, count::Integer, op, comm::Comm), + Scan!(sendbuf, recvbuf, op, comm), false) +@deprecate(Scan!(sendrecvbuf::AbstractArray, count::Integer, op, comm::Comm), + Scan!(view(sendrecvbuf, 1:count), op, comm), false) +@deprecate(Scan!(sendrecvbuf, count::Integer, op, comm::Comm), + Scan!(sendrecvbuf, op, comm), false) + +@deprecate(Exscan!(sendbuf::AbstractArray, recvbuf::AbstractArray, count::Integer, op, comm::Comm), + Exscan!(view(sendbuf, 1:count), view(recvbuf, 1:count), op, comm), false) +@deprecate(Exscan!(sendbuf, recvbuf::AbstractArray, count::Integer, op, comm::Comm), + Exscan!(sendbuf, view(recvbuf, 1:count), op, comm), false) +@deprecate(Exscan!(sendbuf::AbstractArray, recvbuf, count::Integer, op, comm::Comm), + Exscan!(view(sendbuf, 1:count), recvbuf, op, comm), false) +@deprecate(Exscan!(sendbuf, recvbuf, count::Integer, op, comm::Comm), + Exscan!(sendbuf, recvbuf, op, comm), false) +@deprecate(Exscan!(sendrecvbuf::AbstractArray, count::Integer, op, comm::Comm), + Exscan!(view(sendrecvbuf, 1:count), op, comm), false) +@deprecate(Exscan!(sendrecvbuf, count::Integer, op, comm::Comm), + Exscan!(sendrecvbuf, op, comm), false) diff --git a/src/pointtopoint.jl b/src/pointtopoint.jl index 3673b1e06..42d03b7b0 100644 --- a/src/pointtopoint.jl +++ b/src/pointtopoint.jl @@ -215,7 +215,7 @@ end Starts a nonblocking send of `data` to MPI rank `dest` of communicator `comm` using with the message tag `tag`. -`data` can be a `Buffer`, or any object for which [`Buffer_send`](@ref) is defined. +`data` can be a [`Buffer`](@ref), or any object for which [`Buffer_send`](@ref) is defined. Returns the [`Request`](@ref) object for the nonblocking send. diff --git a/test/test_allgather.jl b/test/test_allgather.jl index 6c1e5b20c..ceece8501 100644 --- a/test/test_allgather.jl +++ b/test/test_allgather.jl @@ -32,23 +32,28 @@ for T in Base.uniontypes(MPI.MPIDatatype) # Test passing output buffer with set size A = ArrayType(T[val]) + C = ArrayType{T}(undef, size) - MPI.Allgather!(A, C, length(A), comm) + MPI.Allgather!(A, C, comm) # implied size + @test C == ArrayType{T}(1:size) + + C = ArrayType{T}(undef, size) + MPI.Allgather!(A, UBuffer(C,1), comm) @test C == ArrayType{T}(1:size) # Test assertion error C = ArrayType{T}(undef, size-1) - @test_throws AssertionError MPI.Allgather!(A, C, length(A), comm) + @test_throws AssertionError MPI.Allgather!(A, C, comm) # Test explicit IN_PLACE C = ArrayType{T}([i == rank ? i : size + 1 for i = 0:size-1]) - MPI.Allgather!(MPI.IN_PLACE, C, 1, comm) + MPI.Allgather!(MPI.IN_PLACE, UBuffer(C, 1), comm) @test C isa ArrayType{T,1} @test C == ArrayType{T}(0:size-1) # Test IN_PLACE C = ArrayType{T}([i == rank ? i : size + 1 for i = 0:size-1]) - MPI.Allgather!(C, 1, comm) + MPI.Allgather!(UBuffer(C, 1), comm) @test C isa ArrayType{T,1} @test C == ArrayType{T}(0:size-1) end diff --git a/test/test_allgatherv.jl b/test/test_allgatherv.jl index 2b1f6388f..e8fa2c5fd 100644 --- a/test/test_allgatherv.jl +++ b/test/test_allgatherv.jl @@ -19,22 +19,19 @@ check = collect(Iterators.flatten([fill(r, counts[r+1]) for r = 0:size-1])) for T in Base.uniontypes(MPI.MPIDatatype) A = ArrayType{T}(fill(T(rank), counts[rank+1])) - B = MPI.Allgatherv(A, counts, comm) - @test B isa ArrayType{T} - @test B == ArrayType{T}(check) # Test passing the output buffer B = ArrayType{T}(undef, sum(counts)) - MPI.Allgatherv!(A, B, counts, comm) + MPI.Allgatherv!(A, VBuffer(B, counts), comm) @test B == ArrayType{T}(check) # Test assertion when output size is too small B = ArrayType{T}(undef, sum(counts)-1) - @test_throws AssertionError MPI.Allgatherv!(A, B, counts, comm) + @test_throws AssertionError MPI.Allgatherv!(A, VBuffer(B, counts), comm) # Test explicit MPI_IN_PLACE B = ArrayType(fill(T(rank), sum(counts))) - MPI.Allgatherv!(MPI.IN_PLACE, B, counts, comm) + MPI.Allgatherv!(MPI.IN_PLACE, VBuffer(B, counts), comm) @test B == ArrayType{T}(check) end diff --git a/test/test_allreduce.jl b/test/test_allreduce.jl index e882d969b..dd252d06e 100644 --- a/test/test_allreduce.jl +++ b/test/test_allreduce.jl @@ -36,7 +36,6 @@ for T = [Int] # Assertions when output buffer too small recv_arr = ArrayType{T}(undef, size(send_arr).-1) @test_throws AssertionError MPI.Allreduce!(send_arr, recv_arr, - length(send_arr), op, MPI.COMM_WORLD) # IN_PLACE recv_arr = copy(send_arr) diff --git a/test/test_alltoall.jl b/test/test_alltoall.jl index 673468e74..875fcabee 100644 --- a/test/test_alltoall.jl +++ b/test/test_alltoall.jl @@ -18,19 +18,19 @@ for T in Base.uniontypes(MPI.MPIDatatype) # Allocating version a = ArrayType(fill(T(rank), size)) - b = MPI.Alltoall(a, 1, comm) + b = MPI.Alltoall(UBuffer(a, 1), comm) @test b isa ArrayType{T} @test b == ArrayType{T}(0:size-1) # Non Allocating version a = ArrayType(fill(T(rank),size)) b = ArrayType{T}(undef, size*1) - MPI.Alltoall!(a, b, 1, comm) + MPI.Alltoall!(UBuffer(a,1), UBuffer(b,1), comm) @test b == ArrayType{T}(0:size-1) # IN_PLACE version a = ArrayType{T}(fill(T(rank),size)) - MPI.Alltoall!(MPI.IN_PLACE, a, 1, comm) + MPI.Alltoall!(MPI.IN_PLACE, UBuffer(a,1), comm) @test a == ArrayType{T}(0:size-1) end diff --git a/test/test_alltoallv.jl b/test/test_alltoallv.jl index 3360c467d..c7996385a 100644 --- a/test/test_alltoallv.jl +++ b/test/test_alltoallv.jl @@ -24,19 +24,14 @@ for T in Base.uniontypes(MPI.MPIDatatype) A = ArrayType{T}(send_vals) - # Allocating version - B = MPI.Alltoallv(A, send_counts, recv_counts, comm) - @test B isa ArrayType{T} - @test B == ArrayType{T}(recv_vals) - # Non Allocating version C = ArrayType{T}(undef, sum(recv_counts)) - MPI.Alltoallv!(A, C, send_counts, recv_counts, comm) + MPI.Alltoallv!(VBuffer(A,send_counts), VBuffer(C,recv_counts), comm) @test C == ArrayType{T}(recv_vals) # Test assertion on wrong output buffer length C = ArrayType{T}(undef, sum(recv_counts)-1) - @test_throws AssertionError MPI.Alltoallv!(A, C, send_counts, recv_counts, comm) + @test_throws AssertionError MPI.Alltoallv!(VBuffer(A,send_counts), VBuffer(C,recv_counts), comm) end MPI.Finalize() diff --git a/test/test_exscan.jl b/test/test_exscan.jl index 7431e8933..a3b675447 100644 --- a/test/test_exscan.jl +++ b/test/test_exscan.jl @@ -27,7 +27,7 @@ for T in setdiff([Base.uniontypes(MPI.MPIDatatype)...], [Char, Int8, UInt8]) B = MPI.Exscan(A, *, comm) @test B isa ArrayType{T} - MPI.Exscan!(A, length(A), *, comm) + MPI.Exscan!(A, *, comm) if rank > 0 @test A == ArrayType{T}(fill(T(prodrank), 4)) end diff --git a/test/test_gather.jl b/test/test_gather.jl index 0673f9434..96825a014 100644 --- a/test/test_gather.jl +++ b/test/test_gather.jl @@ -31,14 +31,18 @@ for T in Base.uniontypes(MPI.MPIDatatype) if isroot @test C isa Vector{T} @test C == Vector{T}(1:sz) + else + @test C === nothing end - # Allocating, explicit length - A = ArrayType{T}([MPI.Comm_rank(comm) + 1, 0]) - C = MPI.Gather(A, 1, root, comm) + # Allocating, array + A = ArrayType{T}([MPI.Comm_rank(comm) + 1]) + C = MPI.Gather(A, root, comm) if isroot @test C isa ArrayType{T} @test C == ArrayType{T}(1:MPI.Comm_size(comm)) + else + @test C === nothing end # Allocating, view @@ -52,7 +56,11 @@ for T in Base.uniontypes(MPI.MPIDatatype) # Non Allocating A = ArrayType{T}(fill(T(rank+1), 4)) C = ArrayType{T}(undef, 4sz) - MPI.Gather!(A, C, length(A), root, comm) + MPI.Gather!(A, C, root, comm) + if isroot + @test C == ArrayType{T}(repeat(1:sz, inner=4)) + end + MPI.Gather!(A, UBuffer(C,4), root, comm) if isroot @test C == ArrayType{T}(repeat(1:sz, inner=4)) end @@ -63,9 +71,9 @@ for T in Base.uniontypes(MPI.MPIDatatype) A = ArrayType{T}(fill(T(rank+1), 4*sz)) end if root == MPI.Comm_rank(comm) - MPI.Gather!(nothing, A, 4, root, comm) + MPI.Gather!(MPI.IN_PLACE, UBuffer(A,4), root, comm) else - MPI.Gather!(A, nothing, 4, root, comm) + MPI.Gather!(A, nothing, root, comm) end if isroot @test A == ArrayType{T}(repeat(1:sz, inner=4)) diff --git a/test/test_gatherv.jl b/test/test_gatherv.jl index e44c8b3cd..8cff5afa7 100644 --- a/test/test_gatherv.jl +++ b/test/test_gatherv.jl @@ -21,15 +21,10 @@ check = collect(Iterators.flatten([fill(r, counts[r+1]) for r = 0:size-1])) for T in Base.uniontypes(MPI.MPIDatatype) A = ArrayType(fill(T(rank), mod(rank,2) + 1)) - B = MPI.Gatherv(A, counts, root, comm) - if isroot - @test B isa ArrayType{T} - @test B == ArrayType{T}(check) - end # Test passing the output buffer B = ArrayType{T}(undef, sum(counts)) - MPI.Gatherv!(A, B, counts, root, comm) + MPI.Gatherv!(A, isroot ? VBuffer(B, counts) : nothing, root, comm) if isroot @test B == ArrayType{T}(check) end @@ -37,15 +32,16 @@ for T in Base.uniontypes(MPI.MPIDatatype) # Test assertion when output size is too small B = ArrayType{T}(undef, sum(counts)-1) if isroot - @test_throws AssertionError MPI.Gatherv!(A, B, counts, root, comm) + @test_throws AssertionError MPI.Gatherv!(A, VBuffer(B, counts), root, comm) end # Test explicit MPI_IN_PLACE - B = ArrayType(fill(T(rank), sum(counts))) - if root == MPI.Comm_rank(comm) - MPI.Gatherv!(nothing, B, counts, root, comm) + if isroot + B = ArrayType(fill(T(rank), sum(counts))) + MPI.Gatherv!(MPI.IN_PLACE, VBuffer(B, counts), root, comm) else - MPI.Gatherv!(B, nothing, counts, root, comm) + B = ArrayType(fill(T(rank), counts[rank+1])) + MPI.Gatherv!(B, nothing, root, comm) end if isroot @test B == ArrayType{T}(check) diff --git a/test/test_reduce.jl b/test/test_reduce.jl index acc1c5eda..b05dfa36f 100644 --- a/test/test_reduce.jl +++ b/test/test_reduce.jl @@ -66,7 +66,6 @@ for T = [Int] recv_arr = ArrayType{T}(undef, size(send_arr).-1) if isroot @test_throws AssertionError MPI.Reduce!(send_arr, recv_arr, - length(send_arr), op, root, MPI.COMM_WORLD) end @@ -79,9 +78,9 @@ for T = [Int] end # Allocating version - val = MPI.Reduce(2, op, root, MPI.COMM_WORLD) + r = MPI.Reduce(2, op, root, MPI.COMM_WORLD) if isroot - @test val == sz*2 + @test r == sz*2 end recv_arr = MPI.Reduce(send_arr, op, root, MPI.COMM_WORLD) diff --git a/test/test_scatter.jl b/test/test_scatter.jl index 3b81a4845..e2ec4e5d1 100644 --- a/test/test_scatter.jl +++ b/test/test_scatter.jl @@ -18,29 +18,29 @@ isroot = rank == root for T in Base.uniontypes(MPI.MPIDatatype) - # Allocating version A = isroot ? ArrayType{T}(1:size) : ArrayType{T}(undef, 1) - B = MPI.Scatter(A, 1, root, comm) - @test B isa ArrayType{T} - @test Array(B)[1] == T(rank+1) # Non Allocating version B = ArrayType{T}(undef, 1) - MPI.Scatter!(A, B, 1, root, comm) + MPI.Scatter!(A, B, root, comm) @test Array(B)[1] == T(rank+1) # In place version B = isroot ? copy(A) : ArrayType{T}(undef, 1) if root == MPI.Comm_rank(comm) - MPI.Scatter!(B, nothing, 1, root, comm) + MPI.Scatter!(UBuffer(B,1), MPI.IN_PLACE, root, comm) else - MPI.Scatter!(nothing, B, 1, root, comm) + MPI.Scatter!(nothing, B, root, comm) end @test Array(B)[1] == T(rank+1) # Test throwing - B = ArrayType{T}(undef, 0) - @test_throws AssertionError MPI.Scatter!(A, B, 1, root, comm) + if isroot + B = ArrayType{T}(undef, 0) + @test_throws DivideError MPI.Scatter!(A, B, root, comm) + B = ArrayType{T}(undef, 8) + @test_throws AssertionError MPI.Scatter!(A, B, root, comm) + end end MPI.Finalize() diff --git a/test/test_scatterv.jl b/test/test_scatterv.jl index 39f3d5d72..07465b625 100644 --- a/test/test_scatterv.jl +++ b/test/test_scatterv.jl @@ -20,24 +20,19 @@ counts = Cint[mod(i,2) + 1 for i in 0:(size-1)] ref = collect(Iterators.flatten([fill(r, counts[r+1]) for r = 0:size-1])) for T in Base.uniontypes(MPI.MPIDatatype) - A = isroot ? ArrayType{T}(ref) : ArrayType{T}(undef, 1) - - # Allocating - B = MPI.Scatterv(A, counts, root, comm) - @test B isa ArrayType{T} - @test B == ArrayType{T}(fill(rank, counts[rank+1])) + A = isroot ? ArrayType{T}(ref) : nothing # Non Allocating B = ArrayType{T}(undef, counts[rank+1]) - MPI.Scatterv!(A, B, counts, root, comm) + MPI.Scatterv!(isroot ? VBuffer(A, counts) : nothing, B, root, comm) @test B == ArrayType{T}(fill(rank, counts[rank+1])) # IN_PLACE B = isroot ? copy(A) : ArrayType{T}(undef, counts[rank+1]) if root == MPI.Comm_rank(comm) - MPI.Scatterv!(B, nothing, counts, root, comm) + MPI.Scatterv!(VBuffer(B, counts), MPI.IN_PLACE, root, comm) else - MPI.Scatterv!(nothing, B, counts, root, comm) + MPI.Scatterv!(nothing, B, root, comm) end if isroot @test B == A