Ordering: expose METIS nseps as nd_nseps

ext/SparseDirectSolverMetisExt.jl

@@ -11,13 +11,14 @@ using SparseDirectSolver
 using CliqueTrees
 using Metis
 
-# METIS_NodeND with `ufactor = nd_ubfactor` (-1: METIS default). `nd_nlevels` is
-# cuDSS' *minimum* number of dissection levels; METIS_NodeND dissects recursively
+# METIS_NodeND with `ufactor = nd_ubfactor` and `nseps = nd_nseps` (-1: METIS
+# defaults; values < 1 fall back to the default). `nd_nlevels` is cuDSS'
+# *minimum* number of dissection levels; METIS_NodeND dissects recursively
 # until the parts are small, so the minimum holds whenever the graph is large
 # enough to be split that often. CliqueTrees' level-capped `ND{S}` (AMD below the
 # cap) was 4–13× slower and gave 15–80% more fill on Laplacians (T05 Report).
-function nd_algorithm(nlevels::Integer, ubfactor::Integer)
-    return CliqueTrees.METIS(; ufactor = ubfactor)
+function nd_algorithm(nlevels::Integer, ubfactor::Integer, nseps::Integer)
+    return CliqueTrees.METIS(; ufactor = ubfactor, nseps = nseps >= 1 ? nseps : -1)
 end
 
 function __init__()

src/options.jl

@@ -17,7 +17,7 @@ const CONFIG_PARAMETERS = ("reordering_alg", "factorization_alg", "solve_alg",
                            "hybrid_device_memory_limit", "use_cuda_register_memory", "host_nthreads",
                            "hybrid_execute_mode", "pivot_epsilon_alg", "nd_nlevels", "ubatch_size",
                            "ubatch_index", "use_superpanels", "device_count", "device_indices",
-                           "schur_mode", "deterministic_mode", "nd_ubfactor")
+                           "schur_mode", "deterministic_mode", "nd_ubfactor", "nd_nseps")
 
 """
     DATA_PARAMETERS
@@ -127,6 +127,7 @@ defaults below. Keyword arguments are applied through [`setparam!`](@ref), so
 | `schur_mode` | `0` | off |
 | `deterministic_mode` | `0` | atomic forward solve allowed |
 | `nd_ubfactor` | `-1` | ordering library default |
+| `nd_nseps` | `-1` | ordering library default (1 separator per split) |
 | `ir_mode` | `IR_PLAIN` | plain iterative refinement |
 | `factor_precision` | `nothing` | factors in the input precision |
 | `amalgamation` | `(max_width = 32, zero_fraction = 0.25, min_width = 8)` | |
@@ -176,6 +177,7 @@ mutable struct Options
     schur_mode::Int
     deterministic_mode::Int
     nd_ubfactor::Int
+    nd_nseps::Int
     # configuration parameters beyond cuDSS (PLAN §1.7)
     ir_mode::IRMode
     factor_precision::Union{Nothing, DataType}
@@ -202,7 +204,7 @@ mutable struct Options
         opts = new(
             REORDERING_DEFAULT, FACTORIZATION_DEFAULT, SOLVE_DEFAULT, MATCHING_NONE,
             0, 0, 0.0, PIVOT_AUTO, 0.01, nothing, PIVOT_EPSILON_DEFAULT, -1,
-            0, 0, 1, 0, 0, 10, 0, -1, 1, 0, 0, -1,
+            0, 0, 1, 0, 0, 10, 0, -1, 1, 0, 0, -1, -1,
             IR_PLAIN, nothing, DEFAULT_AMALGAMATION, SCHEDULE_AUTO, PIVOT_PAIRS_DEFAULT, PIVOT_PAIR_TOLERANCE,
             64, 512, copy(DEFAULT_SUBTREE_BUDGETS), SUBTREE_PARALLELISM, 0, -1,
             nothing, nothing, nothing, nothing, nothing, nothing,
@@ -286,6 +288,7 @@ const PARAMETER_SPECS = Dict{String, ParameterSpec}(
     "schur_mode" => ParameterSpec(:config, :port, :schur_mode),
     "deterministic_mode" => ParameterSpec(:config, :port, :deterministic_mode),
     "nd_ubfactor" => ParameterSpec(:config, :port, :nd_ubfactor),
+    "nd_nseps" => ParameterSpec(:config, :port, :nd_nseps),
     # CUDSS.jl data parameters
     "info" => ParameterSpec(:data, :solver, :none),
     "lu_nnz" => ParameterSpec(:data, :output, :none),
@@ -415,6 +418,7 @@ for (field, lo, hi, expected) in (
         (:schur_mode, 0, 1, "0 or 1"),
         (:deterministic_mode, 0, 1, "0 or 1"),
         (:nd_ubfactor, -1, typemax(Int), "-1 (library default) or an integer ≥ 0"),
+        (:nd_nseps, -1, typemax(Int), "-1 (library default) or an integer ≥ 1"),
     )
     @eval _parse_option(::Val{$(QuoteNode(field))}, value, _) =
         _parse_int($(String(field)), value, $lo, $hi, $expected)

src/symbolic/ordering.jl

@@ -8,7 +8,7 @@
 """
     ND_PROVIDER
 
