Fix population cleanup code in back-end (#845)

* add tests for #844

* Do not apply transform to mcounts_from_preserved_nodes in
final_population_cleanup.

Fixes #844

doc/misc/changelog.md

@@ -10,6 +10,11 @@ Bug fixes
 * Fix error initializing the founder genome of {class}`fwdpy11.DiploidPopulation`
   Issue {issue}`836`
   PR {pr}`838`
+* Fix bug in handling mutation counts when final generation was recorded as "ancient samples".
+  This bug resulted in an exception being raised.
+  Thus, previous exception-free results were not affected.
+  Issue {issue}`844`
+  PR {pr}`845`
 
 Dependencies
 

fwdpy11/src/evolve_population/evolvets.cc

@@ -80,24 +80,26 @@ final_population_cleanup(
     bool suppress_edge_table_indexing, bool preserve_selected_fixations,
     bool remove_extinct_mutations_at_finish, bool simulating_neutral_variants,
     bool reset_treeseqs_to_alive_nodes_after_simplification,
-    std::uint32_t last_preserved_generation,
-    const std::vector<std::uint32_t> &last_preserved_generation_counts,
+    std::uint32_t /*last_preserved_generation*/,
+    const std::vector<std::uint32_t> & /*last_preserved_generation_counts*/,
     fwdpy11::DiploidPopulation &pop)
 {
     index_and_count_mutations(suppress_edge_table_indexing, simulating_neutral_variants,
                               reset_treeseqs_to_alive_nodes_after_simplification, pop);
     check_mutation_table_consistency_with_count_vectors(pop, __FILE__, __LINE__);
-    if (pop.generation == last_preserved_generation
-        && !reset_treeseqs_to_alive_nodes_after_simplification
-        && !simulating_neutral_variants)
-        {
-            std::transform(begin(pop.mcounts_from_preserved_nodes),
-                           end(pop.mcounts_from_preserved_nodes),
-                           begin(last_preserved_generation_counts),
-                           begin(pop.mcounts_from_preserved_nodes),
-                           std::minus<std::uint32_t>());
-        }
-    check_mutation_table_consistency_with_count_vectors(pop, __FILE__, __LINE__);
+    // The following block was commented out in GitHub PR 845
+    // in order to fix GitHub issue 844
+    // if (pop.generation == last_preserved_generation
+    //     && !reset_treeseqs_to_alive_nodes_after_simplification
+    //     && !simulating_neutral_variants)
+    //     {
+    //         std::transform(begin(pop.mcounts_from_preserved_nodes),
+    //                        end(pop.mcounts_from_preserved_nodes),
+    //                        begin(last_preserved_generation_counts),
+    //                        begin(pop.mcounts_from_preserved_nodes),
+    //                        std::minus<std::uint32_t>());
+    //     }
+    // check_mutation_table_consistency_with_count_vectors(pop, __FILE__, __LINE__);
     if (!preserve_selected_fixations)
         {
             auto itr = std::remove_if(

tests/test_stopping_criterion.py

@@ -1,6 +1,7 @@
 import unittest
 
 import fwdpy11
+import numpy as np
 
 
 class test_stopping_criterion_DiploidPopulation(unittest.TestCase):
@@ -32,5 +33,34 @@ def generation(pop, simplified):
         self.assertEqual(self.pop.generation, 50)
 
 
+def test_stopping_criterion_with_ancient_samples():
+    """
+    This triggers github issue 844
+    """
+
+    def preserve(pop, sampler):
+        if pop.generation == 10:
+            sampler.assign(np.arange(pop.N, dtype=np.uint32))
+
+    def stopper(pop, _):
+        if pop.generation == 10:
+            return True
+        return False
+
+    pdict = {
+        "sregions": [fwdpy11.ExpS(0, 1, 1, -0.05, 1)],
+        "gvalue": fwdpy11.Multiplicative(2.0),
+        "rates": (0, 0.1, 0),
+        "simlen": 20,
+    }
+    params = fwdpy11.ModelParams(**pdict)
+    pop = fwdpy11.DiploidPopulation(100, 1.0)
+    rng = fwdpy11.GSLrng(90210)
+
+    fwdpy11.evolvets(
+        rng, pop, params, 100, recorder=preserve, stopping_criterion=stopper
+    )
+
+
 if __name__ == "__main__":
     unittest.main()

tests/test_tree_sequences.py

@@ -1084,7 +1084,7 @@ def testQtraitSim(self):
         MD = namedtuple("MD", ["g", "e", "w", "label"])
 
         class Recorder(object):
-            """ Records entire pop every 100 generations """
+            """Records entire pop every 100 generations"""
 
             def __init__(self):
                 self.data = []
@@ -1526,5 +1526,19 @@ def test_table_indexing_during_sim(test_table_indexing_during_sim_recorder):
     assert test_table_indexing_during_sim_recorder.called > 0
 
 
+def test_only_preserve_final_generation():
+    """
+    This test also triggers GitHub issue 844
+    """
+    params, rng, pop = set_up_quant_trait_model(0.1)
+
+    def preserver(pop, sampler):
+        if pop.generation == preserver.when:
+            sampler.assign(np.arange(pop.N, dtype=np.uint32))
+
+    preserver.when = params.simlen
+    fwdpy11.evolvets(rng, pop, params, 100, recorder=preserver)
+
+
 if __name__ == "__main__":
     unittest.main()