fixes for non-local machines

README.md

@@ -32,26 +32,31 @@ or equivalently:
 julia> Pkg.add("MPI", "SparseBase", "MatrixMarket")
 ```
 
+The binaries provided with Julia may not work correctly on your machine (this is currently affecting at least NERSC Perlmutter). 
+If this is the case you may provide your own library with:
+
+```julia
+SuperLUDIST.set_libraries!(; libsuperlu_dist_Int32 = <PATH TO 32 BIT LIBRARY>, libsuperlu_dist_Int64 = <PATH TO 64 BIT LIBRARY>)
+``````
+
+These libraries should be built with the proper integer size (32 bit or 64 bit), 
+and should be linked against the MPI used by MPI.jl (on HPC clusters this should be your system MPI).
+
+You may provide one or the other, or both, but if the 32 or 64 bit integer libraries are unavailable 32 bit or 64 bit matrices will be unavailable respectively.
+
 ## Usage
 Examples for running in replicated and distributed mode are provided in the examples directory. To run these examples
 follow the instructions provided [here](https://juliaparallel.org/MPI.jl/latest/configuration/) to set up your MPI correctly. In particular `mpiexecjl` should be in your path (it is typically found in `~/.julia/bin`)
 
 Then you can invoke a driver example:
 ```
-mpiexecjl -n <nprocs> julia --project examples/pdrive.jl
+mpiexecjl -n <nprocs> julia --project examples/basic_example.jl
 ```
 
 This command does a couple things. The first part launches MPI with `<nprocs>` ranks. The second part `julia --project` activates the current folder. If your installation is in the global environment you may omit `--project`, or provide a path to the environment you would like to use: `--project=~/MySuperLUProject`. Finally `examples/pdrive.jl` is the Julia script to be executed under MPI. 
 
 Driver routine examples are found in the `examples` folder of this repository. The examples currently cover basic usage of `pgssvx` and `pgssvx_ABglobal` driver routines. These routines load matrix market files and generate right hand sides, which is not a typical use-case. Instead users will often build a submatrix on each rank, and construct a `DistributedSuperMatrix`.
 
-More examples will be added over the coming days to illustrate:
-- DistributedSuperMatrix and ReplicatedSuperMatrix construction.
-- Re-using LU factorizations.
-- Use of the computational routines.
-- 3D distribution routines.
-- Extensions for `SparseMatricesCSR.jl`, `LinearSolve.jl`, `PartitionedArrays.jl` and several other packages.
-
 ## Known Issues
 - CUDA support is currently disabled.
 - OpenMP will often oversubscribe single node setups. In this case you should use `SuperLUDIST.superlu_set_num_threads(Int32, <num>)` or `SuperLUDIST.superlu_set_num_threads(Int64, <num>)` to the number of OpenMP threads desired on each rank.

examples/highlevel.jl → examples/basic_example.jl

@@ -12,15 +12,15 @@ MPI.Init()
 nprow, npcol, nrhs = 2, 2, 4
 root = 0
 comm = MPI.COMM_WORLD
-grid = Grid{Int32}(nprow, npcol, comm)
+grid = Grid{Int64}(nprow, npcol, comm)
 iam = grid.iam
 isroot = iam == root
 
 # Utility function for reading a .mtx file and generating suitable
 # rhs and x for testing. 
 # coo is held only on root, b and xtrue are replicated on each rank.
 coo, b, xtrue = SuperLUDIST.mmread_and_generatesolution(
-    Float64, Int32, nrhs, joinpath(@__DIR__, "add32.mtx"), grid; root
+    Float64, Int64, nrhs, joinpath(@__DIR__, "add32.mtx"), grid; root
 )
 
 # A second rhs and xtrue for testing different sized b.
@@ -31,19 +31,26 @@ csr = isroot ? convert(SparseBase.CSRStore, coo) : nothing
 chunksizes = isroot ? distribute_evenly(size(csr, 1), nprow * npcol) : nothing
 
 # on single nodes this will help prevent oversubscription of threads.
-SuperLUDIST.superlu_set_num_threads(Int64, 2)
+# SuperLUDIST.superlu_set_num_threads(Int64, 1)
 
 # If constructing from existing per-node data the following constructors will help:
 # A = DistributedSuperMatrix(store::CSRStore, firstrow, globalsize::NTuple{2, Int})
 # store = CSRStore(ptrs, indices, values, localsize::NTuple{2, Int})
 # @show iam csr
 A = Communication.scatterstore!(
-    DistributedSuperMatrix{Float64, Int32}(grid), csr, chunksizes; root);
+    DistributedSuperMatrix{Float64, Int64}(grid), csr, chunksizes; root);
 
 b_local = b[A.first_row : A.first_row + localsize(A, 1) - 1, :] # shrink b
 xtrue_local = xtrue[A.first_row : A.first_row + localsize(A, 1) - 1, :] # shrink xtrue
