Line data    Source code

       1             : #include "Compartment.hpp"
       2             : #include "CollisionDetector.hpp"
       3             : #include "CollisionHandler.hpp"
       4             : #include "Disc.hpp"
       5             : #include "MathUtils.hpp"
       6             : #include "ReactionEngine.hpp"
       7             : 
       8             : namespace cell
       9             : {
      10             : 
      11          19 : Compartment::Compartment(Compartment* parent, Membrane membrane, SimulationContext simulationContext)
      12          19 :     : parent_(parent)
      13          19 :     , membrane_(std::move(membrane))
      14          19 :     , simulationContext_(std::move(simulationContext))
      15             : {
      16          19 :     membrane_.setCompartment(this);
      17          19 : }
      18             : 
      19          19 : Compartment::~Compartment() = default;
      20             : 
      21      362925 : const Membrane& Compartment::getMembrane() const
      22             : {
      23      362925 :     return membrane_;
      24             : }
      25             : 
      26           0 : void Compartment::setDiscs(std::vector<Disc>&& discs)
      27             : {
      28           0 :     discs_ = std::move(discs);
      29           0 : }
      30             : 
      31         319 : void Compartment::addDisc(Disc disc)
      32             : {
      33         319 :     discs_.push_back(std::move(disc));
      34         319 : }
      35             : 
      36          30 : const std::vector<Disc>& Compartment::getDiscs() const
      37             : {
      38          30 :     return discs_;
      39             : }
      40             : 
      41           2 : void Compartment::addIntrudingDisc(Disc& disc)
      42             : {
      43           2 :     intrudingDiscs_.push_back(&disc);
      44           2 : }
      45             : 
      46      462360 : std::vector<std::unique_ptr<Compartment>>& Compartment::getCompartments()
      47             : {
      48      462360 :     return compartments_;
      49             : }
      50             : 
      51          10 : const std::vector<std::unique_ptr<Compartment>>& Compartment::getCompartments() const
      52             : {
      53          10 :     return compartments_;
      54             : }
      55             : 
      56           3 : const Compartment* Compartment::getParent() const
      57             : {
      58           3 :     return parent_;
      59             : }
      60             : 
      61           9 : auto Compartment::detectCollisions()
      62             : {
      63           9 :     simulationContext_.collisionDetector.buildEntries(discs_, membranes_, intrudingDiscs_);
      64           9 :     return simulationContext_.collisionDetector.detectCollisions(
      65           9 :         CollisionDetector::Params{.discs = &discs_,
      66           9 :                                   .membranes = &membranes_,
      67           9 :                                   .intrudingDiscs = &intrudingDiscs_,
      68          18 :                                   .containingMembrane = &membrane_});
      69             : }
      70             : 
      71           9 : void Compartment::update(double dt)
      72             : {
      73           9 :     auto detectedCollisions = detectCollisions();
      74             : 
      75           9 :     simulationContext_.reactionEngine.applyBimolecularReactions(detectedCollisions);
      76           9 :     moveDiscsIntoChildCompartments(detectedCollisions);
      77           9 :     moveDiscsIntoParentCompartment(detectedCollisions);
      78           9 :     simulationContext_.collisionHandler.resolveCollisions(detectedCollisions);
      79           9 :     updateChildCompartments(dt);
      80           9 :     moveDiscsAndApplyUnimolecularReactions(dt); // Removes destroyed discs
      81           9 :     intrudingDiscs_.clear();
      82           9 : }
      83             : 
      84          10 : Compartment* Compartment::createSubCompartment(Membrane membrane)
      85             : {
      86          10 :     auto compartment = std::make_unique<Compartment>(this, std::move(membrane), simulationContext_);
      87          10 :     membranes_.push_back(compartment->getMembrane());
      88          10 :     compartments_.push_back(std::move(compartment));
      89             : 
      90          20 :     return compartments_.back().get();
      91          10 : }
      92             : 
      93           9 : void Compartment::moveDiscsAndApplyUnimolecularReactions(double dt)
      94             : {
      95           9 :     std::vector<Disc> newDiscs;
      96             : 
      97          27 :     for (std::size_t i = 0; i < discs_.size(); ++i)
      98             :     {
      99          18 :         auto& disc = discs_[i];
     100          18 :         if (disc.isMarkedDestroyed())
     101             :         {
     102             :             // This function also removes destroyed discs and thus needs to be called first
     103           3 :             discs_[i] = std::move(discs_.back());
     104           3 :             discs_.pop_back();
     105           3 :             --i;
     106           3 :             continue;
     107             :         }
     108             : 
     109          15 :         disc.move(disc.getVelocity() * dt);
     110             : 
     111             :         // A -> B returns nothing, A -> B + C returns 1 new disc
     112          15 :         if (auto newDisc = simulationContext_.reactionEngine.applyUnimolecularReactions(disc, dt))
     113           1 :             newDiscs.push_back(std::move(*newDisc));
     114             :     }
     115             : 
     116           9 :     if (!newDiscs.empty())
     117           1 :         discs_.insert(discs_.end(), newDiscs.begin(), newDiscs.end());
     118           9 : }
     119             : 
     120           9 : void Compartment::moveDiscsIntoChildCompartments(const CollisionDetector::DetectedCollisions& detectedCollisions)
     121             : {
     122           9 :     const auto& indexes = detectedCollisions.indexes.find(CollisionDetector::CollisionType::DiscChildMembrane);
     123           9 :     if (indexes == detectedCollisions.indexes.end())
     124           7 :         return;
     125             : 
     126           8 :     for (std::size_t index : indexes->second)
     127             :     {
     128           6 :         const auto& collision = detectedCollisions.collisions[index];
     129           6 :         if (collision.isInvalidatedByDestroyedDisc())
     130           4 :             continue;
     131             : 
     132           5 :         const auto permeability = simulationContext_.membraneTypeRegistry.getByID(collision.membrane->getTypeID())
     133           5 :                                       .getPermeabilityFor(collision.disc->getTypeID());
     134           5 :         if (permeability == MembraneType::Permeability::None || permeability == MembraneType::Permeability::Outward)
     135           3 :             continue;
     136             : 
     137           2 :         const auto& discRadius = simulationContext_.discTypeRegistry.getByID(collision.disc->getTypeID()).getRadius();
     138             :         const auto& membraneRadius =
     139           2 :             simulationContext_.membraneTypeRegistry.getByID(collision.membrane->getTypeID()).getRadius();
     140           2 :         if (mathutils::circleIsFullyContainedByCircle(collision.disc->getPosition(), discRadius,
     141           2 :                                                       collision.membrane->getPosition(), membraneRadius))
     142             :         {
     143           0 :             collision.membrane->getCompartment()->addDisc(*collision.disc);
     144           0 :             collision.disc->markDestroyed();
     145             :         }
     146             :         else
     147           2 :             collision.membrane->getCompartment()->addIntrudingDisc(*collision.disc);
     148             :     }
     149             : }
     150             : 
     151           9 : void Compartment::moveDiscsIntoParentCompartment(const CollisionDetector::DetectedCollisions& detectedCollisions)
     152             : {
     153             :     // Only the case for the containing cell
     154           9 :     if (!parent_)
     155           9 :         return;
     156             : 
     157           2 :     const auto& indexes = detectedCollisions.indexes.find(CollisionDetector::CollisionType::DiscContainingMembrane);
     158           2 :     if (indexes == detectedCollisions.indexes.end())
     159           2 :         return;
     160             : 
     161           0 :     const auto& membraneType = simulationContext_.membraneTypeRegistry.getByID(membrane_.getTypeID());
     162           0 :     const auto M = membrane_.getPosition();
     163           0 :     const auto R = membraneType.getRadius();
     164             : 
     165           0 :     for (std::size_t index : indexes->second)
     166             :     {
     167           0 :         const auto& collision = detectedCollisions.collisions[index];
     168           0 :         if (collision.isInvalidatedByDestroyedDisc())
     169           0 :             continue;
     170             : 
     171           0 :         const auto permeability = simulationContext_.membraneTypeRegistry.getByID(collision.membrane->getTypeID())
     172           0 :                                       .getPermeabilityFor(collision.disc->getTypeID());
     173           0 :         if (permeability == MembraneType::Permeability::None || permeability == MembraneType::Permeability::Inward)
     174           0 :             continue;
     175             : 
     176           0 :         const auto discRadius = simulationContext_.discTypeRegistry.getByID(collision.disc->getTypeID()).getRadius();
     177           0 :         if (!mathutils::circlesOverlap(collision.disc->getPosition(), discRadius, M, R))
     178             :         {
     179           0 :             parent_->addDisc(*collision.disc);
     180           0 :             collision.disc->markDestroyed();
     181             :         }
     182             :         else
     183           0 :             parent_->addIntrudingDisc(*collision.disc);
     184             :     }
     185             : }
     186             : 
     187           9 : void Compartment::updateChildCompartments(double dt)
     188             : {
     189          11 :     for (auto& compartment : compartments_)
     190           2 :         compartment->update(dt);
     191           9 : }
     192             : 
     193             : } // namespace cell