Merge branch 'master' into multitrait

Author: bhaller

EidosScribe/EidosHelpFunctions.rtf

@@ -6686,6 +6686,24 @@ The task keys recognized, and the tasks they govern, are:\
 \f3\fs20 \cf2 Typically, a dictionary of task keys and thread counts is read from a file and set up with this function at initialization time, but it is also possible to change new task thread counts dynamically.  If Eidos is not configured to run multithreaded, this function has no effect.\
 \pard\pardeftab397\li720\fi-446\ri720\sb180\sa60\partightenfactor0
 
+\f1\fs18 \cf2 (string$)readLine(void)\
+\pard\pardeftab397\li547\ri720\sb60\sa60\partightenfactor0
+
+\f0\b\fs20 \cf2 Reads a line of input from the \'93standard input\'94 file.
+\f3\b0   This function is intended to allow for communication between a running Eidos script (such as a SLiM simulation) and an external process that is sending input to it.  It reads one line from the standard input and returns it as a singleton 
+\f1\fs18 string
+\f3\fs20 .  This is done with the C++ function 
+\f1\fs18 std::getline()
+\f3\fs20 , using the 
+\f1\fs18 std::cin
+\f3\fs20  file.  If that call returns false (typically because the end-of-file was reached or there was a read error of some kind), the empty string 
+\f1\fs18 ""
+\f3\fs20  will be returned.  The 
+\f1\fs18 system()
+\f3\fs20  function can also be helpful for related situations.\
+This function will raise an error if called in a GUI environment such as EidosScribe, SLiMgui, or SLiMguiLegacy, since there is no standard input set up for the running script in those environments.\
+\pard\pardeftab397\li720\fi-446\ri720\sb180\sa60\partightenfactor0
+
 \f1\fs18 \cf2 (void)rm([Ns\'a0variableNames\'a0=\'a0NULL])\
 \pard\pardeftab397\li547\ri720\sb60\sa60\partightenfactor0
 

QtSLiM/help/EidosHelpFunctions.html

@@ -648,6 +648,9 @@
 "UNIQUE_MUTRUNS"<span class="Apple-tab-span">	</span></span>uniquing mutation runs (internal bookkeeping)<span class="s2"><br>
 "SURVIVAL"<span class="Apple-tab-span">	</span></span>survival evaluation (no callbacks)</p>
 <p class="p5">Typically, a dictionary of task keys and thread counts is read from a file and set up with this function at initialization time, but it is also possible to change new task thread counts dynamically.<span class="Apple-converted-space">  </span>If Eidos is not configured to run multithreaded, this function has no effect.</p>
+<p class="p4">(string$)readLine(void)</p>
+<p class="p5"><b>Reads a line of input from the “standard input” file.</b><span class="Apple-converted-space">  </span>This function is intended to allow for communication between a running Eidos script (such as a SLiM simulation) and an external process that is sending input to it.<span class="Apple-converted-space">  </span>It reads one line from the standard input and returns it as a singleton <span class="s2">string</span>.<span class="Apple-converted-space">  </span>This is done with the C++ function <span class="s2">std::getline()</span>, using the <span class="s2">std::cin</span> file.<span class="Apple-converted-space">  </span>If that call returns false (typically because the end-of-file was reached or there was a read error of some kind), the empty string <span class="s2">""</span> will be returned.<span class="Apple-converted-space">  </span>The <span class="s2">system()</span> function can also be helpful for related situations.</p>
+<p class="p5">This function will raise an error if called in a GUI environment such as EidosScribe, SLiMgui, or SLiMguiLegacy, since there is no standard input set up for the running script in those environments.</p>
 <p class="p4">(void)rm([Ns variableNames = NULL])</p>
 <p class="p5"><b>Removes variables</b> from the Eidos namespace; in other words, it causes the variables to become undefined.<span class="Apple-converted-space">  </span>Variables are specified by their <span class="s2">string</span> name in the <span class="s2">variableNames</span> parameter.<span class="Apple-converted-space">  </span>If the optional <span class="s2">variableNames</span> parameter is <span class="s2">NULL</span> (the default), <i>all</i> variables will be removed (be careful!).</p>
 <p class="p5">In SLiM 3, there was an optional parameter <span class="s2">removeConstants</span> that, if <span class="s2">T</span>, allowed you to remove defined constants (and then potentially redefine them to have a different value).<span class="Apple-converted-space">  </span>The <span class="s2">removeConstants</span> parameter was removed in SLiM 4, since the <span class="s2">defineGlobal()</span> function now provides the ability to define (and redefine) global variables that are not constant.</p>

VERSIONS

@@ -18,7 +18,8 @@ development head (in the master branch):
 	add the PCG random number generator, switch to pcg32_fast and pcg64_fast, remove all use of the old taus2 and MT19937-64 generators; note this completely breaks backward reproducibility
 	fix various bugs involving conflicts between defined constants and other symbols, including #573 and #574; this sets new definition rules that could break some existing scripts (but is unlikely to)
 	fix #575, QtSLiM terminates early when single-stepping with a rescheduled script block
-	----- last merge into multitrait was done here
+	fix #579, crash in models with (a) tree-seq recording, (b) multiple chromosomes, AND (c) rejection of proposed offspring in modifyChild()
+	add readLine() function for obtaining stdin input to running models, thanks to Chris Talbot
 
 
 multitrait branch:

core/species.cpp

@@ -5842,8 +5842,62 @@ void Species::RetractNewIndividual()
 	// around the code since it seems to keep coming back...
 	//current_new_individual_ = nullptr;
 
-	for (TreeSeqInfo &tsinfo : treeseq_)
-		tsk_table_collection_truncate(&tsinfo.tables_, &tsinfo.table_position_);
+	size_t trees_count = treeseq_.size();
+	
+	if (trees_count > 0)
+	{
+		TreeSeqInfo &tsinfo_0 = treeseq_[0];
+		
+		// BCH 12/1/2025: The base table collection can restore its bookmarked position directly;
+		// that will reset the bookmarked positions in all of the shared tables as well.
+		tsk_table_collection_truncate(&tsinfo_0.tables_, &tsinfo_0.table_position_);
+		
+		// BCH 12/1/2025: In the multichrom case we need to protect against a segfault inside
+		// tsk_table_collection_truncate() for the secondary table collections.  This is because
+		// they have NULL for their various column pointers, and tsk_table_collection_truncate()
+		// accesses index 0 of every offset column to get the offset for row 0.  (It is always
+		// for row 0 in the shared tables because they are zeroed out; their num_rows was zero
+		// in RecordTablePosition().)  See https://github.com/MesserLab/SLiM/issues/579 for details.
+		// BEWARE: This code will need updating if new shared tables are added, or new columns
+		// are added within the existing shared table.  Any offset column that is accessed in the
+		// ..._truncate() functions for the shared tables needs to be protected here.
+		tsk_size_t zero_value = 0;
+		tsk_size_t *pointer_to_zero_value = &zero_value;
+		
+		for (size_t trees_index = 1; trees_index < trees_count; ++trees_index)
+		{
+			TreeSeqInfo &tsinfo_i = treeseq_[trees_index];
+			
+#if DEBUG
+			// This protection scheme relies upon the bookmarked row being zero for shared tables;
+			// only the zeroth element of each offset column is set up by the hack here.
+			if ((tsinfo_i.table_position_.nodes != 0) ||
+				(tsinfo_i.table_position_.individuals != 0) ||
+				(tsinfo_i.table_position_.populations != 0))
+				EIDOS_TERMINATION << "ERROR (Species::RetractNewIndividual): (internal error) tree sequence bookmark for a shared table in a secondary table collection is non-zero." << EidosTerminate();
+#endif
+			
+			tsinfo_i.tables_.nodes.metadata_offset = pointer_to_zero_value;
+			tsinfo_i.tables_.individuals.location_offset = pointer_to_zero_value;
+			tsinfo_i.tables_.individuals.parents_offset = pointer_to_zero_value;
+			tsinfo_i.tables_.individuals.metadata_offset = pointer_to_zero_value;
+			tsinfo_i.tables_.populations.metadata_offset = pointer_to_zero_value;
+			
+			tsk_table_collection_truncate(&tsinfo_i.tables_, &tsinfo_i.table_position_);
+			
+			tsinfo_i.tables_.nodes.metadata_offset = NULL;
+			tsinfo_i.tables_.individuals.location_offset = NULL;
+			tsinfo_i.tables_.individuals.parents_offset = NULL;
+			tsinfo_i.tables_.individuals.metadata_offset = NULL;
+			tsinfo_i.tables_.populations.metadata_offset = NULL;
+		}
+	}
+	
+	// The above code boils down to this, which was the old code before #579 came along.
+	// It is the same apart from the complicated protection against segfaults:
+	//
+	//	for (TreeSeqInfo &tsinfo : treeseq_)
+	//		tsk_table_collection_truncate(&tsinfo.tables_, &tsinfo.table_position_);
 }
 
 void Species::RecordNewHaplosome(slim_position_t *p_breakpoints, int p_breakpoints_count, Haplosome *p_new_haplosome, 
@@ -7585,7 +7639,7 @@ void Species::WriteTreeSequence(std::string &p_recording_tree_path, bool p_simpl
 		// Add in the mutation.parent information; valid tree sequences need parents, but we don't keep them while running
 		ret = tsk_table_collection_build_index(&output_tables, 0);
 		if (ret < 0) handle_error("tsk_table_collection_build_index", ret);
-		ret = tsk_table_collection_compute_mutation_parents(&output_tables, 0);
+		ret = tsk_table_collection_compute_mutation_parents(&output_tables, TSK_NO_CHECK_INTEGRITY);
 		if (ret < 0) handle_error("tsk_table_collection_compute_mutation_parents", ret);
 
 		{
@@ -7987,7 +8041,7 @@ void Species::CrosscheckTreeSeqIntegrity(void)
 				ret = tsk_table_collection_build_index(tables_copy, 0);
 				if (ret < 0) handle_error("tsk_table_collection_build_index", ret);
 
-				ret = tsk_table_collection_compute_mutation_parents(tables_copy, 0);
+				ret = tsk_table_collection_compute_mutation_parents(tables_copy, TSK_NO_CHECK_INTEGRITY);
 				if (ret < 0) handle_error("tsk_table_collection_compute_mutation_parents", ret);
 
 			}

eidos/eidos_functions.cpp

@@ -369,6 +369,7 @@ const std::vector<EidosFunctionSignature_CSP> &EidosInterpreter::BuiltInFunction
 		signatures->emplace_back((EidosFunctionSignature *)(new EidosFunctionSignature("fileExists",		Eidos_ExecuteFunction_fileExists,	kEidosValueMaskLogical | kEidosValueMaskSingleton))->AddString_S(gEidosStr_filePath));
 		signatures->emplace_back((EidosFunctionSignature *)(new EidosFunctionSignature("flushFile",			Eidos_ExecuteFunction_flushFile,	kEidosValueMaskLogical | kEidosValueMaskSingleton))->AddString_S(gEidosStr_filePath));
 		signatures->emplace_back((EidosFunctionSignature *)(new EidosFunctionSignature("readFile",			Eidos_ExecuteFunction_readFile,		kEidosValueMaskString))->AddString_S(gEidosStr_filePath));
+		signatures->emplace_back((EidosFunctionSignature *)(new EidosFunctionSignature("readLine",			Eidos_ExecuteFunction_readLine,		kEidosValueMaskString | kEidosValueMaskSingleton)));
 		signatures->emplace_back((EidosFunctionSignature *)(new EidosFunctionSignature("setwd",				Eidos_ExecuteFunction_setwd,		kEidosValueMaskString | kEidosValueMaskSingleton))->AddString_S("path"));
 		signatures->emplace_back((EidosFunctionSignature *)(new EidosFunctionSignature("tempdir",			Eidos_ExecuteFunction_tempdir,		kEidosValueMaskString | kEidosValueMaskSingleton)));
 		signatures->emplace_back((EidosFunctionSignature *)(new EidosFunctionSignature("writeFile",			Eidos_ExecuteFunction_writeFile,	kEidosValueMaskLogical | kEidosValueMaskSingleton))->AddString_S(gEidosStr_filePath)->AddString("contents")->AddLogical_OS("append", gStaticEidosValue_LogicalF)->AddLogical_OS("compress", gStaticEidosValue_LogicalF));

eidos/eidos_functions.h

@@ -250,6 +250,7 @@ EidosValue_SP Eidos_ExecuteFunction_filesAtPath(const std::vector<EidosValue_SP>
 EidosValue_SP Eidos_ExecuteFunction_flushFile(const std::vector<EidosValue_SP> &p_arguments, EidosInterpreter &p_interpreter);
 EidosValue_SP Eidos_ExecuteFunction_getwd(const std::vector<EidosValue_SP> &p_arguments, EidosInterpreter &p_interpreter);
 EidosValue_SP Eidos_ExecuteFunction_readFile(const std::vector<EidosValue_SP> &p_arguments, EidosInterpreter &p_interpreter);
+EidosValue_SP Eidos_ExecuteFunction_readLine(const std::vector<EidosValue_SP> &p_arguments, EidosInterpreter &p_interpreter);
 EidosValue_SP Eidos_ExecuteFunction_setwd(const std::vector<EidosValue_SP> &p_arguments, EidosInterpreter &p_interpreter);
 EidosValue_SP Eidos_ExecuteFunction_tempdir(const std::vector<EidosValue_SP> &p_arguments, EidosInterpreter &p_interpreter);
 EidosValue_SP Eidos_ExecuteFunction_writeFile(const std::vector<EidosValue_SP> &p_arguments, EidosInterpreter &p_interpreter);

eidos/eidos_functions_files.cpp

@@ -27,7 +27,9 @@
 #include <stdio.h>
 #include <dirent.h>
 #include <fstream>
+#include <iostream>
 #include <sys/stat.h>
+#include <unistd.h>
 
 #include "../eidos_zlib/zlib.h"
 
@@ -220,6 +222,37 @@ EidosValue_SP Eidos_ExecuteFunction_readFile(const std::vector<EidosValue_SP> &p
 	return result_SP;
 }
 
+//	(string$)readLine(void)
+EidosValue_SP Eidos_ExecuteFunction_readLine(const std::vector<EidosValue_SP> &p_arguments, EidosInterpreter &p_interpreter)
+{
+#pragma unused (p_arguments, p_interpreter)
+	
+	// This function was implemented by Chris Talbot 11/19/25 for use in reinforcement learning environments.
+	// Associated with issue #576.
+	
+#ifdef EIDOS_GUI
+	// Cannot be called inside a GUI since there is (presumably) no stdin set up for the script to read.
+	EIDOS_TERMINATION << "ERROR (Eidos_ExecuteFunction_readLine): function readLine() is not available in GUI environments (SLiMgui, SLiMguiLegacy, or EidosScribe)." << EidosTerminate(nullptr);
+#else
+	// Read a single line from stdin (command-line mode with interactive terminal)
+	std::string line;
+	
+	if (std::getline(std::cin, line))
+	{
+		// Control for CRLF vs. LF line endings
+		if (!line.empty() && line[line.size() - 1] == '\r')
+			line.pop_back();
+		
+		return EidosValue_SP(new (gEidosValuePool->AllocateChunk()) EidosValue_String(line));
+	}
+	else
+	{
+		// EOF or error reading from stdin; return empty string
+		return gStaticEidosValue_StringEmpty;
+	}
+#endif
+}
+
 //	(string$)setwd(string$ path)
 EidosValue_SP Eidos_ExecuteFunction_setwd(const std::vector<EidosValue_SP> &p_arguments, __attribute__((unused)) EidosInterpreter &p_interpreter)
 {

treerec/_README

@@ -8,6 +8,8 @@ The code in `tskit/kastore` is from [kastore/c](https://github.com/tskit-dev/kas
 
 **Changelog:**
 
+- *30 November 2025*: updated to tskit 1.0.0
+- *24 September 2025*: updated to tskit C_1.2.0 (=1.0.0b2) (still kastore C_2.1.1)
 - *15 April 2025*: updated to tskit C_1.1.4 (=0.6.2) (still kastore C_2.1.1)
 - *7 November 2023*: updated to tskit C_1.1.2 (=0.5.5) (still kastore C_2.1.1)
 - *2 June 2022*: updated to tskit C_1.0.0 and kastore C_2.1.1

treerec/implementation.md

@@ -222,3 +222,47 @@ To read in a tree sequence, SLiM requires that:
 This is because of how `__TabulateSubpopulationsFromTreeSequence` works;
 probably it could be made more general, but it isn't.
 
+## Multiple chromosomes and multiple tree sequences
+
+In SLiM 5.0 we added support for simulating multiple chromosomes in SLiM, and that added
+some wrinkles to how tree-sequence recording is implemented in SLiM.  In particular, we
+now keep one tree sequence per chromosome being simulated.  In `species.h` we now define
+a struct named `TreeSeqInfo` that keeps the `tsk_table_collection_t`, in particular, as
+well as a couple of other bits of information.  Each Species then keeps its own `std::vector`
+of `TreeSeqInfo` structs, named `treeseq_`.  Where operations used to involve the table
+collection, they now typically involve a loop over the elements of `treeseq_` to perform
+the operation on each table collection in turn.
+
+The main complication here is that all of these table collections share three tskit tables:
+the node, individual, and population tables.  This means that the table collections all
+have a shared structure, and that needs to be preserved across operations like simplify.
+The shared tables are kept by the first chromosome's table collection, in `treeseq_[0]`.
+The other table collections zero-fill their node, individual, and population tables most
+of the time.  That means that, in that state, they are not compliant with tskit
+requirements, and trying to use them will often produce a segfault due to a dereference
+of a NULL pointer.  That is deliberate and useful; it makes it easy to debug situations
+where the table collections are being used when they should not be.  Sometimes we want the
+table collections to actually be usable.  For that, `Species::CopySharedTablesIn()` will
+do a bitwise, shallow copy of the shared tables into a given table collection; it should
+be matched by `DisconnectCopiedSharedTables()` as soon as the operation is done, restoring
+the zero-filled table state.  See https://github.com/tskit-dev/tskit/pull/2665 for
+Jerome's original multi-chromosome parallel simplification example, from which this design
+was derived.
+
+The end goal is that this will allow parallel simplification to happen in SLiM.  The code
+in `Species::SimplifyAllTreeSequences()` now implements the extra bookkeeping needed to
+maintain the shared table structure, so the design is ready to be parallelized when I
+return to the parallelization project.
+
+Single-chromosome simulations still get saved out as a .trees tree sequence file.  With
+multiple chromosomes, we now save out a "trees archive" with one .trees file per chromosome.
+This is not directly supported on the tskit side at the moment; you can just iterate over
+the files in the trees archive and process them as you wish in Python.  An example of this
+is provided in the SimHumanity paper, https://doi.org/10.47248/hpgg2505040006.
+
+There is some new top-level metadata associated with trees archives and multiple chromosomes.
+This is detailed in section 29.1 of the SLiM manual.  In particular, there are new top-level
+keys `this_chromosome` and `chromosomes` that should be provided in every .trees file.  The
+`chromosomes` key provides a table of all of the chromosomes involved in the trees archive.
+The `this_chromosome` key provides information about the particular chromosome represented
+by one particular .trees file in the trees archive.