diff --git a/src/algorithms.hpp b/src/algorithms.hpp
index 10ba44d30a7721e2ca58cf7e98bf399d99c05f24..bc25b98b84518e0a0c444bf6163bd05b8ed37711 100644
--- a/src/algorithms.hpp
+++ b/src/algorithms.hpp
@@ -1,5 +1,7 @@
#pragma once
+#include <iostream>
+
#include <algorithm>
#include <numeric>
#include <type_traits>
@@ -85,14 +87,19 @@ namespace algorithms{
"is_minimal_degree only applies to integral types"
);
+ if(c.size()==0u) {
+ return true;
+ }
+
using value_type = typename C::value_type;
- auto i = value_type(0);
- for(auto v : c) {
- if(i++<v) {
- return false;
- }
+ if(c[0] != value_type(0)) {
+ return false;
}
- return true;
+ auto i = value_type(1);
+ auto it = std::find_if(
+ c.begin()+1, c.end(), [&i](value_type v) { return v>=(i++); }
+ );
+ return it==c.end();
}
template <typename C>
diff --git a/src/fvm.hpp b/src/fvm.hpp
index 9e84119f526c531f3e80c082ea22d8dc7edd997b..d25d4669dfa890ab1a05f5208875e83b0f827f91 100644
--- a/src/fvm.hpp
+++ b/src/fvm.hpp
@@ -149,8 +149,6 @@ fvm_cell<T, I>::fvm_cell(nest::mc::cell const& cell)
face_alpha_[i] = area_face / (c_m * r_L * c.length);
cv_capacitance_[i] = c_m;
- std::cout << "radius " << radius_center << ", c_m " << c_m << ", r_L " << r_L << ", dx " << c.length << "\n";
-
auto halflen = c.length/2;
auto al = math::area_frustrum(halflen, left(c.radius), radius_center);
diff --git a/src/mechanism.hpp b/src/mechanism.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..3432225885848ed06f1f97513772709025757a20
--- /dev/null
+++ b/src/mechanism.hpp
@@ -0,0 +1,67 @@
+#pragma once
+
+#pragma once
+
+#include <memory>
+#include <string>
+
+#include "matrix.hpp"
+#include "parameter_list.hpp"
+#include "util.hpp"
+
+namespace nest {
+namespace mc {
+
+enum class mechanismKind {point, density};
+
+template <typename T, typename I>
+class mechanism {
+ using value_type = T;
+ using size_type = I;
+
+ // define storage types
+ using vector_type = memory::HostVector<value_type>;
+ using view_type = typename vector_type::view_type;
+ using index_type = memory::HostVector<size_type>;
+ using index_view = typename index_type::view_type;
+ //using indexed_view= IndexedView<value_type, size_type, target()>;
+
+ using matrix_type = matrix<value_type, size_type>;
+
+ mechanism(matrix_type *matrix, index_view node_indices)
+ : matrix_(matrix)
+ , node_indices_(node_indices)
+ {}
+
+ std::size_t size() const
+ {
+ return node_indices_.size();
+ }
+
+ virtual void set_params(value_type t_, value_type dt_) = 0;
+ virtual std::string name() const = 0;
+ virtual std::size_t memory() const = 0;
+ virtual void init() = 0;
+ virtual void state() = 0;
+ virtual void current() = 0;
+
+ virtual mechanismKind kind() const = 0;
+
+ matrix_type* matrix_;
+ index_type node_indices_;
+};
+
+template <typename T, typename I>
+using mechanism_ptr = std::unique_ptr<mechanism<T,I>>;
+
+template <typename M>
+mechanism_ptr<typename M::value_type, typename M::size_type>
+make_mechanism(
+ typename M::matrix_type* matrix,
+ typename M::index_view node_indices
+) {
+ return util::make_unique<M>(matrix, node_indices);
+}
+
+} // namespace nest
+} // namespace mc
diff --git a/src/mechanism_interface.cpp b/src/mechanism_interface.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..ba24eb06ce13c3452bb5db6d8847e3b1b1ce3081
--- /dev/null
+++ b/src/mechanism_interface.cpp
@@ -0,0 +1,24 @@
+#include "mechanism_interface.hpp"
+
+/* include the mechanisms
+#include <mechanisms/hh.hpp>
+#include <mechanisms/pas.hpp>
+*/
+
+namespace nest {
+namespace mc {
+namespace mechanisms {
+
+std::map<std::string, mechanism_helper<double, int>> mechanism_map;
+
+void setup_mechanism_helpers() {
+ /* manually insert
+ mechanism_map["hh"] = mechanisms::hh;
+ mechanism_map["pas"] = mechanisms::pas;
+ */
+}
+
+} // namespace mechanisms
+} // namespace nest
+} // namespace mc
+
diff --git a/src/mechanism_interface.hpp b/src/mechanism_interface.hpp
new file mode 100644
index 0000000000000000000000000000000000000000..389ef05997b90c96782dff8593aa5817ca4b2fde
--- /dev/null
+++ b/src/mechanism_interface.hpp
@@ -0,0 +1,32 @@
+#pragma once
+
+#include <map>
+#include <string>
+
+#include "matrix.hpp"
+#include "mechanism.hpp"
+#include "parameter_list.hpp"
+
+namespace nest {
+namespace mc {
+
+///
+template <typename T, typename I>
+struct mechanism_helper {
+ using index_type = memory::HostView<I>;
+ using index_view = typename index_type::view_type;
+ using mechanism_type = mechanism<T, I>;
+ using matrix_type = typename mechanism_type::matrix_type;
+
+ virtual std::string name() const = 0;
+
+ virtual mechanism_type new_mechanism(matrix_type*, index_view) const = 0;
+
+ virtual void set_parameters(mechanism_type&, parameter_list const&) const = 0;
+};
+
+extern std::map<std::string, mechanism_helper<double, int>> mechanism_helpers;
+void setup_mechanism_helpers();
+
+} // namespace nest
+} // namespace mc
diff --git a/src/parameter_list.hpp b/src/parameter_list.hpp
index e7dd1bd1cb1da945d48d8aa11fca8c07674a030f..56848466f58f36171e44b2d6da2c371a4edd84a9 100644
--- a/src/parameter_list.hpp
+++ b/src/parameter_list.hpp
@@ -146,8 +146,8 @@ namespace mc {
membrane_parameters()
: base("membrane")
{
- base::add_parameter({"c_m", 0.01, {0., 1e9}}); // typically 10 nF/mm^2 == 0.01 F/m2
- base::add_parameter({"r_L", 180.00, {0., 1e9}}); // Ohm.cm
+ base::add_parameter({"r_L", 0.01, {0., 1e9}}); // typically 10 nF/mm^2 == 0.01 F/m2
+ base::add_parameter({"c_m", 180.00, {0., 1e9}}); // Ohm.cm
}
};
diff --git a/tests/test_algorithms.cpp b/tests/test_algorithms.cpp
index 9328b315f992a9b6b12dd9e88b0e13ef6ace5d22..cfe058a3c8b92d930c739fbbde3f5a3bad6f634c 100644
--- a/tests/test_algorithms.cpp
+++ b/tests/test_algorithms.cpp
@@ -47,11 +47,6 @@ TEST(algorithms, minimal_degree)
EXPECT_TRUE(nest::mc::algorithms::is_minimal_degree(v));
}
- {
- std::vector<int> v = {0, 1, 2, 3, 4};
- EXPECT_TRUE(nest::mc::algorithms::is_minimal_degree(v));
- }
-
{
std::vector<int> v = {0, 0, 1, 2, 3, 4};
EXPECT_TRUE(nest::mc::algorithms::is_minimal_degree(v));
@@ -76,6 +71,11 @@ TEST(algorithms, minimal_degree)
std::vector<int> v = {0, 2};
EXPECT_FALSE(nest::mc::algorithms::is_minimal_degree(v));
}
+
+ {
+ std::vector<int> v = {0, 1, 2};
+ EXPECT_FALSE(nest::mc::algorithms::is_minimal_degree(v));
+ }
}
TEST(algorithms, is_strictly_monotonic_increasing)