Filippo-style simulations

examples/auto_charge.cpp

@@ -0,0 +1,21 @@
+#include <dsf.hpp>
+#include <iostream>
+
+int main() {
+  using namespace dsf::mobility;
+  try {
+    TrafficSimulator sim;
+    sim.setName("auto_charge_example");
+    // Import a small example network (user must provide path)
+    sim.importRoadNetwork("examples/data/edges.csv", "");
+    sim.setTimeFrame(0, 600);                  // 600 steps
+    sim.saveData(60, true, true, true, true);  // save every 60 steps
+    // Run auto-charge: base 5 agents every 30 steps, save interval 60, maxSteps=600
+    sim.runAutoCharge(
+        5, 30, 60, 600, std::nullopt, 120, 120, 2, true, "auto_charge_events.csv");
+  } catch (const std::exception& e) {
+    std::cerr << "Error: " << e.what() << std::endl;
+    return 1;
+  }
+  return 0;
+}

examples/auto_charge_py_example.py

@@ -0,0 +1,22 @@
+from dsf_cpp import mobility
+
+sim = mobility.TrafficSimulator()
+sim.setName("auto_charge_py")
+# User should import a valid network first
+# sim.importRoadNetwork('examples/data/edges.csv')
+sim.setTimeFrame(0, 600)
+# save every 60 steps: avg_stats etc.
+sim.saveData(60, True, True, True, True)
+# Run AutoCharge via bindings
+sim.runAutoCharge(
+    5,  # baseAgentCount
+    30,  # dtAgent
+    60,  # saveIntervalSeconds
+    None,  # maxSteps
+    None,  # stopMeanDensityVpk
+    120,  # stabilityHoldSeconds
+    120,  # stabilityCooldownSeconds
+    2,  # chargeIncrement
+    True,  # injectOnCharge
+    "auto_charge_events.csv",  # stabilityLogFile
+)

src/dsf/bindings.cpp

@@ -1241,6 +1241,62 @@ PYBIND11_MODULE(dsf_cpp, m) {
           checkDeltaT:           Interval in seconds between occupancy checks.
           agentIncrement:        How many agents to add to the target at each check (default 1).
       )doc")
+      .def(
+          "runAutoCharge",
+          [](dsf::mobility::TrafficSimulator& self,
+             std::size_t baseAgentCount,
+             std::size_t dtAgent,
+             std::size_t saveIntervalSeconds,
+             pybind11::object maxSteps,
+             pybind11::object stopMeanDensityVpk,
+             std::size_t stabilityHoldSeconds,
+             std::size_t stabilityCooldownSeconds,
+             std::size_t chargeIncrement,
+             bool injectOnCharge,
+             std::string stabilityLogFile) {
+            std::optional<std::size_t> optMaxSteps = std::nullopt;
+            if (!maxSteps.is_none()) {
+              optMaxSteps = static_cast<std::size_t>(pybind11::cast<std::uint64_t>(maxSteps));
+            }
+            std::optional<double> optStopMeanDensity = std::nullopt;
+            if (!stopMeanDensityVpk.is_none()) {
+              optStopMeanDensity = pybind11::cast<double>(stopMeanDensityVpk);
+            }
+            self.runAutoCharge(baseAgentCount,
+                               dtAgent,
+                               saveIntervalSeconds,
+                               optMaxSteps,
+                               optStopMeanDensity,
+                               stabilityHoldSeconds,
+                               stabilityCooldownSeconds,
+                               chargeIncrement,
+                               injectOnCharge,
+                               stabilityLogFile);
+          },
+          pybind11::arg("baseAgentCount"),
+          pybind11::arg("dtAgent"),
+          pybind11::arg("saveIntervalSeconds"),
+          pybind11::arg("maxSteps") = pybind11::none(),
+          pybind11::arg("stopMeanDensityVpk") = pybind11::none(),
+          pybind11::arg("stabilityHoldSeconds") = 60,
+          pybind11::arg("stabilityCooldownSeconds") = 0,
+          pybind11::arg("chargeIncrement") = 1,
+          pybind11::arg("injectOnCharge") = true,
+          pybind11::arg("stabilityLogFile") = std::string(),
+          R"doc(Run the simulation in auto-charge mode.
+
+      Args:
+          baseAgentCount: Number of agents to insert every dtAgent steps.
+          dtAgent: Interval in steps between agent insertions.
+          saveIntervalSeconds: Interval in steps between save/check operations.
+          maxSteps: Optional maximum number of steps to run.
+          stopMeanDensityVpk: Optional threshold to stop when mean density reaches this value.
+          stabilityHoldSeconds: Seconds of stability required before auto-charging.
+          stabilityCooldownSeconds: Minimum seconds since last charge before charging again.
+          chargeIncrement: Number of agents to add when charging triggers.
+          injectOnCharge: Whether to inject agents immediately when charging.
+          stabilityLogFile: Optional CSV file path to append stability events.
+      )doc")
       .def(
           "database",
           [](dsf::mobility::TrafficSimulator& self) { return self.database(); },

src/dsf/mobility/TrafficSimulator.cpp

@@ -434,6 +434,139 @@ namespace dsf::mobility {
     }
   }
 
