(Semi-)Automatic Sparsity Detection

This section describes how to enable Sparsity Detection. For a detailed tutorial on how to use this in an actual problem, see this tutorial on Efficiently Solving Large Sparse Ill-Conditioned Nonlinear Systems.

Notation wise we are trying to solve for x such that nlfunc(x) = 0.

Big Table for Determining Sparsity Detection and Coloring Algorithms

f.sparsityf.jac_prototypef.colorvecSparsity DetectionColoring Algorithm
AnyNoSparsityDetector()NoColoringAlgorithm()
Not StructuredAnyNoSparsityDetector()NoColoringAlgorithm()
StructuredKnownJacobianSparsityDetector(f.jac_prototype)GreedyColoringAlgorithm(LargestFirst())
StructuredKnownJacobianSparsityDetector(f.jac_prototype)GreedyColoringAlgorithm(LargestFirst())
-----
AbstractMatrixKnownJacobianSparsityDetector(f.sparsity)ConstantColoringAlgorithm(f.colorvec)
AbstractMatrixKnownJacobianSparsityDetector(f.sparsity)GreedyColoringAlgorithm(LargestFirst())
AbstractMatrixNot StructuredKnownJacobianSparsityDetector(f.sparsity)ConstantColoringAlgorithm(f.colorvec)
AbstractMatrixNot StructuredKnownJacobianSparsityDetector(f.sparsity)GreedyColoringAlgorithm(LargestFirst())
AbstractMatrixStructuredAny🔴🔴
-----
AbstractSparsityDetectorAnyf.sparsityGreedyColoringAlgorithm(LargestFirst())
AbstractSparsityDetectorNot Structuredf.sparsityConstantColoringAlgorithm(f.colorvec)
AbstractSparsityDetectorNot Structuredf.sparsityGreedyColoringAlgorithm(LargestFirst())
AbstractSparsityDetectorStructuredKnownJacobianSparsityDetector(f.jac_prototype)ConstantColoringAlgorithm(f.colorvec)
AbstractSparsityDetectorStructuredKnownJacobianSparsityDetector(f.jac_prototype)GreedyColoringAlgorithm(LargestFirst())
  1. Structured means either a AbstractSparseMatrix or ArrayInterface.isstructured(x) is true.
  2. ❌ means not provided (default)
  3. ✅ means provided
  4. 🔴 means an error will be thrown
  5. Providing a colorvec without specifying either sparsity or jac_prototype with a sparse or structured matrix will cause us to ignore the colorvec.
  6. The function calls demonstrated above are simply pseudo-code to show the general idea.

Case I: Sparse Jacobian Prototype is Provided

Let's say you have a Sparse Jacobian Prototype jac_prototype, in this case you can create your problem as follows:

prob = NonlinearProblem(NonlinearFunction(nlfunc; jac_prototype = jac_prototype), x0)

NonlinearSolve will automatically perform matrix coloring and use sparse differentiation.

Now you can help the solver further by providing the color vector. This can be done by

prob = NonlinearProblem(
    NonlinearFunction(nlfunc; jac_prototype = jac_prototype, colorvec = colorvec), x0)

If the colorvec is not provided, then it is computed on demand.

Note

One thing to be careful about in this case is that colorvec is dependent on the autodiff backend used. ADTypes.mode(ad) isa ADTypes.ForwardMode will assume that the colorvec is the column colorvec, otherwise we will assume that the colorvec is the row colorvec.

Case II: Sparsity Detection algorithm is provided

If you don't have a Sparse Jacobian Prototype, but you know the which sparsity detection algorithm you want to use, then you can create your problem as follows:

prob = NonlinearProblem(
    NonlinearFunction(nlfunc; sparsity = SymbolicsSparsityDetector()), x0)  # Remember to have Symbolics.jl loaded
# OR
prob = NonlinearProblem(
    NonlinearFunction(nlfunc; sparsity = TracerSparsityDetector()), x0) # From SparseConnectivityTracer.jl

Refer to the documentation of DifferentiationInterface.jl and SparseConnectivityTracer.jl for more information on sparsity detection algorithms.