Particles: Emit particles over time

This adds a basic internal emitter for every particle simulation.
The emitter cannot be controlled by the user yet. That will
come next.

10 files changed, 491 insertions, 42 deletions

diff --git a/source/blender/simulation/CMakeLists.txt b/source/blender/simulation/CMakeLists.txt
index a19e96e1a91..243b056db74 100644
--- a/source/blender/simulation/CMakeLists.txt
+++ b/source/blender/simulation/CMakeLists.txt
@@ -41,6 +41,7 @@ set(SRC
   intern/hair_volume.cpp
   intern/implicit_blender.c
   intern/implicit_eigen.cpp
+  intern/particle_allocator.cc
   intern/particle_function.cc
   intern/simulation_collect_influences.cc
   intern/simulation_solver.cc
@@ -49,11 +50,13 @@ set(SRC
   intern/ConstrainedConjugateGradient.h
   intern/eigen_utils.h
   intern/implicit.h
+  intern/particle_allocator.hh
+  intern/particle_function.hh
   intern/simulation_collect_influences.hh
   intern/simulation_solver.hh
+  intern/time_interval.hh
 
   SIM_mass_spring.h
-  SIM_particle_function.hh
   SIM_simulation_update.hh
 )
 
diff --git a/source/blender/simulation/intern/particle_allocator.cc b/source/blender/simulation/intern/particle_allocator.cc
new file mode 100644
index 00000000000..b65c0197c76
--- /dev/null
+++ b/source/blender/simulation/intern/particle_allocator.cc
@@ -0,0 +1,77 @@
+/*
+ * This program is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU General Public License
+ * as published by the Free Software Foundation; either version 2
+ * of the License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
+ */
+
+#include "particle_allocator.hh"
+
+namespace blender::sim {
+
+AttributesAllocator::~AttributesAllocator()
+{
+  for (std::unique_ptr<AttributesBlock> &block : allocated_blocks_) {
+    for (uint i : attributes_info_.index_range()) {
+      const fn::CPPType &type = attributes_info_.type_of(i);
+      type.destruct_n(block->buffers[i], block->size);
+      MEM_freeN(block->buffers[i]);
+    }
+  }
+}
+
+fn::MutableAttributesRef AttributesAllocator::allocate_uninitialized(uint size)
+{
+  std::unique_ptr<AttributesBlock> block = std::make_unique<AttributesBlock>();
+  block->buffers = Array<void *>(attributes_info_.size(), nullptr);
+  block->size = size;
+
+  for (uint i : attributes_info_.index_range()) {
+    const fn::CPPType &type = attributes_info_.type_of(i);
+    void *buffer = MEM_mallocN_aligned(size * type.size(), type.alignment(), AT);
+    block->buffers[i] = buffer;
+  }
+
+  fn::MutableAttributesRef attributes{attributes_info_, block->buffers, size};
+
+  {
+    std::lock_guard lock{mutex_};
+    allocated_blocks_.append(std::move(block));
+    allocated_attributes_.append(attributes);
+    total_allocated_ += size;
+  }
+
+  return attributes;
+}
+
+fn::MutableAttributesRef ParticleAllocator::allocate(uint size)
+{
+  const fn::AttributesInfo &info = attributes_allocator_.attributes_info();
+  fn::MutableAttributesRef attributes = attributes_allocator_.allocate_uninitialized(size);
+  for (uint i : info.index_range()) {
+    const fn::CPPType &type = info.type_of(i);
+    StringRef name = info.name_of(i);
+    if (name == "ID") {
+      uint start_id = next_id_.fetch_add(size);
+      MutableSpan<int> ids = attributes.get<int>("ID");
+      for (uint pindex : IndexRange(size)) {
+        ids[pindex] = start_id + pindex;
+      }
+    }
+    else {
+      type.fill_uninitialized(info.default_of(i), attributes.get(i).buffer(), size);
+    }
+  }
+  return attributes;
+}
+
+}  // namespace blender::sim
diff --git a/source/blender/simulation/intern/particle_allocator.hh b/source/blender/simulation/intern/particle_allocator.hh
new file mode 100644
index 00000000000..f854413c9aa
--- /dev/null
+++ b/source/blender/simulation/intern/particle_allocator.hh
@@ -0,0 +1,95 @@
+/*
+ * This program is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU General Public License
+ * as published by the Free Software Foundation; either version 2
+ * of the License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
+ */
+
+#ifndef __SIM_PARTICLE_ALLOCATOR_HH__
+#define __SIM_PARTICLE_ALLOCATOR_HH__
+
+#include "BLI_array.hh"
+#include "BLI_vector.hh"
+
+#include "FN_attributes_ref.hh"
+
+#include <atomic>
+#include <mutex>
+
+namespace blender::sim {
+
+class AttributesAllocator : NonCopyable, NonMovable {
+ private:
+  struct AttributesBlock {
+    Array<void *> buffers;
+    uint size;
+  };
+
+  const fn::AttributesInfo &attributes_info_;
+  Vector<std::unique_ptr<AttributesBlock>> allocated_blocks_;
+  Vector<fn::MutableAttributesRef> allocated_attributes_;
+  uint total_allocated_ = 0;
+  std::mutex mutex_;
+
+ public:
+  AttributesAllocator(const fn::AttributesInfo &attributes_info)
+      : attributes_info_(attributes_info)
+  {
+  }
+
+  ~AttributesAllocator();
+
+  Span<fn::MutableAttributesRef> get_allocations() const
+  {
+    return allocated_attributes_;
+  }
+
+  uint total_allocated() const
+  {
+    return total_allocated_;
+  }
+
+  const fn::AttributesInfo &attributes_info() const
+  {
+    return attributes_info_;
+  }
+
+  fn::MutableAttributesRef allocate_uninitialized(uint size);
+};
+
+class ParticleAllocator : NonCopyable, NonMovable {
+ private:
+  AttributesAllocator attributes_allocator_;
+  std::atomic<uint> next_id_;
+
+ public:
+  ParticleAllocator(const fn::AttributesInfo &attributes_info, uint next_id)
+      : attributes_allocator_(attributes_info), next_id_(next_id)
+  {
+  }
+
+  Span<fn::MutableAttributesRef> get_allocations() const
+  {
+    return attributes_allocator_.get_allocations();
+  }
+
+  uint total_allocated() const
+  {
+    return attributes_allocator_.total_allocated();
+  }
+
+  fn::MutableAttributesRef allocate(uint size);
+};
+
+}  // namespace blender::sim
+
+#endif /* __SIM_PARTICLE_ALLOCATOR_HH__ */
diff --git a/source/blender/simulation/intern/particle_function.cc b/source/blender/simulation/intern/particle_function.cc
index 7a0c9ccdb13..3788fd17e36 100644
--- a/source/blender/simulation/intern/particle_function.cc
+++ b/source/blender/simulation/intern/particle_function.cc
@@ -14,7 +14,7 @@
  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
  */
 
-#include "SIM_particle_function.hh"
+#include "particle_function.hh"
 
 namespace blender::sim {
 
@@ -49,12 +49,11 @@ ParticleFunction::ParticleFunction(const fn::MultiFunction *global_fn,
   }
 }
 
-ParticleFunctionEvaluator::ParticleFunctionEvaluator(const ParticleFunction &particle_fn,
-                                                     IndexMask mask,
-                                                     fn::AttributesRef particle_attributes)
+ParticleFunctionEvaluator::ParticleFunctionEvaluator(
+    const ParticleFunction &particle_fn, const ParticleChunkContext &particle_chunk_context)
     : particle_fn_(particle_fn),
-      mask_(mask),
-      particle_attributes_(particle_attributes),
+      particle_chunk_context_(particle_chunk_context),
+      mask_(particle_chunk_context_.index_mask()),
       outputs_(particle_fn_.output_types_.size(), nullptr)
 {
 }
@@ -112,7 +111,7 @@ void ParticleFunctionEvaluator::compute_globals()
 
   /* Add input parameters. */
   for (const ParticleFunctionInput *input : particle_fn_.global_inputs_) {
-    input->add_input(particle_attributes_, params, resources_);
+    input->add_input(particle_chunk_context_.attributes(), params, resources_);
   }
 
   /* Add output parameters. */
@@ -139,7 +138,7 @@ void ParticleFunctionEvaluator::compute_per_particle()
 
   /* Add input parameters. */
   for (const ParticleFunctionInput *input : particle_fn_.per_particle_inputs_) {
-    input->add_input(particle_attributes_, params, resources_);
+    input->add_input(particle_chunk_context_.attributes(), params, resources_);
   }
 
   /* Add output parameters. */
diff --git a/source/blender/simulation/SIM_particle_function.hh b/source/blender/simulation/intern/particle_function.hh
index eae61b1ae11..abed9063bae 100644
--- a/source/blender/simulation/SIM_particle_function.hh
+++ b/source/blender/simulation/intern/particle_function.hh
@@ -22,6 +22,8 @@
 
 #include "BLI_resource_collector.hh"
 
+#include "simulation_solver.hh"
+
 namespace blender::sim {
 
 class ParticleFunctionInput {
@@ -58,17 +60,16 @@ class ParticleFunctionEvaluator {
  private:
   ResourceCollector resources_;
   const ParticleFunction &particle_fn_;
+  const ParticleChunkContext &particle_chunk_context_;
   IndexMask mask_;
   fn::MFContextBuilder global_context_;
   fn::MFContextBuilder per_particle_context_;
-  fn::AttributesRef particle_attributes_;
   Vector<void *> outputs_;
   bool is_computed_ = false;
 
  public:
   ParticleFunctionEvaluator(const ParticleFunction &particle_fn,
-                            IndexMask mask,
-                            fn::AttributesRef particle_attributes);
+                            const ParticleChunkContext &particle_chunk_context);
   ~ParticleFunctionEvaluator();
 
   void compute();
diff --git a/source/blender/simulation/intern/simulation_collect_influences.cc b/source/blender/simulation/intern/simulation_collect_influences.cc
index 3e35b040dbe..ec1e7db0f72 100644
--- a/source/blender/simulation/intern/simulation_collect_influences.cc
+++ b/source/blender/simulation/intern/simulation_collect_influences.cc
@@ -15,7 +15,7 @@
  */
 
 #include "simulation_collect_influences.hh"
-#include "SIM_particle_function.hh"
+#include "particle_function.hh"
 
 #include "FN_attributes_ref.hh"
 #include "FN_multi_function_network_evaluation.hh"
@@ -23,6 +23,8 @@
 
 #include "NOD_node_tree_multi_function.hh"
 
+#include "BLI_rand.hh"
+
 namespace blender::sim {
 
 struct DummyDataSources {
@@ -191,12 +193,15 @@ class ParticleFunctionForce : public ParticleForce {
   {
   }
 
-  void add_force(fn::AttributesRef attributes, MutableSpan<float3> r_combined_force) const override
+  void add_force(ParticleForceContext &context) const override
   {
-    IndexMask mask = IndexRange(attributes.size());
-    ParticleFunctionEvaluator evaluator{particle_fn_, mask, attributes};
+    IndexMask mask = context.particle_chunk().index_mask();
+    MutableSpan<float3> r_combined_force = context.force_dst();
+
+    ParticleFunctionEvaluator evaluator{particle_fn_, context.particle_chunk()};
     evaluator.compute();
     fn::VSpan<float3> forces = evaluator.get<float3>(0, "Force");
+
     for (uint i : mask) {
       r_combined_force[i] += forces[i];
     }
@@ -247,6 +252,48 @@ static void collect_forces(nodes::MFNetworkTreeMap &network_map,
   }
 }
 
+class MyBasicEmitter : public ParticleEmitter {
+ private:
+  std::string name_;
+
+ public:
+  MyBasicEmitter(std::string name) : name_(std::move(name))
+  {
+  }
+
+  void emit(ParticleEmitterContext &context) const override
+  {
+    ParticleAllocator *allocator = context.try_get_particle_allocator(name_);
+    if (allocator == nullptr) {
+      return;
+    }
+
+    fn::MutableAttributesRef attributes = allocator->allocate(10);
+    RandomNumberGenerator rng{(uint)context.simulation_time_interval().start() ^
+                              DefaultHash<std::string>{}(name_)};
+
+    MutableSpan<float3> positions = attributes.get<float3>("Position");
+    MutableSpan<float3> velocities = attributes.get<float3>("Velocity");
+
+    for (uint i : IndexRange(attributes.size())) {
+      positions[i] = rng.get_unit_float3();
+      velocities[i] = rng.get_unit_float3();
+    }
+  }
+};
+
+static void collect_emitters(nodes::MFNetworkTreeMap &network_map,
+                             ResourceCollector &resources,
+                             SimulationInfluences &r_influences)
+{
+  for (const nodes::DNode *dnode :
+       network_map.tree().nodes_by_type("SimulationNodeParticleSimulation")) {
+    std::string name = dnode_to_path(*dnode);
+    ParticleEmitter &emitter = resources.construct<MyBasicEmitter>(AT, name);
+    r_influences.particle_emitters.append(&emitter);
+  }
+}
+
 void collect_simulation_influences(Simulation &simulation,
                                    ResourceCollector &resources,
                                    SimulationInfluences &r_influences,
@@ -267,6 +314,7 @@ void collect_simulation_influences(Simulation &simulation,
   // WM_clipboard_text_set(network.to_dot().c_str(), false);
 
   collect_forces(network_map, resources, data_sources, r_influences);
+  collect_emitters(network_map, resources, r_influences);
 
   for (const nodes::DNode *dnode : tree.nodes_by_type("SimulationNodeParticleSimulation")) {
     r_states_info.particle_simulation_names.add(dnode_to_path(*dnode));
diff --git a/source/blender/simulation/intern/simulation_solver.cc b/source/blender/simulation/intern/simulation_solver.cc
index c51e8e1a37f..f158bd6e8bf 100644
--- a/source/blender/simulation/intern/simulation_solver.cc
+++ b/source/blender/simulation/intern/simulation_solver.cc
@@ -26,6 +26,10 @@ ParticleForce::~ParticleForce()
 {
 }
 
+ParticleEmitter::~ParticleEmitter()
+{
+}
+
 class CustomDataAttributesRef {
  private:
   Vector<void *> buffers_;
@@ -84,48 +88,61 @@ void initialize_simulation_states(Simulation &simulation,
                                   Depsgraph &UNUSED(depsgraph),
                                   const SimulationInfluences &UNUSED(influences))
 {
-  RandomNumberGenerator rng;
-
-  LISTBASE_FOREACH (ParticleSimulationState *, state, &simulation.states) {
-    state->tot_particles = 1000;
-    CustomData_realloc(&state->attributes, state->tot_particles);
-    ensure_attributes_exist(state);
-
-    CustomDataAttributesRef custom_data_attributes{state->attributes, (uint)state->tot_particles};
-
-    fn::MutableAttributesRef attributes = custom_data_attributes;
-    MutableSpan<float3> positions = attributes.get<float3>("Position");
-    MutableSpan<float3> velocities = attributes.get<float3>("Velocity");
-    MutableSpan<int32_t> ids = attributes.get<int32_t>("ID");
-
-    for (uint i : positions.index_range()) {
-      positions[i] = {i / 100.0f, 0, 0};
-      velocities[i] = {0, rng.get_float() - 0.5f, rng.get_float() - 0.5f};
-      ids[i] = i;
-    }
-  }
+  simulation.current_simulation_time = 0.0f;
 }
 
 void solve_simulation_time_step(Simulation &simulation,
-                                Depsgraph &UNUSED(depsgraph),
+                                Depsgraph &depsgraph,
                                 const SimulationInfluences &influences,
                                 float time_step)
 {
+  SimulationSolveContext solve_context{simulation, depsgraph, influences};
+  TimeInterval simulation_time_interval{simulation.current_simulation_time, time_step};
+
+  Map<std::string, std::unique_ptr<fn::AttributesInfo>> attribute_infos;
+  Map<std::string, std::unique_ptr<ParticleAllocator>> particle_allocators;
   LISTBASE_FOREACH (ParticleSimulationState *, state, &simulation.states) {
     ensure_attributes_exist(state);
 
-    CustomDataAttributesRef custom_data_attributes{state->attributes, (uint)state->tot_particles};
+    fn::AttributesInfoBuilder builder;
+    CustomData &custom_data = state->attributes;
+    for (const CustomDataLayer &layer : Span(custom_data.layers, custom_data.totlayer)) {
+      switch (layer.type) {
+        case CD_PROP_INT32: {
+          builder.add<int32_t>(layer.name, 0);
+          break;
+        }
+        case CD_PROP_FLOAT3: {
+          builder.add<float3>(layer.name, {0, 0, 0});
+          break;
+        }
+      }
+    }
+    auto info = std::make_unique<fn::AttributesInfo>(builder);
+    particle_allocators.add_new(
+        state->head.name, std::make_unique<ParticleAllocator>(*info, state->next_particle_id));
+    attribute_infos.add_new(state->head.name, std::move(info));
+  }
 
+  LISTBASE_FOREACH (ParticleSimulationState *, state, &simulation.states) {
+    CustomDataAttributesRef custom_data_attributes{state->attributes, (uint)state->tot_particles};
     fn::MutableAttributesRef attributes = custom_data_attributes;
+
     MutableSpan<float3> positions = attributes.get<float3>("Position");
     MutableSpan<float3> velocities = attributes.get<float3>("Velocity");
 
     Array<float3> force_vectors{(uint)state->tot_particles, {0, 0, 0}};
     const Vector<const ParticleForce *> *forces = influences.particle_forces.lookup_ptr(
         state->head.name);
+
     if (forces != nullptr) {
+      ParticleChunkContext particle_chunk_context{IndexMask((uint)state->tot_particles),
+                                                  attributes};
+      ParticleForceContext particle_force_context{
+          solve_context, particle_chunk_context, force_vectors};
+
       for (const ParticleForce *force : *forces) {
-        force->add_force(attributes, force_vectors);
+        force->add_force(particle_force_context);
       }
     }
 
@@ -134,6 +151,44 @@ void solve_simulation_time_step(Simulation &simulation,
       positions[i] += velocities[i] * time_step;
     }
   }
+
+  for (const ParticleEmitter *emitter : influences.particle_emitters) {
+    ParticleEmitterContext emitter_context{
+        solve_context, particle_allocators, simulation_time_interval};
+    emitter->emit(emitter_context);
+  }
+
+  LISTBASE_FOREACH (ParticleSimulationState *, state, &simulation.states) {
+    ParticleAllocator &allocator = *particle_allocators.lookup_as(state->head.name);
+
+    const uint emitted_particle_amount = allocator.total_allocated();
+    const uint old_particle_amount = state->tot_particles;
+    const uint new_particle_amount = old_particle_amount + emitted_particle_amount;
+
+    CustomData_realloc(&state->attributes, new_particle_amount);
+
+    CustomDataAttributesRef custom_data_attributes{state->attributes, new_particle_amount};
+    fn::MutableAttributesRef attributes = custom_data_attributes;
+
+    uint offset = old_particle_amount;
+    for (fn::MutableAttributesRef emitted_attributes : allocator.get_allocations()) {
+      fn::MutableAttributesRef dst_attributes = attributes.slice(
+          IndexRange(offset, emitted_attributes.size()));
+      for (uint attribute_index : attributes.info().index_range()) {
+        fn::GMutableSpan emitted_data = emitted_attributes.get(attribute_index);
+        fn::GMutableSpan dst = dst_attributes.get(attribute_index);
+        const fn::CPPType &type = dst.type();
+        type.copy_to_uninitialized_n(
+            emitted_data.buffer(), dst.buffer(), emitted_attributes.size());
+      }
+      offset += emitted_attributes.size();
+    }
+
+    state->tot_particles = new_particle_amount;
+    state->next_particle_id += emitted_particle_amount;
+  }
+
+  simulation.current_simulation_time = simulation_time_interval.end();
 }
 
 }  // namespace blender::sim
diff --git a/source/blender/simulation/intern/simulation_solver.hh b/source/blender/simulation/intern/simulation_solver.hh
index 818bcf067a3..5a45a676ec7 100644
--- a/source/blender/simulation/intern/simulation_solver.hh
+++ b/source/blender/simulation/intern/simulation_solver.hh
@@ -24,19 +24,133 @@
 
 #include "FN_attributes_ref.hh"
 
+#include "particle_allocator.hh"
+#include "time_interval.hh"
+
 struct Depsgraph;
 
 namespace blender::sim {
 
+class ParticleEmitterContext;
+class ParticleForceContext;
+
+class ParticleEmitter {
+ public:
+  virtual ~ParticleEmitter();
+  virtual void emit(ParticleEmitterContext &context) const = 0;
+};
+
 class ParticleForce {
  public:
   virtual ~ParticleForce();
-  virtual void add_force(fn::AttributesRef attributes,
-                         MutableSpan<float3> r_combined_force) const = 0;
+  virtual void add_force(ParticleForceContext &context) const = 0;
 };
 
 struct SimulationInfluences {
   Map<std::string, Vector<const ParticleForce *>> particle_forces;
+  Vector<const ParticleEmitter *> particle_emitters;
+};
+
+class SimulationSolveContext {
+ private:
+  Simulation &simulation_;
+  Depsgraph &depsgraph_;
+  const SimulationInfluences &influences_;
+
+ public:
+  SimulationSolveContext(Simulation &simulation,
+                         Depsgraph &depsgraph,
+                         const SimulationInfluences &influences)
+      : simulation_(simulation), depsgraph_(depsgraph), influences_(influences)
+  {
+  }
+};
+
+class ParticleChunkContext {
+ private:
+  IndexMask index_mask_;
+  fn::MutableAttributesRef attributes_;
+
+ public:
+  ParticleChunkContext(IndexMask index_mask, fn::MutableAttributesRef attributes)
+      : index_mask_(index_mask), attributes_(attributes)
+  {
+  }
+
+  IndexMask index_mask() const
+  {
+    return index_mask_;
+  }
+
+  fn::MutableAttributesRef attributes()
+  {
+    return attributes_;
+  }
+
+  fn::AttributesRef attributes() const
+  {
+    return attributes_;
+  }
+};
+
+class ParticleEmitterContext {
+ private:
+  SimulationSolveContext &solve_context_;
+  Map<std::string, std::unique_ptr<ParticleAllocator>> &particle_allocators_;
+  TimeInterval simulation_time_interval_;
+
+ public:
+  ParticleEmitterContext(SimulationSolveContext &solve_context,
+                         Map<std::string, std::unique_ptr<ParticleAllocator>> &particle_allocators,
+                         TimeInterval simulation_time_interval)
+      : solve_context_(solve_context),
+        particle_allocators_(particle_allocators),
+        simulation_time_interval_(simulation_time_interval)
+  {
+  }
+
+  ParticleAllocator *try_get_particle_allocator(StringRef particle_simulation_name)
+  {
+    auto *ptr = particle_allocators_.lookup_ptr_as(particle_simulation_name);
+    if (ptr != nullptr) {
+      return ptr->get();
+    }
+    else {
+      return nullptr;
+    }
+  }
+
+  TimeInterval simulation_time_interval() const
+  {
+    return simulation_time_interval_;
+  }
+};
+
+class ParticleForceContext {
+ private:
+  SimulationSolveContext &solve_context_;
+  const ParticleChunkContext &particle_chunk_context_;
+  MutableSpan<float3> force_dst_;
+
+ public:
+  ParticleForceContext(SimulationSolveContext &solve_context,
+                       const ParticleChunkContext &particle_chunk_context,
+                       MutableSpan<float3> force_dst)
+      : solve_context_(solve_context),
+        particle_chunk_context_(particle_chunk_context),
+        force_dst_(force_dst)
+  {
+  }
+
+  const ParticleChunkContext &particle_chunk() const
+  {
+    return particle_chunk_context_;
+  }
+
+  MutableSpan<float3> force_dst()
+  {
+    return force_dst_;
+  }
 };
 
 void initialize_simulation_states(Simulation &simulation,
diff --git a/source/blender/simulation/intern/simulation_update.cc b/source/blender/simulation/intern/simulation_update.cc
index d8b7e1df19a..5328304863b 100644
--- a/source/blender/simulation/intern/simulation_update.cc
+++ b/source/blender/simulation/intern/simulation_update.cc
@@ -14,7 +14,6 @@
  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
  */
 
-#include "SIM_particle_function.hh"
 #include "SIM_simulation_update.hh"
 
 #include "BKE_customdata.h"
@@ -32,6 +31,7 @@
 #include "BLI_rand.h"
 #include "BLI_vector.hh"
 
+#include "particle_function.hh"
 #include "simulation_collect_influences.hh"
 #include "simulation_solver.hh"
 
diff --git a/source/blender/simulation/intern/time_interval.hh b/source/blender/simulation/intern/time_interval.hh
new file mode 100644
index 00000000000..49600dd10de
--- /dev/null
+++ b/source/blender/simulation/intern/time_interval.hh
@@ -0,0 +1,57 @@
+/*
+ * This program is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU General Public License
+ * as published by the Free Software Foundation; either version 2
+ * of the License, or (at your option) any later version.
+ *
+ * This program is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+ * GNU General Public License for more details.
+ *
+ * You should have received a copy of the GNU General Public License
+ * along with this program; if not, write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
+ */
+
+#ifndef __SIM_TIME_INTERVAL_HH__
+#define __SIM_TIME_INTERVAL_HH__
+
+#include "BLI_utildefines.h"
+
+namespace blender::sim {
+
+/**
+ * The start time is inclusive and the end time is exclusive. The duration is zero, the interval
+ * describes a single point in time.
+ */
+class TimeInterval {
+ private:
+  float start_;
+  float duration_;
+
+ public:
+  TimeInterval(float start, float duration) : start_(start), duration_(duration)
+  {
+    BLI_assert(duration_ >= 0.0f);
+  }
+
+  float start() const
+  {
+    return start_;
+  }
+
+  float end() const
+  {
+    return start_ + duration_;
+  }
+
+  float duration() const
+  {
+    return duration_;
+  }
+};
+
+}  // namespace blender::sim
+
+#endif /* __SIM_TIME_INTERVAL_HH__ */