+
+# Form a factorization object.
+# Alternatively, to solve and factor simultaneously call:
+# b_local, F = pgssvx!(A, b_local)
 F = lu!(A);
 
+# Out of place solve
+# Alternatively for in-place call: for in-place.
+# ldiv!(A, b_local) 
 b_local = F \ b_local
 
 if !(iam == root) || (nprow * npcol == 1)
@@ -54,6 +61,7 @@ end
 b2_local = b2[A.first_row : A.first_row + localsize(A, 1) - 1, :] # shrink b2
 xtrue2_local = x2[A.first_row : A.first_row + localsize(A, 1) - 1, :] # shrink xtrue2
 
+# Solve again, reusing the initial facctorization:
 b2_local = F \ b2_local
 
 if !(iam == root) || (nprow * npcol == 1)

lib/libsuperlu_dist32.jl

@@ -1,9 +1,9 @@
 module SuperLU_Int32
-using SuperLU_DIST_jll
 import MPI: MPI_Comm, MPI_Request, MPI_Datatype, MPI_Errhandler
 using ..SuperLUDIST_Common
+import ..SuperLUDIST: libsuperlu_dist_Int32
 using SuperLUBase.Common
-const libsuperlu_dist_Int32 = SuperLU_DIST_jll.libsuperlu_dist_Int32
+
 function superlu_abort_and_exit_dist(arg1)
     @ccall libsuperlu_dist_Int32.superlu_abort_and_exit_dist(arg1::Ptr{Cchar})::Cvoid
 end
@@ -3289,4 +3289,4 @@ const DEG_TREE = 2
 
 const MAX_LOOKAHEADS = 50
 
-end
+end

lib/libsuperlu_dist64.jl

@@ -2,8 +2,8 @@ module SuperLU_Int64
 import MPI: MPI_Comm, MPI_Request, MPI_Datatype, MPI_Errhandler
 using SuperLU_DIST_jll
 using ..SuperLUDIST_Common
+import ..SuperLUDIST: libsuperlu_dist_Int64
 using SuperLUBase.Common
-const libsuperlu_dist_Int64 = SuperLU_DIST_jll.libsuperlu_dist_Int64
 
 function superlu_abort_and_exit_dist(arg1)
     @ccall libsuperlu_dist_Int64.superlu_abort_and_exit_dist(arg1::Ptr{Cchar})::Cvoid
@@ -3290,4 +3290,4 @@ const THR_N = BLK_N ÷ DIM_Y
 
 const DEG_TREE = 2
 
-end
+end

src/SuperLUDIST.jl

@@ -1,7 +1,7 @@
 module SuperLUDIST
 using SparseBase
 using MPI
-using SuperLU_DIST_jll
+include("libs.jl")
 using CIndices: CIndex
 using DocStringExtensions
 using MatrixMarket
@@ -17,7 +17,7 @@ function __init__()
         LinearAlgebra.BLAS.lbt_forward(OpenBLAS32_jll.libopenblas_path)
       end
     end
-    LinearAlgebra.BLAS.set_num_threads(1)
+    # LinearAlgebra.BLAS.set_num_threads(1)
     #superlu_set_num_threads(Int32, 1)
     #superlu_set_num_threads(Int64, 1)
 end

src/libs.jl

