addOperatorSparsity(J, mat) implementation

Applies the same connectivity in `mat` to all pairs of variables in
`J`. Passes unit tests.

Author: bendudson

include/bout/petsc_jacobian.hxx

@@ -27,6 +27,13 @@
 /// @param in_var  Column variable index in [0, nvars).
 void addOperatorSparsity(Mat Jfd, Mat sub, int out_var, int in_var);
 
+/// @brief Insert the nonzero pattern of @p sub into every variable block of
+///        @p Jfd.
+///
+/// Equivalent to calling addOperatorSparsity(Jfd, sub, out_var, in_var) for
+/// every combination of @p out_var and @p in_var in [0, nvars), where
+/// @c nvars is inferred as @c Jfd_global / @c sub_global.
+///
 /// @param Jfd     The Jacobian matrix to populate. Must be preallocated.
 /// @param sub     Evolving-cell submatrix providing the nonzero pattern.
 void addOperatorSparsity(Mat Jfd, Mat sub);

src/mesh/petsc_jacobian.cxx

@@ -45,6 +45,21 @@ void addOperatorSparsity(Mat Jfd, Mat sub, int out_var, int in_var) {
   }
 }
 
-void addOperatorSparsity(Mat Jfd, Mat sub) {}
+void addOperatorSparsity(Mat Jfd, Mat sub) {
+  PetscInt jfd_global{0}, sub_global{0};
+  MatGetSize(Jfd, &jfd_global, nullptr);
+  MatGetSize(sub, &sub_global, nullptr);
+
+  ASSERT1(sub_global > 0);
+  ASSERT1(jfd_global % sub_global == 0);
+
+  const int nvars = static_cast<int>(jfd_global / sub_global);
+
+  for (int out_var = 0; out_var < nvars; ++out_var) {
+    for (int in_var = 0; in_var < nvars; ++in_var) {
+      addOperatorSparsity(Jfd, sub, out_var, in_var);
+    }
+  }
+}
 
 #endif // BOUT_HAS_PETSC