+  void TrafficSimulator::m_runAutoCharge(std::size_t baseAgentCount,
+                                         std::size_t dtAgent,
+                                         std::size_t saveIntervalSeconds,
+                                         std::optional<std::size_t> maxSteps,
+                                         std::optional<double> stopMeanDensityVpk,
+                                         std::size_t stabilityHoldSeconds,
+                                         std::size_t stabilityCooldownSeconds,
+                                         std::size_t chargeIncrement,
+                                         bool injectOnCharge,
+                                         std::string const& stabilityLogFile) {
+    if (m_endTime < m_initTime) {
+      throw std::runtime_error(
+          "End time must be greater than or equal to initial time for the simulation.");
+    }
+    if (dtAgent == 0) {
+      throw std::invalid_argument("dtAgent must be > 0");
+    }
+
+    auto const totalTimeSteps = static_cast<std::size_t>(m_endTime - m_initTime);
+
+    m_preparePersistence();
+
+    spdlog::info("Starting auto-charge run from {} to {} ({} time steps).",
+                 m_timeToStr(m_initTime),
+                 m_timeToStr(m_endTime),
+                 totalTimeSteps);
+    auto pbar = dsf::utility::default_progress_bar("Running simulation", totalTimeSteps);
+
+    std::size_t base = baseAgentCount;
+    std::size_t last_charge_time_step = 0;
+    std::optional<std::size_t> stable_since_time_step = std::nullopt;
+    // Prepare moving window for stability detection
+    std::deque<double> densityWindow;
+    auto windowSize = std::size_t(1);
+    if (saveIntervalSeconds > 0 && stabilityHoldSeconds > 0) {
+      windowSize = std::max<std::size_t>(1, stabilityHoldSeconds / saveIntervalSeconds);
+    }
+    const double stabilityRelStdThreshold = 0.01;  // relative std/mean threshold
+
+    std::optional<CSVWriter> stabilityWriter;
+    if (!stabilityLogFile.empty()) {
+      bool fileExists = std::filesystem::exists(stabilityLogFile);
+      stabilityWriter.emplace(stabilityLogFile, ';');
+      if (!fileExists) {
+        stabilityWriter->writeHeader("time_step", "base_agent_count", "mean_density");
+      }
+    }
+
+    for (std::size_t i = 0; i < totalTimeSteps; ++i) {
+      if ((m_updatePathDeltaT > 0 && i % m_updatePathDeltaT == 0) || (i == 0)) {
+        m_dynamics->updatePaths();
+      }
+
+      if (maxSteps.has_value() && i >= maxSteps.value()) {
+        spdlog::info("Reached max steps ({}). Stopping.", maxSteps.value());
+        break;
+      }
+
+      if (i % dtAgent == 0) {
+        if (base > 0) {
+          m_dynamics->addAgents(base, AgentInsertionMethod::ODS);
+        }
+      }
+
+      bool const shouldSave = (saveIntervalSeconds > 0 && i % saveIntervalSeconds == 0);
+      auto stepData = m_dynamics->evolve(shouldSave ? StepDataRequest{m_saveAverageStats,
+                                                                      m_saveStreetData,
+                                                                      m_saveTravelData,
+                                                                      m_saveAgentData}
+                                                    : StepDataRequest{});
+
+      if (shouldSave) {
+        if (stepData.averageStats.has_value()) {
+          double meanDensity = stepData.averageStats->meanDensity;
+
+          if (stopMeanDensityVpk.has_value() &&
+              meanDensity >= stopMeanDensityVpk.value()) {
+            spdlog::info(
+                "Stopping: mean_density_vpk={} reached threshold {} at time_step {}.",
+                meanDensity,
+                stopMeanDensityVpk.value(),
+                stepData.timeStep);
+            break;
+          }
+
+          // Push density into moving window
+          densityWindow.push_back(meanDensity);
+          if (densityWindow.size() > windowSize)
+            densityWindow.pop_front();
+          if (densityWindow.size() >= windowSize) {
+            // compute mean and std
+            double sum = 0.;
+            for (auto v : densityWindow)
+              sum += v;
+            double mean = sum / static_cast<double>(densityWindow.size());
+            double sq = 0.;
+            for (auto v : densityWindow)
+              sq += (v - mean) * (v - mean);
+            double std = std::sqrt(sq / static_cast<double>(densityWindow.size()));
+            double relStd = (mean > 0.) ? std / mean : std;
+            if (relStd <= stabilityRelStdThreshold) {
+              if (!stable_since_time_step.has_value()) {
+                stable_since_time_step = i;
+              } else if ((i - stable_since_time_step.value() >= stabilityHoldSeconds) &&
+                         (i - last_charge_time_step >= stabilityCooldownSeconds)) {
+                // Trigger charge
+                base += chargeIncrement;
+                last_charge_time_step = i;
+                stable_since_time_step = std::nullopt;
+                if (stabilityWriter.has_value()) {
+                  stabilityWriter->writeRow(i, base, meanDensity);
+                }
+                spdlog::info(
+                    "[auto-charge] New BASE_AGENT_COUNT: {} at step {}.", base, i);
+                if (injectOnCharge && chargeIncrement > 0) {
+                  m_dynamics->addAgents(chargeIncrement, AgentInsertionMethod::ODS);
+                }
+              }
+            } else {
+              stable_since_time_step = std::nullopt;
+            }
+          }
+        }
+
+        m_pendingStepData.push_back(std::move(stepData));
+        m_flushStepData(std::move(m_pendingStepData.back()));
+        m_pendingStepData.clear();
+      }
+
+      pbar->update();
+    }
+  }
+
   void TrafficSimulator::connectDataBase(std::string_view const dbPath,
                                          std::string_view const queries) {
     m_database = std::make_unique<SQLite::Database>(

src/dsf/mobility/TrafficSimulator.hpp

@@ -191,6 +191,16 @@ namespace dsf::mobility {
                          std::time_t const agentInsertionDeltaT,
                          std::time_t const checkDeltaT,
                          std::size_t const agentIncrement = 1);
+    void m_runAutoCharge(std::size_t baseAgentCount,
+                         std::size_t dtAgent,
+                         std::size_t saveIntervalSeconds,
+                         std::optional<std::size_t> maxSteps,
+                         std::optional<double> stopMeanDensityVpk,
+                         std::size_t stabilityHoldSeconds,
+                         std::size_t stabilityCooldownSeconds,
+                         std::size_t chargeIncrement,
+                         bool injectOnCharge,
+                         std::string const& stabilityLogFile);
 
   public:
     /// @brief Construct a new TrafficSimulator with a generated simulation id.
@@ -270,6 +280,34 @@ namespace dsf::mobility {
       m_runSlowCharge(nInitialAgents, agentInsertionDeltaT, checkDeltaT, agentIncrement);
     }
 
+    inline void runAutoCharge(std::size_t baseAgentCount,
+                              std::size_t dtAgent,
+                              std::size_t saveIntervalSeconds,
+                              std::optional<std::size_t> maxSteps = std::nullopt,
+                              std::optional<double> stopMeanDensityVpk = std::nullopt,
+                              std::size_t stabilityHoldSeconds = 60,
+                              std::size_t stabilityCooldownSeconds = 0,
+                              std::size_t chargeIncrement = 1,
+                              bool injectOnCharge = true,
+                              std::string const& stabilityLogFile = std::string()) {
+      if (m_dynamics == nullptr) {
+        throw std::runtime_error(
+            "Cannot run the simulation without imported road network dynamics.");
+      }
+      m_dynamics->prepareNetwork();
+      m_preparePersistence();
+      m_runAutoCharge(baseAgentCount,
+                      dtAgent,
+                      saveIntervalSeconds,
+                      maxSteps,
+                      stopMeanDensityVpk,
+                      stabilityHoldSeconds,
+                      stabilityCooldownSeconds,
+                      chargeIncrement,
+                      injectOnCharge,
+                      stabilityLogFile);
+    }
+
     /// @brief Get the database connection (const version)
     /// @return const SQLite::Database const*, The database connection
     inline auto* database() const { return m_database.get(); }