Fix SMP serial-unstable filters: SurfaceNets3D BoundaryLabels + LengthDistribution RNG (#176)

Co-authored-by: Claude Opus 4.8 <noreply@anthropic.com>

Author: akaszynski

Common/Math/vtkReservoirSampler.h

@@ -42,6 +42,21 @@ template <typename Integer, bool Monotonic = true>
 class vtkReservoirSampler : public vtkReservoirSamplerBase
 {
 public:
+  vtkReservoirSampler() = default;
+
+  /// Construct a sampler with a fixed seed for reproducible sampling.
+  ///
+  /// By default a sampler draws its seed from std::random_device, so results
+  /// vary run-to-run (and, when the same instance is shared across threads,
+  /// from thread to thread). Passing an explicit seed makes every draw this
+  /// instance performs deterministic, which callers that need reproducible
+  /// output (and byte-stable serial-vs-parallel behavior) should prefer.
+  explicit vtkReservoirSampler(SeedType seed)
+    : Seed(seed)
+    , HasSeed(true)
+  {
+  }
+
   /// Choose kk items from a sequence of (0, nn - 1).
   ///
   /// This will throw an exception if kk <= 0.
@@ -98,7 +113,7 @@ class vtkReservoirSampler : public vtkReservoirSamplerBase
       return;
     }
 
-    std::mt19937 generator(vtkReservoirSampler::RandomSeed());
+    std::mt19937 generator(this->HasSeed ? this->Seed : vtkReservoirSampler::RandomSeed());
     std::uniform_real_distribution<> unitUniform(0., 1.);
     std::uniform_int_distribution<Integer> randomIndex(0, kk - 1);
     double w = exp(log(unitUniform(generator)) / kk);
@@ -136,6 +151,9 @@ class vtkReservoirSampler : public vtkReservoirSamplerBase
       std::sort(data.begin(), data.end());
     }
   }
+
+  SeedType Seed = 0;
+  bool HasSeed = false;
 };
 
 VTK_ABI_NAMESPACE_END

Filters/Core/vtkSurfaceNets3D.cxx

@@ -20,6 +20,7 @@
 #include "vtkStreamingDemandDrivenPipeline.h"
 #include "vtkTriangle.h"
 
+#include <algorithm>
 #include <cstdint>
 #include <memory>
 
@@ -1411,7 +1412,16 @@ void SurfaceNets<TArray>::ConfigureOutput(
 
     // Scalars, which are of type T and 2-components
     newScalars->SetNumberOfTuples(outputEMD.NumQuads);
-    this->NewScalars = vtk::DataArrayValueRange<2>(newScalars).begin();
+    // The quad count is a padded upper bound: an isolated boundary case can
+    // reserve a scalar tuple that GenerateQuads never writes (the connectivity
+    // array zero-fills on ResizeExact, but SetNumberOfTuples leaves the scalar
+    // storage uninitialized). Those never-written BoundaryLabels tuples read
+    // back as garbage -- run-to-run nondeterministic under threading, and even
+    // serially unstable on a dirtied heap. Zero-init so any padding slot is
+    // deterministic; every genuinely-emitted quad overwrites its own tuple.
+    auto scalarRange = vtk::DataArrayValueRange<2>(newScalars);
+    std::fill(scalarRange.begin(), scalarRange.end(), 0);
+    this->NewScalars = scalarRange.begin();
 
     // Smoothing stencils, which are represented by a vtkCellArray
     stencils->ResizeExact(outputEMD.NumPoints, outputEMD.NumStencilEdges);

Filters/Statistics/vtkLengthDistribution.cxx

@@ -22,6 +22,7 @@ void vtkLengthDistribution::PrintSelf(ostream& os, vtkIndent indent)
   this->Superclass::PrintSelf(os, indent);
   os << indent << "SampleSize: " << this->SampleSize << "\n";
   os << indent << "SortSample: " << (this->SortSample ? "T" : "F") << "\n";
+  os << indent << "Seed: " << this->Seed << "\n";
 }
 
 double vtkLengthDistribution::GetLengthQuantile(double qq)
