From 8c284a007cb49eb216e23a812d7c98b96f101325 Mon Sep 17 00:00:00 2001
From: Kyle <kyle.c.klenk@gmail.com>
Date: Mon, 26 Jun 2023 23:11:01 +0000
Subject: [PATCH] Added C++ actor model to docker file
Created basic actor code that can run s1 python file
---
.gitignore | 6 +
CMakeLists.txt | 15 +++
environment_setup/Dockerfile | 19 +++
find_s1.py | 111 ++++++++++++++++
find_s2.py | 131 +++++++++++++++++++
main.cpp | 75 +++++++++++
runS1S2.py | 158 +----------------------
run_with_actors.py | 15 ++-
basic.par => simulations/input/basic.par | 0
test.py | 13 ++
10 files changed, 383 insertions(+), 160 deletions(-)
create mode 100644 .gitignore
create mode 100644 CMakeLists.txt
create mode 100644 find_s1.py
create mode 100644 find_s2.py
create mode 100644 main.cpp
rename basic.par => simulations/input/basic.par (100%)
create mode 100644 test.py
diff --git a/.gitignore b/.gitignore
new file mode 100644
index 0000000..276ff8a
--- /dev/null
+++ b/.gitignore
@@ -0,0 +1,6 @@
+my_program
+output7.txt
+simulations/output/*
+build/
+meshes/
+output8.txt
diff --git a/CMakeLists.txt b/CMakeLists.txt
new file mode 100644
index 0000000..562f1f9
--- /dev/null
+++ b/CMakeLists.txt
@@ -0,0 +1,15 @@
+cmake_minimum_required(VERSION 3.0)
+project(my_project)
+
+find_package(CAF COMPONENTS core io REQUIRED)
+find_package (PythonLibs REQUIRED)
+include_directories (${PYTHON_INCLUDE_DIRS})
+add_executable(my_program main.cpp)
+# Add this line to set the C++ 17 standard
+set_property(TARGET my_program PROPERTY CXX_STANDARD 17)
+
+target_link_libraries(my_program CAF::core CAF::io)
+target_link_libraries(my_program ${PYTHON_LIBRARIES})
+
+# Move the executable up one directory
+set_target_properties(my_program PROPERTIES RUNTIME_OUTPUT_DIRECTORY ..)
diff --git a/environment_setup/Dockerfile b/environment_setup/Dockerfile
index fc3b33c..70c80fe 100755
--- a/environment_setup/Dockerfile
+++ b/environment_setup/Dockerfile
@@ -18,6 +18,15 @@ RUN apt-get update -y && apt-get upgrade -y && \
python3-testresources python-is-python3 \
bc
+# Install network tools - (ping, traceroute, etc.)
+RUN apt-get -y install iputils-ping libssl-dev
+
+
+# Install performance tools - (top, htop, etc.)
+RUN apt-get -y install htop
+
+# Install development tools - (git, gcc, etc.)
+RUN apt-get -y install gcc build-essential
# Install PETSc
WORKDIR /opt
@@ -41,6 +50,16 @@ ENV PETSC_DIR=/usr/local/petsc
ENV PATH=$PATH:$PETSC_DIR/bin
ENV PATH=$PATH:/root/.local/bin
+# Install common applications - (C++ Actor Framework)
+WORKDIR /opt
+WORKDIR /opt/c++_actor_framework
+
+RUN git clone https://github.com/actor-framework/actor-framework.git
+WORKDIR /opt/c++_actor_framework/actor-framework/build
+RUN cmake ..
+RUN make
+RUN make install
+
# Install OpenCARP
WORKDIR /opt
RUN git clone https://git.opencarp.org/openCARP/openCARP.git
diff --git a/find_s1.py b/find_s1.py
new file mode 100644
index 0000000..95b5eff
--- /dev/null
+++ b/find_s1.py
@@ -0,0 +1,111 @@
+#!/usr/bin/env python
+
+import os
+import sys
+NAME = 'Find_S1'
+AUTHORS = 'Joyce Reimer, Kyle Klenk, Raymond Spiteri'
+GUIinclude = True
+
+from datetime import date
+from carputils import settings
+from carputils import tools
+from carputils import mesh
+from carputils.carpio import txt
+import numpy as np
+from numpy import array as nplist
+
+def parser_for_s1():
+ parser = tools.standard_parser()
+ group = parser.add_argument_group('Experiment for finding S1 threshold')
+ group.add_argument('--duration',
+ type = float,
+ default = 20.,
+ help = 'Duration of simulation in [ms] (default: 20.)')
+ group.add_argument('--S1-strength',
+ type = float,
+ default = 20.,
+ help = 'pick transmembrane current stimulus strength in [uA/cm^2] = [pA/pF] considering fixed Cm (default: 20.)')
+ group.add_argument('--S1-dur',
+ type = float,
+ default = 2.,
+ help = 'pick transmembrane current stimulus duration in [ms] (default: 2.)')
+ return parser
+
+def jobID(args):
+ """
+ Generate name of top level output directory.
+ """
+ today = date.today()
+ return 'simulations/output/s1_{}_basic_{}'.format(today.isoformat(), args.S1_strength)
+
+@tools.carpexample(parser_for_s1, jobID)
+def run(args, job):
+ #===========================================
+ # 1: Define the mesh: Units are mm
+ #===========================================
+
+ # Create a block with a thin z dimension
+ x = 10.0 #(mm)
+ y = 1.0 #(mm)
+ z = 0.2 #(mm)
+ geom = mesh.Block(size=(x,y,z))
+
+ # Set fibre angle to 0, sheet angle to 0
+ geom.set_fibres(0, 0, 90, 90) #endo, epi, sheet endo, sheet epicard
+
+ # Finilaize the creation of the mesh
+ meshname = mesh.generate(geom)
+
+ #===========================================
+ # 2: Define the ionic models & conductivities
+ #===========================================
+
+ # Query for element labels
+ _, etags, _ = txt.read(meshname + '.elem')
+ etags = np.unique(etags)
+
+ IntraTags = etags[etags != 0].tolist()
+ ExtraTags = etags.tolist().copy()
+
+ #===========================================
+ # 3: Configure the simulation
+ #===========================================
+
+ # Define the geometry of the stimulus at one end of the block
+ # Units are um
+ stim = ['-num_stim', 1,
+ '-stimulus[0].name', 'S1',
+ '-stimulus[0].stimtype', 0,
+ '-stimulus[0].strength', args.S1_strength,
+ '-stimulus[0].duration', args.S1_dur ]
+
+ #electrode = mesh.block_boundary_condition(geom, 'stimulus', 0, 'x', True) #True refers to a BC between tissue and bath
+ E1_lower_bound = nplist([ 0.0, 0.0, 0.0])
+ E1_upper_bound = nplist([ 2.0, 0.3, 0.2])
+ electrode = mesh.block_region(geom, 'stimulus', 0, E1_lower_bound, E1_upper_bound, False) # introudce a stimulus in the region defined by the lower and upper bounds
+
+ #electrode = mesh.block_boundary_condition(geom, 'stimulus', 0, 'x', True) #True refers to a BC between tissue and bath
+ E1_lower_bound = nplist([ 0.0, 0.0, 0.0])
+ E1_upper_bound = nplist([ 2.0, 0.3, 0.2])
+ electrode = mesh.block_region(geom, 'stimulus', 0, E1_lower_bound, E1_upper_bound, False) # introudce a stimulus in the region defined by the lower and upper bounds
+
+ #===========================================
+ # 4: Defining simulator options:
+ #===========================================
+
+ # Get basic command line, including solver options
+ cmd = tools.carp_cmd(os.path.join('simulations/input/basic.par')) #cmd is a carp command
+
+ cmd += ['-simID', job.ID,
+ '-meshname', meshname,
+ '-dt', 25,
+ '-tend', args.duration]
+ cmd += tools.gen_physics_opts(ExtraTags=ExtraTags, IntraTags=IntraTags)
+ cmd += stim + electrode
+
+ # Run simulation
+ job.carp(cmd) #cmd includes the params
+
+
+
+run()
diff --git a/find_s2.py b/find_s2.py
new file mode 100644
index 0000000..c522555
--- /dev/null
+++ b/find_s2.py
@@ -0,0 +1,131 @@
+#!/usr/bin/env python
+
+import os
+NAME = 'Find_S2'
+AUTHORS = 'Joyce Reimer, Kyle Klenk, Raymond Spiteri'
+GUIinclude = True
+
+from datetime import date
+from carputils import settings
+from carputils import tools
+from carputils import mesh
+from carputils.carpio import txt
+import numpy as np
+from numpy import array as nplist
+
+def parser_for_s2(s1_strength):
+ parser = tools.standard_parser()
+ group = parser.add_argument_group('experiment specific options')
+ group.add_argument('--duration',
+ type = float,
+ default = 300.,
+ help = 'Duration of simulation in [ms] (default: 20.)')
+ group.add_argument('--S1-strength',
+ type = float,
+ default = s1_strength,
+ help = 'pick S1 transmembrane current stimulus strength in [uA/cm^2] = [pA/pF] considering fixed Cm (default: 20.)')
+ group.add_argument('--S1-dur',
+ type = float,
+ default = 2.,
+ help = 'pick transmembrane current stimulus duration in [ms] (default: 2.)')
+ group.add_argument('--S1-start',
+ type = float,
+ default = 0.0,
+ help = 'pick S1 start time')
+ group.add_argument('--S2-strength',
+ type = float,
+ default = 200.,
+ help = 'pick S2 transmembrane current stimulus strength in [uA/cm^2] = [pA/pF] considering fixed Cm (default: 20.)')
+ group.add_argument('--S2-dur',
+ type = float,
+ default = 2.,
+ help = 'pick transmembrane current stimulus duration in [ms] (default: 2.)')
+ group.add_argument('--S2-start',
+ type = float,
+ default = 150.0,
+ help = 'pick S2 start time')
+ return parser
+
+
+def jobID(args):
+ """
+ Generate name of top level output directory.
+ """
+ today = date.today()
+ return 'output/s2_{}_basic_{}'.format(today.isoformat(), args.S2_strength)
+
+@tools.carpexample(parser_for_s2, jobID)
+def run(args, job):
+ #===========================================
+ # 1: Define the mesh: Units are mm
+ #===========================================
+
+ # Create a block with a thin z dimension
+ x = 10.0 #(mm)
+ y = 1.0 #(mm)
+ z = 0.2 #(mm)
+ geom = mesh.Block(size=(x,y,z))
+
+ # Set fibre angle to 0, sheet angle to 0
+ geom.set_fibres(0, 0, 90, 90) #endo, epi, sheet endo, sheet epicard
+
+ # Finilaize the creation of the mesh
+ meshname = mesh.generate(geom)
+
+ #===========================================
+ # 2: Define the ionic models & conductivities
+ #===========================================
+
+ # Query for element labels
+ _, etags, _ = txt.read(meshname + '.elem')
+ etags = np.unique(etags)
+
+ IntraTags = etags[etags != 0].tolist()
+ ExtraTags = etags.tolist().copy()
+
+ #===========================================
+ # 3: Configure the simulation
+ #===========================================
+
+ # Define the geometry of the stimulus at one end of the block
+ # Units are um
+ stim = ['-num_stim', 2,
+ '-stimulus[0].name', 'S1', #index 0 because only one type of electrode
+ '-stimulus[0].stimtype', 0, #0 means intracellular stim
+ '-stimulus[0].strength', args.S1_strength, #uA/Cm^2
+ '-stimulus[0].duration', args.S1_dur,
+ '-stimulus[0].start', args.S1_start,
+ '-stimulus[1].name', 'S2',
+ '-stimulus[1].stimtype', 0,
+ '-stimulus[1].strength', args.S2_strength,
+ '-stimulus[1].duration', args.S2_dur,
+ '-stimulus[1].start', args.S2_start]
+
+ #electrode = mesh.block_boundary_condition(geom, 'stimulus', 0, 'x', True) #True refers to a BC between tissue and bath
+ E1_lower_bound = nplist([ 0.0, 0.0, 0.0])
+ E1_upper_bound = nplist([ 2.0, 0.3, 0.2])
+ electrode = mesh.block_region(geom, 'stimulus', 0, E1_lower_bound, E1_upper_bound, False) # introudce a stimulus in the region defined by the lower and upper bounds
+
+ #electrode = mesh.block_boundary_condition(geom, 'stimulus', 0, 'x', True) #True refers to a BC between tissue and bath
+ E1_lower_bound = nplist([ 0.0, 0.0, 0.0])
+ E1_upper_bound = nplist([ 2.0, 0.3, 0.2])
+ electrode = mesh.block_region(geom, 'stimulus', 0, E1_lower_bound, E1_upper_bound, False) # introudce a stimulus in the region defined by the lower and upper bounds
+
+ #===========================================
+ # 4: Defining simulator options:
+ #===========================================
+
+ # Get basic command line, including solver options
+ cmd = tools.carp_cmd(os.path.join('simulations/input/basic.par')) #cmd is a carp command
+
+ cmd += ['-simID', job.ID,
+ '-meshname', meshname,
+ '-dt', 25,
+ '-tend', args.duration]
+ cmd += tools.gen_physics_opts(ExtraTags=ExtraTags, IntraTags=IntraTags)
+ cmd += stim + electrode
+
+ # Run simulation
+ job.carp(cmd) #cmd includes the params
+
+run()
\ No newline at end of file
diff --git a/main.cpp b/main.cpp
new file mode 100644
index 0000000..dc5d3fe
--- /dev/null
+++ b/main.cpp
@@ -0,0 +1,75 @@
+#include <iostream>
+#include <fstream>
+#include <chrono>
+#include <thread>
+#include <caf/all.hpp>
+#include "caf/io/all.hpp"
+#include <cstdio>
+#include <sstream>
+
+
+
+
+using namespace caf;
+
+
+behavior mirror(event_based_actor* self) {
+
+
+ return {
+ [=](float duration, float s1_strength, float s1_duration) {
+ /*
+ Run the S1 Script
+ */
+ std::stringstream ss;
+ ss << "python3 /code/find_s1.py"
+ << " --duration " << duration
+ << " --S1-strength " << s1_strength
+ << " --S1-dur " << s1_duration;
+
+ FILE* pipe = popen(ss.str().c_str(), "r"); // create a pipe and run the
+
+ if (!pipe) std::cerr << "popen() failed!" << std::endl;
+ char buffer[128];
+ std::string result = "";
+ while (!feof(pipe)) {
+ if (fgets(buffer, 128, pipe) != NULL)
+ result += buffer;
+ }
+ pclose(pipe); // close the pipe
+
+ // Log output to a txt file
+ std::stringstream ss_out;
+ ss_out << "/code/output" << self->id() << ".txt";
+
+ std::ofstream out(ss_out.str().c_str());
+ out << result;
+ out.close();
+ }
+ };
+}
+
+
+
+
+void caf_main(actor_system& system) {
+ scoped_actor self{system};
+
+
+
+ auto actor_1 = self->spawn(mirror);
+ auto actor_2 = self->spawn(mirror);
+
+ float duration = 20.0;
+ float s1_strength_min = 100.0;
+ float s1_duration = 2.0;
+ float s1_strength_max = 200.0;
+
+
+ self->send(actor_1, duration, s1_strength_min, s1_duration);
+ self->send(actor_2, duration, s1_strength_max, s1_duration);
+
+}
+
+
+CAF_MAIN(io::middleman)
diff --git a/runS1S2.py b/runS1S2.py
index 01c5e89..d4d3874 100755
--- a/runS1S2.py
+++ b/runS1S2.py
@@ -1,154 +1,8 @@
#!/usr/bin/env python
-"""
-Overview
-========
-
-.. code-block:: bash
-
- cd ${TUTORIALS}/02_EP_tissue/01_basic_usage
-
-This example introduces to the basics of using the openCARP executable
-for simulating EP at the tissue and organ scale.
-Details on how to inquire input parameters are provided in the :ref:`openCARP <openCARP>` section
-of the manual.
-
-In absence of pacemaker currents electrical stimulus currents must be supplied in order to
-depolarize tissue above the threshold to elicit a propagated action potential.
-The simplest way of initiating a propagating action potential is the use of a transmembrane stimulus, :math:`I_{\mathrm{tr}}`.
-According to the definition in openCARP supplying a positive :math:`I_{\mathrm{tr}}` depolarizes the tissue,
-a negative :math:`I_{\mathrm{tr}}` hyperpolarizes the tissue.
-
-
-Experimental setup
-==================
-
-To introduce definition and usage of transmembrane stimulus currents in openCARP a simple experiment
-using a thin monolayer preparation was defined. The following input parameters are available to steer the experiment:
-
-.. code-block:: bash
-
- --duration
- Duration of simulation (ms)
-
- --S1-strength
- pick transmembrane current stimulus strength in
- [uA/cm^2] = [pA/pF] considering fixed Cm (default: 20.)
-
- --S1-dur
- pick transmembrane current stimulus duration in [ms
-
-
-.. _tutorial-subthreshold-stimuation:
-
-Modeling subthreshold conditions
---------------------------------
-
-.. _exp01-ep-ti-basic:
-
-**Experiment exp01 (subthreshold depolarizing stimulus)**
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-
-We start with a longer stimulus lasting for 15 ms and gradually increase the stimulus strength
-until we achieve capture. This is equivalent to measuring one data point of a `strength-duration curve <https://en.wikipedia.org/wiki/Rheobase>`_
-to determine chronaxie and rheobase of our preparation.
-
-.. code-block:: bash
-
- ./run.py --duration 20 --S1-strength 20. --S1-dur 15 --visualize
-
-This weak current slowly depolarizes the membrane, but does not reach threshold.
-After the end of stimulation, the membrane returns to its resting state, as shown in :numref:`fig-tutorial-subthreshold-depol-stim`.
-
-.. _fig-tutorial-subthreshold-depol-stim:
-
-.. figure:: /images/02_01_subthreshold_stim_response.png
- :width: 50%
- :align: center
-
- Membrane response to a subthresold *depolarizing* transmembrane stimulus.
-
-Inspect the variation in membrane depolarization as a function of distance to the stimulus location
-by varying the vertex index in meshalyzer.
-
-To understand how the stimulus definition is passed on to the simulator we inspect the command line generated by the script.
-For this sake run
-
-.. code-block:: bash
-
- ./run.py --duration 20 --S1-strength 20. --S1-dur 15 --dry
-
-The stimulus section of the command line is then
-
-.. code-block:: bash
-
- -num_stim 1 \ # use one stimulus electrode
- -stimulus[0].name S1 \
- \ # pulse and type definition
- -stimulus[0].stimtype 0 \
- -stimulus[0].strength 20.0 \
- -stimulus[0].duration 15.0 \
- \ # electrode definition (block)
- -stimulus[0].x0 -5050.0 \
- -stimulus[0].xd 100.0 \
- -stimulus[0].y0 -550.0 \
- -stimulus[0].yd 1100.0 \
- -stimulus[0].z0 -150.0 \
- -stimulus[0].zd 300.0
-
-
-
-**Experiment exp02 (subthreshold hyperpolarizing stimulus)**
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-
-We repeat exp01 with a negative transmembrane stimulus current.
-
-.. code-block:: bash
-
- ./run.py --duration 20 --S1-strength -20 --S1-dur 15 --visualize
-
-
-As the membrane behaves essentially only passively in the subthreshold regime
-we observe the same response symmetric to the resting potential, see :numref:`fig-tutorial-subthreshold-hyper-stim`
-
-.. _fig-tutorial-subthreshold-hyper-stim:
-
-.. figure:: /images/02_01_subthreshold_stim_response_hyper.png
- :width: 50%
- :align: center
-
- Membrane response to a subthresold *hyperpolarizing* transmembrane stimulus.
-
-
-.. _tutorial-suprathreshold-stimuation:
-
-Initiating action potential propagation
-^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
-
-We repeat exp01 with a stronger transmembrane stimulus current of shorter duration.
-Peak sodium current is known to be around :math:`250 \mu A/cm^2` and last for less than a millisecond.
-We chose a strength of :math:`250 \mu A/cm^2` of 1 millisecond duration to achieve suprathreshold depolarization.
-Keep in mind that *suprathreshold depolarization* at the stimulus site does not automatically lead to a propagated response.
-The local membrane response at the stimulus location is shown in :numref:`fig-tutorial-suprathreshold-stim`.
-
-.. code-block:: bash
-
- ./run.py --duration 20 --S1-strength 250. --S1-dur 2 --visualize
-
-
-.. _fig-tutorial-suprathreshold-stim:
-
-.. figure:: /images/02_01_suprathreshold_stim_response.png
- :width: 50%
- :align: center
-
- Membrane response to a suprathresold *depolarizing* transmembrane stimulus.
-
-"""
import os
-EXAMPLE_DESCRIPTIVE_NAME = 'Basics of usage '
-EXAMPLE_AUTHOR = 'Gernot Plank <gernot.plank@medunigraz.at>'
-EXAMPLE_DIR = os.path.dirname(__file__)
+NAME = 'S1_S2_with_actors'
+AUTHORS = 'Joyce Reimer, Kyle Klenk, Raymond Spiteri'
GUIinclude = True
from datetime import date
@@ -198,11 +52,10 @@ def jobID(args):
Generate name of top level output directory.
"""
today = date.today()
- return '{}_basic_{}'.format(today.isoformat(), args.S2_strength)
+ return '{}_basic_{}/test'.format(today.isoformat(), args.S2_strength)
@tools.carpexample(parser, jobID)
def run(args, job):
-
#===========================================
# 1: Define the mesh: Units are mm
#===========================================
@@ -211,7 +64,7 @@ def run(args, job):
x = 10.0 #(mm)
y = 1.0 #(mm)
z = 0.2 #(mm)
- geom = mesh.Block(size=(x,y,z))
+ geom = mesh.Block(size=(x,y,z))
# Set fibre angle to 0, sheet angle to 0
geom.set_fibres(0, 0, 90, 90) #endo, epi, sheet endo, sheet epicard
@@ -219,7 +72,6 @@ def run(args, job):
# Finilaize the creation of the mesh
meshname = mesh.generate(geom)
-
#===========================================
# 2: Define the ionic models & conductivities
#===========================================
@@ -259,7 +111,7 @@ def run(args, job):
#===========================================
# Get basic command line, including solver options
- cmd = tools.carp_cmd(os.path.join(EXAMPLE_DIR, 'basic.par')) #cmd is a carp command
+ cmd = tools.carp_cmd(os.path.join('simulations/input/basic.par')) #cmd is a carp command
cmd += ['-simID', job.ID,
'-meshname', meshname,
diff --git a/run_with_actors.py b/run_with_actors.py
index 06d2f07..f83560a 100644
--- a/run_with_actors.py
+++ b/run_with_actors.py
@@ -1,11 +1,14 @@
from thespian.actors import *
-import runS1S2
+# import runS1S2
import datetime
import os
import re
from subprocess import Popen, PIPE
import subprocess
import numpy as np
+from carputils import tools
+
+
class output_actor(Actor):
@@ -27,7 +30,7 @@ class simulation_actor(Actor):
class Simulation_Manager(Actor):
- def __init__(self):
+ def __init__(self)
# ===========================================
# Define the initial parameters for the S1S2 protocol
# ===========================================
@@ -105,14 +108,12 @@ class Simulation_Manager(Actor):
+
def actor_system_start():
print("Starting actor system")
-
- asys = ActorSystem('multiprocTCPBase')
-
-
+ asys = ActorSystem()
first = asys.createActor(Simulation_Manager)
- asys.tell(first, "start")
+ # asys.tell(first, "start")
if __name__ == '__main__':
actor_system_start()
diff --git a/basic.par b/simulations/input/basic.par
similarity index 100%
rename from basic.par
rename to simulations/input/basic.par
diff --git a/test.py b/test.py
new file mode 100644
index 0000000..61224c8
--- /dev/null
+++ b/test.py
@@ -0,0 +1,13 @@
+from thespian.actors import *
+
+class MyActor(Actor):
+ def __init__(self, name):
+ self.name = name
+
+ def receiveMessage(self, message, sender):
+ print(f"My name is {self.name} and I received {message} from {sender}")
+
+system = ActorSystem()
+
+actor = system.createActor(MyActor.__init__("test", 5))
+system.tell(actor, "Hello")
\ No newline at end of file
--
GitLab