@@ -0,0 +1,39 @@
+using Libdl: dlopen, dlclose, dlpath, dlsym, RTLD_DEEPBIND, RTLD_LAZY
+using Preferences
+
+const _PREFERENCE_SLUINT32 = @load_preference("libsuperlu_dist_Int32", nothing)
+const _PREFERENCE_SLUINT64 = @load_preference("libsuperlu_dist_Int64", nothing)
+
+flags = RTLD_DEEPBIND | RTLD_LAZY
+if _PREFERENCE_SLUINT32 === nothing && _PREFERENCE_SLUINT64 === nothing
+    using SuperLU_DIST_jll
+elseif _PREFERENCE_SLUINT32 === nothing
+    libsuperlu_dist_Int64 = _PREFERENCE_SLUINT64
+    dlopen(libsuperlu_dist_Int64, flags; throw_error=true)
+    libsuperlu_dist_Int32 = nothing
+elseif _PREFERENCE_SLUINT32 === nothing
+    libsuperlu_dist_Int32 = _PREFERENCE_SLUINT32
+    dlopen(libsuperlu_dist_Int32, flags; throw_error=true)
+    libsuperlu_dist_Int64 = nothing
+else
+    libsuperlu_dist_Int64 = _PREFERENCE_SLUINT64
+    dlopen(libsuperlu_dist_Int64, flags; throw_error=true)
+    libsuperlu_dist_Int32 = _PREFERENCE_SLUINT32
+    dlopen(libsuperlu_dist_Int32, flags; throw_error=true)
+end
+
+function set_libraries!(;libsuperlu_dist_Int32 = nothing, libsuperlu_dist_Int64 = nothing)
+    if isnothing(libsuperlu_dist_Int32)
+        @delete_preferences!("libsuperlu_dist_Int32")
+    else
+        isfile(libsuperlu_dist_Int32) || throw(ArgumentError("$libsuperlu_dist_Int32 is not a file that exists."))
+        @set_preferences!("libsuperlu_dist_Int32" => libsuperlu_dist_Int32)
+    end
+    if isnothing(libsuperlu_dist_Int64)
+        @delete_preferences!("libsuperlu_dist_Int64")
+    else
+        isfile(libsuperlu_dist_Int64) || throw(ArgumentError("$libsuperlu_dist_Int64 is not a file that exists."))
+        @set_preferences!("libsuperlu_dist_Int64" => libsuperlu_dist_Int64)
+    end
+    @info "Please restart Julia and reload SuperLUDIST.jl for the library changes to take effect."
+end

src/matrixmarket.jl

@@ -106,7 +106,7 @@ function mmread_and_generatesolution(Tv, Ti, nrhs, path, grid; root = 0)
   xtrue = Matrix{Tv}(undef, n, nrhs)
   b = Matrix{Tv}(undef, m, nrhs)
   if iam == root
-      SuperLUDIST.GenXtrue_dist!(xtrue, Int32)
+      SuperLUDIST.GenXtrue_dist!(xtrue, Ti)
       SuperLUDIST.FillRHS_dist!(b, Acsc, xtrue)
   end
   MPI.Bcast!(b, root, comm)

src/structs.jl

@@ -2,7 +2,6 @@
 ##############################
 # const Grid{I} = SuperLUDIST_Common.gridinfo_t{I}
 # const Grid3D{I} = SuperLUDIST_Common.gridinfo3d_t{I}
-
 struct Grid{I}
     comm::MPI.Comm
     gridinfo::SuperLUDIST_Common.gridinfo_t{I}
@@ -30,6 +29,11 @@ function gridinit!(g::Grid{$I}, nprow, npcol)
     return g
 end
 function Grid{$I}(nprow::Integer, npcol::Integer, comm = MPI.COMM_WORLD; batch = false, usermap = nothing)
+    $I == Int32 && libsuperlu_dist_Int32 === nothing && 
+        (throw(ArgumentError("libsuperlu_dist_Int32 is not loaded.")))
+    $I == Int64 && libsuperlu_dist_Int64 === nothing && 
+        (throw(ArgumentError("libsuperlu_dist_Int64 is not loaded.")))
+
     !MPI.Initialized() && MPI.Init()
     g = Grid{$I}(
         comm,