Open ChrisRackauckas opened 2 years ago
Hi @ChrisRackauckas,
PSparseMatrix(::Vector{Int64}, ::Vector{Int64}, ::Vector{Float64}, ::Vector{Int64}, ::Vector{Int64})
This method error is working as expected (idem for the one with PVector). The construciton of theese objects needs a bit more work. E.g., define how the rows are partitioned over the parts.
Here you have a quick example for PVector.
# mycode.jl
using PartitionedArrays
using LinearAlgebra
np = 4
parts = get_part_ids(sequential,np)
# Once everything works on a standard julia repl,
# change sequential->mpi and launch your code as
# $ mpiexecjl --project=. -n 4 julia mycode.jl
n = 10
# This range knows that we are partitioning 10 rows in 4 parts
rows = PRange(parts,n)
display(rows)
display(map_parts(k->k.lid_to_gid,rows.partition))
# Generate some values at each part
values = map_parts(rows.partition) do lids
rand(num_lids(lids))
end
display(values)
# Now, we have the values
# and we know how the rows are partitioned.
# We can construct a PVector
v = PVector(values,rows)
# Work with it
display(norm(v))
# Note that scalar indexing is not implemented (it would be VERY ineficient anyway)
v[3]
lids.lid_to_gid gets the global i
A more detailed example, including the creation of a matrix and solvers
using PartitionedArrays
using SparseArrays
using LinearAlgebra
# Create the "parallel" environment
np = 3
parts = get_part_ids(sequential,np)
a = map_parts(i->2*i,parts)
# Create a vector
n = 10
rows = PRange(parts,n)
values = map_parts(rows.partition) do lids
10.0*lids.lid_to_gid
end
v = PVector(values,rows)
# Use the vector
w = 2 .* v
z = w .+ v
t = similar(w,Int)
# Create a sparse matrix
# (a diagonal one for demo purposes)
Avalues = map_parts(rows.partition) do lids
n = length(lids.lid_to_gid)
I = collect(1:n)
J = collect(1:n)
V = ones(n)
sparse(I,J,V)
end
cols = rows
A = PSparseMatrix(Avalues,rows,cols)
# Solvers
x = A\w
factors = lu(A)
y = similar(w,Float64,cols)
ldiv!(y,factors,w)
lu!(factors,A)
ldiv!(y,factors,w)
In order to launch previous example with MPI, just replace the top lines with:
using PartitionedArrays
using SparseArrays
using LinearAlgebra
using MPI
MPI.Init()
# Create the "parallel" environment
np = 3
parts = get_part_ids(mpi,np)
and run it with $ ~/.julia/bin/mpiexecjl --project=. -n 3 julia demo.jl
I'm trying to make a version of the DifferentialEquations.jl Solving Large Stiff Equations PDE time stepping tutorial that is distributed, and am running into some issues figuring out how to construct PArrays. MWE with the method errors are below:
Then the solving code is:
But since the construction doesn't run right now those will all fail of course.