richardposner · wshlavacek · May 23, 2026 · May 23, 2026
diff --git a/CHANGELOG.md b/CHANGELOG.md
@@ -5,6 +5,28 @@ All notable changes to this project are documented in this file.
 The format is based on [Keep a Changelog](https://keepachangelog.com/en/1.1.0/),
 and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0.html).
 
+## [3.2.1] — 2026-05-23
+
+### Fixed
+
+- **Exact-species pattern methods are now component-order-insensitive.**
+  The `#9` session methods `get_species_count` / `remove_species` /
+  `set_species_count` matched a species only when its components were
+  written in RuleMonkey's own canonical order; a semantically identical
+  pattern with the components swapped (e.g. `X(p~0,y)` for the canonical
+  `X(y,p~0)`) silently matched nothing and returned 0. The `add_species`
+  path, by contrast, already placed components by molecule-type
+  declaration index, so it *did* canonicalize — leaving
+  `set_species_count` with a non-canonical pattern diffing against a
+  wrong baseline of 0 and overshooting its target. The fix routes the
+  match path through the same declaration-order placement
+  (`pattern_to_complex_graph` now orders each molecule's components by
+  `comp_type_index` and remaps bond endpoints accordingly), matching
+  `extract_complex` and NFsim's order-insensitive exact-species matcher.
+  Closes [#13](https://github.com/richardposner/RuleMonkey/issues/13);
+  follow-up to
+  [#9](https://github.com/richardposner/RuleMonkey/issues/9).
+
 ## [3.2.0] — 2026-05-18
 
 ### Added

diff --git a/CMakeLists.txt b/CMakeLists.txt
@@ -1,6 +1,6 @@
 cmake_minimum_required(VERSION 3.20)
 
-project(RuleMonkey VERSION 3.2.0 LANGUAGES CXX)
+project(RuleMonkey VERSION 3.2.1 LANGUAGES CXX)
 
 if(CMAKE_VERSION VERSION_LESS 3.21)
     if(CMAKE_SOURCE_DIR STREQUAL PROJECT_SOURCE_DIR)

diff --git a/cpp/rulemonkey/engine.cpp b/cpp/rulemonkey/engine.cpp
@@ -7774,29 +7774,43 @@ canonical::ComplexGraph extract_complex(const AgentPool& pool, const Model& mode
 // Build a ComplexGraph straight from a parsed species Pattern (issue #9
 // §1) — the parse-side analogue of extract_complex.  `pat` is an exact,
 // fully-specified species from parse_species_pattern, so each pattern
-// molecule carries every component (in listed order) with a concrete
-// state, and pat.bonds gives the edges by flat component index.  The
-// graph is isomorphic to the one extract_complex builds for the same
-// physical species, so canonical_label yields the same string — that
-// is what makes pattern-keyed species_count a byte-equal lookup.
+// molecule carries every component (once) with a concrete state, and
+// pat.bonds gives the edges by flat component index.  The graph is
+// isomorphic to the one extract_complex builds for the same physical
+// species, so canonical_label yields the same string — that is what
+// makes pattern-keyed species_count a byte-equal lookup.
+//
+// Components are placed in molecule-type DECLARATION order (by
+// comp_type_index), exactly as extract_complex feeds the pool's
+// components.  canonical_label preserves the input's component-name
+// slot layout (canonical.cpp §3.3 "component-order contract"), so a
+// pattern listed in a non-canonical component order — e.g. "X(p~0,y)"
+// for the canonical "X(y,p~0)" — must be reordered here, or it yields a
+// different label and silently matches nothing while the add path
+// (instantiate_pattern_complex, which already keys by comp_type_index)
+// does match.  That divergence was issue #13.
 canonical::ComplexGraph pattern_to_complex_graph(const Pattern& pat) {
   canonical::ComplexGraph g;
   for (const auto& pm : pat.molecules) {
-    std::vector<std::pair<std::string, std::string>> comps;
-    comps.reserve(pm.components.size());
+    // Slot i carries the molecule type's i-th declared component.  The
+    // pattern is fully specified (every type component listed exactly
+    // once), so comp_type_index is a bijection onto [0, n) and every
+    // slot is filled.
+    std::vector<std::pair<std::string, std::string>> comps(pm.components.size());
     for (const auto& pc : pm.components)
-      comps.emplace_back(pc.name, pc.required_state); // "" for a stateless component
+      comps[pc.comp_type_index] = {pc.name, pc.required_state}; // "" if stateless
     g.add_molecule(pm.type_name, comps);
   }
-  // pat.bonds endpoints are flat component indices; map each back to its
-  // (molecule index, local component index) — components were fed to
-  // add_molecule in the same listed order pat.flat_index() counts in.
+  // pat.bonds endpoints are flat indices in the pattern's *listed* order;
+  // map each back to (molecule index, declaration-order local index) so
+  // the bond lands on the reordered component above.
   const auto flat_to_loc = [&pat](int flat) -> std::pair<int, int> {
     int acc = 0;
     for (int mi = 0; mi < static_cast<int>(pat.molecules.size()); ++mi) {
-      const int n = static_cast<int>(pat.molecules[mi].components.size());
+      const auto& mol = pat.molecules[mi];
+      const int n = static_cast<int>(mol.components.size());
       if (flat < acc + n)
-        return {mi, flat - acc};
+        return {mi, mol.components[flat - acc].comp_type_index};
       acc += n;
     }
     throw std::runtime_error("pattern_to_complex_graph: flat-index overflow");

diff --git a/pyproject.toml b/pyproject.toml
@@ -1,6 +1,6 @@
 [project]
 name = "rulemonkey-harness"
-version = "3.2.0"
+version = "3.2.1"
 description = "Test and benchmarking harness for the RuleMonkey C++ engine"
 readme = "README.md"
 requires-python = ">=3.10"

diff --git a/tests/CMakeLists.txt b/tests/CMakeLists.txt
@@ -134,6 +134,7 @@ add_test(
     NAME species_methods_test
     COMMAND species_methods_test
         ${_rulemonkey_test_dir}/reference/nfsim/xml/A_plus_A.xml
+        ${_rulemonkey_test_dir}/cpp/species_order_model.xml
 )
 
 # Species enumeration + `.species` output (issue #9 §2): exercises the

diff --git a/tests/cpp/species_methods_test.cpp b/tests/cpp/species_methods_test.cpp
@@ -126,6 +126,45 @@ void test_ssa_runs_after_mutation(rulemonkey::RuleMonkeySimulator& sim) {
   check(sim.get_species_count(kMonomer) >= 0, "get_species_count is callable post-simulate");
 }
 
+// Issue #13: get/match must be component-order-insensitive, consistent
+// with the add path.  The species_order model declares X with a
+// stateful `p~0~1` and a stateless `y`, seeded with 5000 X(p~0,y).  An
+// exact pattern written with the components in either order denotes the
+// same species, so get_species_count must return 5000 for both — and
+// set_species_count must reach an exact target regardless of the order
+// the caller wrote the pattern in (the silent wrong-baseline bug).
+void test_component_order_insensitive(const std::string& order_xml) {
+  rulemonkey::RuleMonkeySimulator sim(order_xml);
+  sim.initialize(/*seed=*/11);
+
+  const int n = sim.get_species_count("X(p~0,y)");
+  check(n == 5000, "X(p~0,y) seeds 5000 copies");
+  check(sim.get_species_count("X(y,p~0)") == n,
+        "get_species_count is component-order-insensitive: X(y,p~0) == X(p~0,y)");
+
+  // Mutate the *other* species (X(p~1,y)) through both orders and
+  // confirm the match path tracks it either way.
+  sim.add_species("X(y,p~1)", 25); // non-canonical order on the add path
+  check(sim.get_species_count("X(p~1,y)") == 25,
+        "add_species with swapped components creates the species the canonical order matches");
+  check(sim.get_species_count("X(y,p~1)") == 25, "...and the swapped-order query agrees");
+
+  // set_species_count with a NON-canonical order must reach exactly N,
+  // not diff against a wrong-order baseline of 0 and overshoot.
+  sim.set_species_count("X(p~1,y)", 100); // canonical order
+  check(sim.get_species_count("X(y,p~1)") == 100, "set to 100 via canonical order lands on 100");
+  sim.set_species_count("X(y,p~1)", 40); // non-canonical order, drives down
+  check(sim.get_species_count("X(p~1,y)") == 40,
+        "set_species_count with swapped components reaches the exact target");
+
+  // remove_species via the swapped order targets the right species.
+  sim.remove_species("X(y,p~1)", 40);
+  check(sim.get_species_count("X(p~1,y)") == 0,
+        "remove_species with swapped components removes the matched species");
+
+  sim.destroy_session();
+}
+
 void test_error_surface(const std::string& xml) {
   rulemonkey::RuleMonkeySimulator sim(xml);
 
@@ -154,11 +193,12 @@ void test_error_surface(const std::string& xml) {
 } // namespace
 
 int main(int argc, char* argv[]) {
-  if (argc < 2) {
-    std::fprintf(stderr, "Usage: species_methods_test <A_plus_A.xml>\n");
+  if (argc < 3) {
+    std::fprintf(stderr, "Usage: species_methods_test <A_plus_A.xml> <species_order_model.xml>\n");
     return 2;
   }
   const std::string xml = argv[1];
+  const std::string order_xml = argv[2];
 
   try {
     rulemonkey::RuleMonkeySimulator sim(xml);
@@ -171,6 +211,7 @@ int main(int argc, char* argv[]) {
     sim.destroy_session();
 
     test_error_surface(xml);
+    test_component_order_insensitive(order_xml);
   } catch (const std::exception& e) {
     std::fprintf(stderr, "EXCEPTION: %s\n", e.what());
     return 2;

diff --git a/tests/cpp/species_order_model.xml b/tests/cpp/species_order_model.xml
@@ -0,0 +1,98 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<!-- Hand-authored regression model for issue #13: a molecule X with a
+     stateful component p~0~1 and a stateless component y, so an exact
+     species pattern can be written in two component orders that denote
+     the same species ("X(p~0,y)" vs "X(y,p~0)").  Seeded with X(p~0,y);
+     a single state-change rule X(p~0) -> X(p~1) keeps the model loadable
+     and lets add/set produce the X(p~1,y) species too. -->
+<sbml xmlns="http://www.sbml.org/sbml/level3" level="3" version="1">
+  <model id="species_order">
+    <ListOfParameters>
+      <Parameter id="X_tot" type="Constant" value="5000" expr="5000"/>
+      <Parameter id="k" type="Constant" value="0.1" expr="0.1"/>
+    </ListOfParameters>
+    <ListOfMoleculeTypes>
+      <MoleculeType id="X">
+        <ListOfComponentTypes>
+          <ComponentType id="p">
+            <ListOfAllowedStates>
+              <AllowedState id="0"/>
+              <AllowedState id="1"/>
+            </ListOfAllowedStates>
+          </ComponentType>
+          <ComponentType id="y"/>
+        </ListOfComponentTypes>
+      </MoleculeType>
+    </ListOfMoleculeTypes>
+    <ListOfCompartments>
+    </ListOfCompartments>
+    <ListOfSpecies>
+      <Species id="S1" concentration="X_tot" name="X(p~0,y)">
+        <ListOfMolecules>
+          <Molecule id="S1_M1" name="X">
+            <ListOfComponents>
+              <Component id="S1_M1_C1" name="p" state="0" numberOfBonds="0"/>
+              <Component id="S1_M1_C2" name="y" numberOfBonds="0"/>
+            </ListOfComponents>
+          </Molecule>
+        </ListOfMolecules>
+      </Species>
+    </ListOfSpecies>
+    <ListOfReactionRules>
+      <ReactionRule id="RR1" name="_R1" symmetry_factor="1">
+        <ListOfReactantPatterns>
+          <ReactantPattern id="RR1_RP1">
+            <ListOfMolecules>
+              <Molecule id="RR1_RP1_M1" name="X">
+                <ListOfComponents>
+                  <Component id="RR1_RP1_M1_C1" name="p" state="0" numberOfBonds="0"/>
+                </ListOfComponents>
+              </Molecule>
+            </ListOfMolecules>
+          </ReactantPattern>
+        </ListOfReactantPatterns>
+        <ListOfProductPatterns>
+          <ProductPattern id="RR1_PP1">
+            <ListOfMolecules>
+              <Molecule id="RR1_PP1_M1" name="X">
+                <ListOfComponents>
+                  <Component id="RR1_PP1_M1_C1" name="p" state="1" numberOfBonds="0"/>
+                </ListOfComponents>
+              </Molecule>
+            </ListOfMolecules>
+          </ProductPattern>
+        </ListOfProductPatterns>
+        <RateLaw id="RR1_RateLaw" type="Ele" totalrate="0">
+          <ListOfRateConstants>
+            <RateConstant value="k"/>
+          </ListOfRateConstants>
+        </RateLaw>
+        <Map>
+          <MapItem sourceID="RR1_RP1_M1" targetID="RR1_PP1_M1"/>
+          <MapItem sourceID="RR1_RP1_M1_C1" targetID="RR1_PP1_M1_C1"/>
+        </Map>
+        <ListOfOperations>
+          <StateChange site="RR1_RP1_M1_C1" finalState="1"/>
+        </ListOfOperations>
+      </ReactionRule>
+    </ListOfReactionRules>
+    <ListOfObservables>
+      <Observable id="O1" name="X0" type="Species">
+        <ListOfPatterns>
+          <Pattern id="O1_P1" matchOnce="1">
+            <ListOfMolecules>
+              <Molecule id="O1_P1_M1" name="X">
+                <ListOfComponents>
+                  <Component id="O1_P1_M1_C1" name="p" state="0" numberOfBonds="0"/>
+                  <Component id="O1_P1_M1_C2" name="y" numberOfBonds="0"/>
+                </ListOfComponents>
+              </Molecule>
+            </ListOfMolecules>
+          </Pattern>
+        </ListOfPatterns>
+      </Observable>
+    </ListOfObservables>
+    <ListOfFunctions>
+    </ListOfFunctions>
+  </model>
+</sbml>