-`Ref` holding `nothing` or a function `(nd_nlevels, nd_ubfactor) -> alg` that
+`Ref` holding `nothing` or a function `(nd_nlevels, nd_ubfactor, nd_nseps) -> alg` that
 returns the CliqueTrees nested-dissection algorithm. Set by the Metis extension.
 """
 const ND_PROVIDER = Ref{Any}(nothing)
@@ -25,7 +25,7 @@ function _nd_algorithm(opts::Options)
     nd_available() ||
         throw(NotSupportedError("nested dissection (reordering_alg = \"algo4\") needs Metis.jl; " *
                                 "run `using Metis` to load SparseDirectSolverMetisExt"))
-    return ND_PROVIDER[](opts.nd_nlevels, opts.nd_ubfactor)
+    return ND_PROVIDER[](opts.nd_nlevels, opts.nd_ubfactor, opts.nd_nseps)
 end
 
 """
@@ -142,7 +142,9 @@ Fill-reducing ordering of the pattern `P` (PLAN §2.3 step 2):
 * `reordering_alg = "algo5"`: natural ordering; `"algo3"`: AMD
   (`CliqueTrees.AMD()`); `"algo4"`: nested dissection with Metis (needs
   `using Metis`, otherwise [`NotSupportedError`](@ref)): `METIS_NodeND` with
-  `ufactor = nd_ubfactor` (`-1`: METIS default); `nd_nlevels` is the cuDSS
+  `ufactor = nd_ubfactor` and `nseps = nd_nseps` (`-1`: METIS defaults; on
+  fill-bound problems `nd_nseps = 4` trades a slower ordering for less fill);
+  `nd_nlevels` is the cuDSS
   *minimum* number of dissection levels, which METIS' full recursion meets on
   every graph large enough to be split that often;
   `"algo1"`/`"algo2"`: AMD on the symmetric pattern;

test/test_options.jl

@@ -15,7 +15,7 @@ const CUDSS_JL_CONFIG_PARAMETERS = ("reordering_alg", "factorization_alg", "solv
                                     "hybrid_device_memory_limit", "use_cuda_register_memory", "host_nthreads",
                                     "hybrid_execute_mode", "pivot_epsilon_alg", "nd_nlevels", "ubatch_size",
                                     "ubatch_index", "use_superpanels", "device_count", "device_indices",
-                                    "schur_mode", "deterministic_mode", "nd_ubfactor")
+                                    "schur_mode", "deterministic_mode", "nd_ubfactor", "nd_nseps")
 
 # Configuration parameters PLAN §1.3/§1.7 mark port or reinterpret:
 # name => (valid value, value of the wrong type).
@@ -41,6 +41,7 @@ const PORTED_CONFIG_VALUES = Dict{String, Tuple{Any, Any}}(
     "schur_mode" => (1, "on"),
     "deterministic_mode" => (1, 1.0),
     "nd_ubfactor" => (30, 0.3),
+    "nd_nseps" => (4, 0.25),
     # PLAN §1.7
     "pivot_sign" => (Int8[1, -1, 0, 1], [1.0, -1.0]),
     "ir_mode" => ("fgmres", :fgmres),
@@ -118,7 +119,7 @@ end
                           ("reordering_alg", "algo6"), ("reordering_alg", 6), ("factorization_alg", "algo3"),
                           ("solve_alg", "algo2"), ("matching_alg", "algo7"), ("pivot_epsilon_alg", "algo3"),
                           ("use_superpanels", 2), ("ubatch_index", -2), ("ubatch_size", -1),
-                          ("nd_ubfactor", -2), ("schedule", "fast"), ("ir_mode", "gmres"),
+                          ("nd_ubfactor", -2), ("nd_nseps", -2), ("schedule", "fast"), ("ir_mode", "gmres"),
                           ("amalgamation", (min_width = 64,)), ("amalgamation", (zero_fraction = -0.1,)),
                           ("amalgamation", (foo = 1,)), ("pivot_sign", [2, 0]),
                           ("user_schur_indices", [0, 2]), ("ubatch_mask", [1, -1]))

test/test_symbolic_etree.jl

@@ -106,8 +106,10 @@ end
             end
         end
     end
-    @test SDS.ND_PROVIDER[](10, 50).ufactor == 50
-    @test SDS.ND_PROVIDER[](0, -1).ufactor == -1
+    @test SDS.ND_PROVIDER[](10, 50, -1).ufactor == 50
+    @test SDS.ND_PROVIDER[](0, -1, -1).ufactor == -1
+    @test SDS.ND_PROVIDER[](10, -1, 4).nseps == 4
+    @test SDS.ND_PROVIDER[](10, -1, 0).nseps == -1
 
     @testset "user_perm" begin
         A = laplacian2d(7, 6)