@@ -81,7 +82,10 @@ int vtkLengthDistribution::RequestData(
   vtkNew<vtkIdList> dummyIds;
   dataIn->GetCellPoints(0, dummyIds.GetPointer());
 
-  vtkReservoirSampler<vtkIdType> sampler;
+  // Seed the sampler so both the cell subset and the per-cell connectivity
+  // picks are deterministic -- otherwise vtkReservoirSampler draws from
+  // std::random_device and each thread (and each run) samples differently.
+  vtkReservoirSampler<vtkIdType> sampler(this->Seed);
   std::vector<vtkIdType> ids = sampler(numSamples, dataIn->GetNumberOfCells());
   vtkSMPTools::For(0, static_cast<vtkIdType>(ids.size()),
     [&dataIn, &lengths, &sampler, &ids](vtkIdType begin, vtkIdType end)

Filters/Statistics/vtkLengthDistribution.h

@@ -54,6 +54,14 @@ class VTKFILTERSSTATISTICS_EXPORT vtkLengthDistribution : public vtkTableAlgorit
   vtkSetMacro(SortSample, bool);
   vtkBooleanMacro(SortSample, bool);
 
+  /// Set/get the seed used to randomly sample cells and connectivity entries.
+  ///
+  /// The sampling is seeded so the output is reproducible from run to run
+  /// (and byte-identical whether the sampling loop runs serially or threaded).
+  /// Change the seed to obtain a different random subset. The default is 0.
+  vtkGetMacro(Seed, unsigned int);
+  vtkSetMacro(Seed, unsigned int);
+
   /// Return the length scale at a particular quantile.
   ///
   /// This method must only be invoked after the filter
@@ -77,6 +85,7 @@ class VTKFILTERSSTATISTICS_EXPORT vtkLengthDistribution : public vtkTableAlgorit
 
   vtkIdType SampleSize = 100000;
   bool SortSample = true;
+  unsigned int Seed = 0;
 
 private:
   vtkLengthDistribution(const vtkLengthDistribution&) = delete;

cvista-validate/smpParityCases.cxx

@@ -341,7 +341,7 @@ std::vector<Case> RegisterCases()
   };
   // Tag the just-registered case as a known, documented nondeterminism exception
   // (see Case::knownIssue): reported with its magnitude but not a gate failure.
-  auto markKnown = [&cases](const std::string& reason) {
+  [[maybe_unused]] auto markKnown = [&cases](const std::string& reason) {
     cases.back().knownIssue = true;
     cases.back().knownIssueReason = reason;
   };
@@ -523,22 +523,17 @@ std::vector<Case> RegisterCases()
     f->SetInputData(in.ugrid);
     return vtkSmartPointer<vtkAlgorithm>(f);
   });
-  add(
-    "vtkLengthDistribution", "Filters/Statistics", Risk::Reduce,
-    [](const Inputs& in) {
-      vtkNew<vtkLengthDistribution> f;
-      f->SetInputData(in.ugrid);
-      return vtkSmartPointer<vtkAlgorithm>(f);
-    },
-    /*orderRelaxed=*/true);
-  // NOT a threading defect: vtkLengthDistribution draws its samples from a
-  // vtkReservoirSampler seeded from std::random_device on every call
-  // (Common/Math/vtkReservoirSampler.cxx), and for each sampled cell it randomly
-  // picks which 2 vertices to measure. Its output is therefore nondeterministic
-  // even SINGLE-THREADED -- the serial reference is itself random -- so a
-  // parallel-vs-serial gate cannot apply. Excluded from the pass/fail count.
-  markKnown("unseeded RNG (vtkReservoirSampler <- std::random_device); "
-            "nondeterministic even single-threaded, not a threading defect");
+  add("vtkLengthDistribution", "Filters/Statistics", Risk::Reduce, [](const Inputs& in) {
+    vtkNew<vtkLengthDistribution> f;
+    f->SetInputData(in.ugrid);
+    return vtkSmartPointer<vtkAlgorithm>(f);
+  });
+  // Now GATED byte-exact: vtkLengthDistribution used to seed its
+  // vtkReservoirSampler from std::random_device on every call, making the
+  // output nondeterministic even single-threaded. The sampler is now seedable
+  // and the filter seeds it (default Seed=0), so both the sampled cell subset
+  // and the per-cell vertex picks are deterministic; each sample writes its own
+  // table row, so serial and threaded output are byte-identical.
   add("vtkSelectEnclosedPoints", "Filters/Modeling", Risk::Reduce, [](const Inputs& in) {
     vtkNew<vtkSelectEnclosedPoints> f;
     f->SetInputData(in.cloud);
@@ -650,15 +645,16 @@ std::vector<Case> RegisterCases()
     f->SetClipFunction(plane(1, 1, 0));
     return vtkSmartPointer<vtkAlgorithm>(f);
   });
-  // Expected byte-exact: the GEOMETRY of the default OUTPUT_STYLE_BOUNDARY path is
+  // GATED byte-exact: the GEOMETRY of the default OUTPUT_STYLE_BOUNDARY path is
   // deterministic (integer voxel-center coords, serial prefix-sum offsets,
-  // per-thread label lookup, Jacobi double-buffered smoothing) -- and indeed the
-  // validator measures points + connectivity byte-identical. HOWEVER the
-  // "BoundaryLabels" cell-data comes out as garbage (~9e8, far outside the 0..3
-  // label range) that changes every run -- an uninitialized/wild read that is
-  // nondeterministic even SINGLE-THREADED. The serial-vs-serial pre-check catches
-  // this and reports it as SERIAL-UNSTABLE (a real defect, but orthogonal to
-  // threading), so it is not charged against the parallel-vs-serial gate.
+  // per-thread label lookup, Jacobi double-buffered smoothing) -- points +
+  // connectivity are byte-identical. The "BoundaryLabels" cell-data used to come
+  // out as garbage (~9e8, far outside the 0..3 label range) that changed every
+  // run: the quad count is a padded upper bound, so an isolated boundary case
+  // could reserve a scalar tuple GenerateQuads never wrote, and SetNumberOfTuples
+  // left that storage uninitialized. The scalar buffer is now zero-initialized at
+  // allocation, so any padding tuple is deterministic and every emitted quad
+  // overwrites its own -- BoundaryLabels is now byte-identical serial-vs-parallel.
   add("vtkSurfaceNets3D", "Filters/Core", Risk::Iso, [](const Inputs& in) {
     vtkNew<vtkSurfaceNets3D> f;
     f->SetInputData(in.labelImage);