diff --git a/arbor/backends/gpu/mechanism.cpp b/arbor/backends/gpu/mechanism.cpp index f1c7422514c62ffbab34b889274cdbcdf22e105d..9d32f5ff062e718031032e572bfb9597f7de1fb9 100644 --- a/arbor/backends/gpu/mechanism.cpp +++ b/arbor/backends/gpu/mechanism.cpp @@ -94,6 +94,7 @@ void mechanism::instantiate(unsigned id, ion_view.reversal_potential = oion->eX_.data(); ion_view.internal_concentration = oion->Xi_.data(); ion_view.external_concentration = oion->Xo_.data(); + ion_view.ionic_charge = oion->charge.data(); } event_stream_ptr_ = &shared.deliverable_events; diff --git a/arbor/backends/gpu/mechanism.hpp b/arbor/backends/gpu/mechanism.hpp index 94144bf37ade1b0746229c58fbea858c327a7212..c9c242d243d3e6fecbc49d6efb8d05d332279458 100644 --- a/arbor/backends/gpu/mechanism.hpp +++ b/arbor/backends/gpu/mechanism.hpp @@ -33,8 +33,6 @@ protected: using array = arb::gpu::array; using iarray = arb::gpu::iarray; - using ion_state_vuew = arb::gpu::ion_state_view; - public: std::size_t size() const override { return width_; diff --git a/arbor/backends/gpu/mechanism_ppack_base.hpp b/arbor/backends/gpu/mechanism_ppack_base.hpp index 457f3eeb47b7be4af5db91cc6932fe317bcc99e7..3232ef4d77f0f80cd8740e648efb532fa2a23230 100644 --- a/arbor/backends/gpu/mechanism_ppack_base.hpp +++ b/arbor/backends/gpu/mechanism_ppack_base.hpp @@ -18,6 +18,7 @@ struct ion_state_view { value_type* reversal_potential; value_type* internal_concentration; value_type* external_concentration; + value_type* ionic_charge; }; // Parameter pack base: diff --git a/arbor/backends/gpu/shared_state.cpp b/arbor/backends/gpu/shared_state.cpp index 149f0bd5034740d6ae493d103eb9918d4102896b..5c29d9fb31dafc5989c376651d0cd83d699bdfdc 100644 --- a/arbor/backends/gpu/shared_state.cpp +++ b/arbor/backends/gpu/shared_state.cpp @@ -25,7 +25,7 @@ void init_concentration_impl( void nernst_impl( std::size_t n, fvm_value_type factor, - const fvm_value_type* Xi, const fvm_value_type* Xo, fvm_value_type* eX); + const fvm_value_type* charge, const fvm_value_type* Xi, const fvm_value_type* Xo, fvm_value_type* eX); void update_time_to_impl( std::size_t n, fvm_value_type* time_to, const fvm_value_type* time, @@ -68,7 +68,7 @@ ion_state::ion_state( Xo_(cv.size(), NAN), weight_Xi_(make_const_view(iconc_norm_area)), weight_Xo_(make_const_view(econc_norm_area)), - charge(info.charge), + charge(1u, info.charge), default_int_concentration(info.default_int_concentration), default_ext_concentration(info.default_ext_concentration) { @@ -91,8 +91,8 @@ void ion_state::nernst(fvm_value_type temperature_K) { // 1e3 factor required to scale from V -> mV. constexpr fvm_value_type RF = 1e3*constant::gas_constant/constant::faraday; - fvm_value_type factor = RF*temperature_K/charge; - nernst_impl(Xi_.size(), factor, Xo_.data(), Xi_.data(), eX_.data()); + fvm_value_type factor = RF*temperature_K; + nernst_impl(Xi_.size(), factor, charge.data(), Xo_.data(), Xi_.data(), eX_.data()); } void ion_state::init_concentration() { diff --git a/arbor/backends/gpu/shared_state.cu b/arbor/backends/gpu/shared_state.cu index a5901cf585a3a39c0dd1d7a5d1f1e565545a5c62..34912740283b8fa07764c257b0a5b32f9898bc44 100644 --- a/arbor/backends/gpu/shared_state.cu +++ b/arbor/backends/gpu/shared_state.cu @@ -15,11 +15,11 @@ namespace kernel { template <typename T> __global__ -void nernst_impl(unsigned n, T factor, const T* Xo, const T* Xi, T* eX) { +void nernst_impl(unsigned n, T factor, const T* charge, const T* Xo, const T* Xi, T* eX) { unsigned i = threadIdx.x+blockIdx.x*blockDim.x; if (i<n) { - eX[i] = factor*std::log(Xo[i]/Xi[i]); + eX[i] = factor*std::log(Xo[i]/Xi[i])/charge[0]; } } @@ -93,13 +93,13 @@ using impl::block_count; void nernst_impl( std::size_t n, fvm_value_type factor, - const fvm_value_type* Xo, const fvm_value_type* Xi, fvm_value_type* eX) + const fvm_value_type* charge, const fvm_value_type* Xo, const fvm_value_type* Xi, fvm_value_type* eX) { if (!n) return; constexpr int block_dim = 128; int nblock = block_count(n, block_dim); - kernel::nernst_impl<<<nblock, block_dim>>>(n, factor, Xo, Xi, eX); + kernel::nernst_impl<<<nblock, block_dim>>>(n, factor, charge, Xo, Xi, eX); } void init_concentration_impl( diff --git a/arbor/backends/gpu/shared_state.hpp b/arbor/backends/gpu/shared_state.hpp index 27174810ba3c5f94183e4647020bc7b41bb8e0b2..03c8e44a00267cc63317f3b9b2df087f26b89422 100644 --- a/arbor/backends/gpu/shared_state.hpp +++ b/arbor/backends/gpu/shared_state.hpp @@ -34,7 +34,7 @@ struct ion_state { array weight_Xi_; // (1) concentration weight internal array weight_Xo_; // (1) concentration weight external - int charge; // charge of ionic species + array charge; // charge of ionic species (global, length 1) fvm_value_type default_int_concentration; // (mM) default internal concentration fvm_value_type default_ext_concentration; // (mM) default external concentration diff --git a/arbor/backends/multicore/mechanism.cpp b/arbor/backends/multicore/mechanism.cpp index fde15d3e68911ff18c61e2ab8d19b9eab3a7aa56..ab101ab7512b3a82afa223e7462bad34577f88cf 100644 --- a/arbor/backends/multicore/mechanism.cpp +++ b/arbor/backends/multicore/mechanism.cpp @@ -95,6 +95,7 @@ void mechanism::instantiate(unsigned id, backend::shared_state& shared, const la ion_view.reversal_potential = oion->eX_.data(); ion_view.internal_concentration = oion->Xi_.data(); ion_view.external_concentration = oion->Xo_.data(); + ion_view.ionic_charge = oion->charge.data(); } event_stream_ptr_ = &shared.deliverable_events; diff --git a/arbor/backends/multicore/mechanism.hpp b/arbor/backends/multicore/mechanism.hpp index 8054c6568023a76dbbe615130626621e7a9a55a2..937432ac12ae1f653a64cb9d1a6aab17931dede4 100644 --- a/arbor/backends/multicore/mechanism.hpp +++ b/arbor/backends/multicore/mechanism.hpp @@ -39,6 +39,7 @@ protected: value_type* reversal_potential; value_type* internal_concentration; value_type* external_concentration; + value_type* ionic_charge; }; public: diff --git a/arbor/backends/multicore/shared_state.cpp b/arbor/backends/multicore/shared_state.cpp index c2917d2fd2ce867233ba8d8c286918ae4fc86cb2..73c6605569e57dd08f6a20eaa47bc03052e3dad1 100644 --- a/arbor/backends/multicore/shared_state.cpp +++ b/arbor/backends/multicore/shared_state.cpp @@ -61,7 +61,7 @@ ion_state::ion_state( Xo_(cv.size(), NAN, pad(alignment)), weight_Xi_(iconc_norm_area.begin(), iconc_norm_area.end(), pad(alignment)), weight_Xo_(econc_norm_area.begin(), econc_norm_area.end(), pad(alignment)), - charge(info.charge), + charge(1u, info.charge, pad(alignment)), default_int_concentration(info.default_int_concentration), default_ext_concentration(info.default_ext_concentration) { @@ -84,7 +84,7 @@ void ion_state::nernst(fvm_value_type temperature_K) { // 1e3 factor required to scale from V -> mV. constexpr fvm_value_type RF = 1e3*constant::gas_constant/constant::faraday; - simd_value_type factor = RF*temperature_K/charge; + simd_value_type factor = RF*temperature_K/charge[0]; for (std::size_t i=0; i<Xi_.size(); i+=simd_width) { simd_value_type xi(Xi_.data()+i); simd_value_type xo(Xo_.data()+i); diff --git a/arbor/backends/multicore/shared_state.hpp b/arbor/backends/multicore/shared_state.hpp index 95b0ab136d3d360ac35277d1d3a946968e39a332..63e8b51f05c4315d6b296a436d553db23eef899a 100644 --- a/arbor/backends/multicore/shared_state.hpp +++ b/arbor/backends/multicore/shared_state.hpp @@ -49,7 +49,7 @@ struct ion_state { array weight_Xi_; // (1) concentration weight internal array weight_Xo_; // (1) concentration weight external - int charge; // charge of ionic species + array charge; // charge of ionic species (global value, length 1) fvm_value_type default_int_concentration; // (mM) default internal concentration fvm_value_type default_ext_concentration; // (mM) default external concentration diff --git a/arbor/fvm_layout.cpp b/arbor/fvm_layout.cpp index 144aaf5b6e1a62a0a33c800947e046e94350a933..4e021b2859f2c732da8118ed9c175d955b172134 100644 --- a/arbor/fvm_layout.cpp +++ b/arbor/fvm_layout.cpp @@ -257,7 +257,7 @@ fvm_discretization fvm_discretize(const std::vector<cable_cell>& cells) { // IIb. Density mechanism CVs, parameter values; ion channel default concentration contributions. // IIc. Point mechanism CVs, parameter values, and targets. -fvm_mechanism_data fvm_build_mechanism_data(const mechanism_catalogue& catalogue, const std::vector<cable_cell>& cells, const fvm_discretization& D, bool coalesce_syanpses) { +fvm_mechanism_data fvm_build_mechanism_data(const cable_cell_global_properties& gprop, const std::vector<cable_cell>& cells, const fvm_discretization& D) { using util::assign; using util::sort_by; using util::optional; @@ -269,6 +269,8 @@ fvm_mechanism_data fvm_build_mechanism_data(const mechanism_catalogue& catalogue using string_set = std::unordered_set<std::string>; using string_index_map = std::unordered_map<std::string, size_type>; + const mechanism_catalogue& catalogue = *gprop.catalogue; + fvm_mechanism_data mechdata; // I. Collect segment mechanism info from cells. @@ -332,6 +334,29 @@ fvm_mechanism_data fvm_build_mechanism_data(const mechanism_catalogue& catalogue } }; + auto verify_ion_usage = + [&gprop](const std::string& mech_name, const mechanism_info* info) + { + const auto& global_ions = gprop.ion_default; + + for (const auto& ion: info->ions) { + const auto& ion_name = ion.first; + const auto& ion_dep = ion.second; + + if (!global_ions.count(ion_name)) { + throw cable_cell_error( + "mechanism "+mech_name+" uses ion "+ion_name+ " which is missing in global properties"); + } + + if (ion_dep.verify_ion_charge) { + if (ion_dep.expected_ion_charge!=global_ions.at(ion_name).charge) { + throw cable_cell_error( + "mechanism "+mech_name+" uses ion "+ion_name+ " expecting a different valence"); + } + } + } + }; + auto cell_segment_part = D.cell_segment_part(); size_type target_id = 0; @@ -374,6 +399,14 @@ fvm_mechanism_data fvm_build_mechanism_data(const mechanism_catalogue& catalogue } } + for (auto& entry: density_mech_table) { + verify_ion_usage(entry.first, entry.second.info); + } + + for (auto& entry: point_mech_table) { + verify_ion_usage(entry.first, entry.second.info); + } + // II. Build ion and mechanism configs. // Shared temporary lookup info across mechanism instances, set by build_param_data. @@ -574,7 +607,7 @@ fvm_mechanism_data fvm_build_mechanism_data(const mechanism_catalogue& catalogue // Generate config.param_values: contains parameters of group of synapses that can be coalesced into one instance // Generate multiplicity: contains number of synapses in each coalesced group of synapses // Generate target: contains the synapse target number - if(catalogue[name].linear && coalesce_syanpses) { + if(catalogue[name].linear && gprop.coalesce_synapses) { // cv_param_vec used to lexicographically sort the cv, parameters and target, which are stored in that order std::vector<cv_param> cv_param_vec(cv_order.size()); diff --git a/arbor/fvm_layout.hpp b/arbor/fvm_layout.hpp index e470a7190c1a776c89607a7122c5fb1c0d3d5d21..90c202203e5fa696b781ada08b999fc1f8fc25c7 100644 --- a/arbor/fvm_layout.hpp +++ b/arbor/fvm_layout.hpp @@ -140,6 +140,6 @@ struct fvm_mechanism_data { std::size_t ntarget = 0; }; -fvm_mechanism_data fvm_build_mechanism_data(const mechanism_catalogue& catalogue, const std::vector<cable_cell>& cells, const fvm_discretization& D, bool coalesce = true); +fvm_mechanism_data fvm_build_mechanism_data(const cable_cell_global_properties& gprop, const std::vector<cable_cell>& cells, const fvm_discretization& D); } // namespace arb diff --git a/arbor/fvm_lowered_cell_impl.hpp b/arbor/fvm_lowered_cell_impl.hpp index 3badb637e571c9a2d11f6193b07f2283c2887224..60394aeb06e6f2c20f21f2055dec894ab9115ba1 100644 --- a/arbor/fvm_lowered_cell_impl.hpp +++ b/arbor/fvm_lowered_cell_impl.hpp @@ -385,7 +385,7 @@ void fvm_lowered_cell_impl<B>::initialize( // Discretize mechanism data. - fvm_mechanism_data mech_data = fvm_build_mechanism_data(*catalogue, cells, D, global_props.coalesce_synapses); + fvm_mechanism_data mech_data = fvm_build_mechanism_data(global_props, cells, D); // Discritize and build gap junction info diff --git a/arbor/include/arbor/mechinfo.hpp b/arbor/include/arbor/mechinfo.hpp index 5e23e2ea38b5f9b6c8f7ebaf96ccbd0b82e19256..1badd32c887f16a6f86c5bd5f6ba7c605955444d 100644 --- a/arbor/include/arbor/mechinfo.hpp +++ b/arbor/include/arbor/mechinfo.hpp @@ -32,7 +32,15 @@ struct mechanism_field_spec { struct ion_dependency { bool write_concentration_int; bool write_concentration_ext; + + bool read_reversal_potential; bool write_reversal_potential; + + bool read_ion_charge; + + // Support for NMODL 'VALENCE n' construction. + bool verify_ion_charge; + int expected_ion_charge; }; // A hash of the mechanism dynamics description is used to ensure that offline-compiled diff --git a/mechanisms/mod/test_ca.mod b/mechanisms/mod/test_ca.mod index 09630f4377bc7df67b6e30cc73425a7a47e52f1a..f8621f2a38a51098d4a814db2d0a66a4975bebd9 100644 --- a/mechanisms/mod/test_ca.mod +++ b/mechanisms/mod/test_ca.mod @@ -2,7 +2,7 @@ NEURON { SUFFIX test_ca - USEION ca READ ica WRITE cai + USEION ca READ ica WRITE cai VALENCE 2 } UNITS { diff --git a/modcc/blocks.hpp b/modcc/blocks.hpp index 118077dfbbee7e694c0069ccf07683da881d96a8..b71284ed8f36ce0a4dc4566db0c057e383b31b01 100644 --- a/modcc/blocks.hpp +++ b/modcc/blocks.hpp @@ -14,7 +14,10 @@ struct IonDep { std::string name; // name of ion channel std::vector<Token> read; // name of channels parameters to write std::vector<Token> write; // name of channels parameters to read - std::string valence; + + Token valence_var; // optional variable name following VALENCE + int expected_valence = 0; // optional integer following VALENCE + bool has_valence_expr = false; bool has_variable(std::string const& name) const { return writes_variable(name) || reads_variable(name); @@ -41,6 +44,13 @@ struct IonDep { return writes_variable("e"+name); }; + bool uses_valence() const { + return valence_var.type==tok::identifier; + } + bool verifies_valence() const { + return has_valence_expr && !uses_valence(); + } + bool reads_variable(const std::string& name) const { return std::find_if(read.begin(), read.end(), [&name](const Token& t) {return t.spelling==name;}) != read.end(); diff --git a/modcc/identifier.hpp b/modcc/identifier.hpp index 58aea5691066793ece8e079584751dccdb290861..f01fe61f9fb45523d32308ad975f0765d5f54bc5 100644 --- a/modcc/identifier.hpp +++ b/modcc/identifier.hpp @@ -43,6 +43,7 @@ enum class sourceKind { ion_revpot, ion_iconc, ion_econc, + ion_valence, temperature, no_source }; diff --git a/modcc/module.cpp b/modcc/module.cpp index 3c86d0decd51364b5c45e40449f875f45467723c..dc7b559340853a3b460c1e23b409df648121b8e6 100644 --- a/modcc/module.cpp +++ b/modcc/module.cpp @@ -641,6 +641,12 @@ void Module::add_variables_to_symbols() { for(auto const& var : ion.write) { update_ion_symbols(var, accessKind::write, ion.name); } + + if(ion.uses_valence()) { + Token valence_var = ion.valence_var; + create_indexed_variable(valence_var.spelling, "ion_valence", sourceKind::ion_valence, + tok::eq, accessKind::read, ion.name, valence_var.location); + } } // then GLOBAL variables diff --git a/modcc/parser.cpp b/modcc/parser.cpp index d2fd3a021a18bea25d967ea7fedfea8c90069ba8..aa02f3ec829bfa2a629ff59a0fd5112a7d6585c0 100644 --- a/modcc/parser.cpp +++ b/modcc/parser.cpp @@ -1,5 +1,5 @@ -#include <list> #include <cstring> +#include <string> #include "parser.hpp" #include "perfvisitor.hpp" @@ -328,9 +328,19 @@ void Parser::parse_neuron_block() { } if(token_.type == tok::valence) { - //Consume "Valence" + ion.has_valence_expr = true; + + // consume "Valence" get_token(); - ion.valence == value_literal(); + + // take and consume variable name or signed integer + if(token_.type == tok::identifier) { + ion.valence_var = token_; + get_token(); + } + else { + ion.expected_valence = value_signed_integer(); + } } // add the ion dependency to the NEURON block @@ -634,6 +644,29 @@ std::string Parser::value_literal() { } } +// Parse an integral value with possible preceding unary plus or minus, +// and return as an int. +int Parser::value_signed_integer() { + std::string value; + + if(token_.type==tok::minus) { + value = "-"; + get_token(); + } + else if(token_.type==tok::plus) { + get_token(); + } + if(token_.type != tok::integer) { + error(pprintf("numeric constant not an integer '%'", token_)); + return 0; + } + else { + value += token_.spelling; + get_token(); + return std::stoi(value); + } +} + std::vector<Token> Parser::unit_description() { static const tok legal_tokens[] = {tok::identifier, tok::divide, tok::real, tok::integer}; int startline = location_.line; diff --git a/modcc/parser.hpp b/modcc/parser.hpp index 8b9414237db2948a899eeed990881b08d93c8e6f..8590938c25189ee40c7c04caeaeac7763ab7563d 100644 --- a/modcc/parser.hpp +++ b/modcc/parser.hpp @@ -62,6 +62,7 @@ private: std::vector<Token> comma_separated_identifiers(); std::vector<Token> unit_description(); std::string value_literal(); + int value_signed_integer(); std::pair<Token, Token> range_description(); /// build the identifier list diff --git a/modcc/printer/infoprinter.cpp b/modcc/printer/infoprinter.cpp index 5e9fe0b979142ab8ccc872cadbebb314e81b926f..1326923e8f1548de6ff6372370a906389b914de3 100644 --- a/modcc/printer/infoprinter.cpp +++ b/modcc/printer/infoprinter.cpp @@ -47,7 +47,11 @@ std::ostream& operator<<(std::ostream& out, const ion_dep_info& wrap) { return out << "{\"" << ion.name << "\", {" << boolalpha[ion.writes_concentration_int()] << ", " << boolalpha[ion.writes_concentration_ext()] << ", " - << boolalpha[ion.writes_rev_potential()] << "}}"; + << boolalpha[ion.uses_rev_potential()] << ", " + << boolalpha[ion.writes_rev_potential()] << ", " + << boolalpha[ion.uses_valence()] << ", " + << boolalpha[ion.verifies_valence()] << ", " + << ion.expected_valence << "}}"; } std::string build_info_header(const Module& m, const printer_options& opt) { @@ -58,6 +62,11 @@ std::string build_info_header(const Module& m, const printer_options& opt) { auto ids = public_variable_ids(m); auto ns_components = namespace_components(opt.cpp_namespace); + bool any_fields = + !ids.global_parameter_ids.empty() || + !ids.range_parameter_ids.empty() || + !ids.state_ids.empty(); + io::pfxstringstream out; out << @@ -73,8 +82,12 @@ std::string build_info_header(const Module& m, const printer_options& opt) { "::arb::concrete_mech_ptr<Backend> make_mechanism_" << name << "();\n" "\n" "inline const ::arb::mechanism_info& mechanism_" << name << "_info() {\n" - << indent << - "using spec = ::arb::mechanism_field_spec;\n" + << indent; + + any_fields && out << + "using spec = ::arb::mechanism_field_spec;\n"; + + out << "static mechanism_info info = {\n" << indent << "// globals\n" diff --git a/modcc/printer/printerutil.cpp b/modcc/printer/printerutil.cpp index ca42489ca04faba6547fb3e028770f21b81ab709..e3dbbcf54f2425eac7962128d04b564449ccfd26 100644 --- a/modcc/printer/printerutil.cpp +++ b/modcc/printer/printerutil.cpp @@ -147,6 +147,10 @@ indexed_variable_info decode_indexed_variable(IndexedVariable* sym) { case sourceKind::ion_econc: data_var=ion_pfx+".external_concentration"; break; + case sourceKind::ion_valence: + data_var=ion_pfx+".ionic_charge"; + index_var=""; // scalar global + break; case sourceKind::temperature: data_var="temperature_degC_"; index_var=""; // scalar global diff --git a/test/unit-modcc/test_parser.cpp b/test/unit-modcc/test_parser.cpp index 84592054094f941c2c8e080e9c9d07ccdf9bc4bb..6551eb10bc53d536b403504043c5538b3206b344 100644 --- a/test/unit-modcc/test_parser.cpp +++ b/test/unit-modcc/test_parser.cpp @@ -606,6 +606,89 @@ TEST(Parser, parse_state_block) { } } +static std::vector<IonDep> extract_useion(const std::string& neuron_block) { + Module m(neuron_block, "dummy"); + Parser p(m, false); + p.parse_neuron_block(); + EXPECT_EQ(lexerStatus::happy, p.status()); + verbose_print(expression_ptr{}, p, neuron_block.c_str()); + + return m.neuron_block().ions; +} + +TEST(Parser, parse_neuron_block_useion) { + { + const char* neuron_block = "NEURON { USEION x }"; + IonDep ion = extract_useion(neuron_block).at(0); + + EXPECT_EQ("x", ion.name); + EXPECT_EQ(true, ion.read.empty()); + EXPECT_EQ(true, ion.write.empty()); + EXPECT_EQ(false, ion.uses_valence()); + EXPECT_EQ(false, ion.verifies_valence()); + } + { + const char* neuron_block = "NEURON { USEION x READ ix, xi, xo }"; + IonDep ion = extract_useion(neuron_block).at(0); + + EXPECT_EQ("x", ion.name); + EXPECT_EQ(false, ion.read.empty()); + EXPECT_EQ(true, ion.write.empty()); + EXPECT_EQ(false, ion.uses_valence()); + EXPECT_EQ(false, ion.verifies_valence()); + + EXPECT_EQ(true, ion.uses_current()); + EXPECT_EQ(true, ion.uses_concentration_int()); + EXPECT_EQ(true, ion.uses_concentration_ext()); + EXPECT_EQ(false, ion.writes_concentration_int()); + EXPECT_EQ(false, ion.writes_concentration_ext()); + } + { + const char* neuron_block = "NEURON { USEION x WRITE xi, xo }"; + IonDep ion = extract_useion(neuron_block).at(0); + + EXPECT_EQ("x", ion.name); + EXPECT_EQ(true, ion.read.empty()); + EXPECT_EQ(false, ion.write.empty()); + EXPECT_EQ(false, ion.uses_valence()); + EXPECT_EQ(false, ion.verifies_valence()); + + EXPECT_EQ(false, ion.uses_current()); + EXPECT_EQ(true, ion.uses_concentration_int()); + EXPECT_EQ(true, ion.uses_concentration_ext()); + EXPECT_EQ(true, ion.writes_concentration_int()); + EXPECT_EQ(true, ion.writes_concentration_ext()); + } + { + const char* neuron_block = "NEURON { USEION x WRITE ex VALENCE -2}"; + IonDep ion = extract_useion(neuron_block).at(0); + + EXPECT_EQ("x", ion.name); + EXPECT_EQ(true, ion.read.empty()); + EXPECT_EQ(false, ion.write.empty()); + EXPECT_EQ(false, ion.uses_valence()); + EXPECT_EQ(true, ion.verifies_valence()); + EXPECT_EQ(-2, ion.expected_valence); + + EXPECT_EQ(false, ion.uses_current()); + EXPECT_EQ(false, ion.uses_concentration_int()); + EXPECT_EQ(false, ion.uses_concentration_ext()); + EXPECT_EQ(true, ion.uses_rev_potential()); + EXPECT_EQ(true, ion.writes_rev_potential()); + } + { + const char* neuron_block = "NEURON { USEION x WRITE ex VALENCE zx}"; + IonDep ion = extract_useion(neuron_block).at(0); + + EXPECT_EQ("x", ion.name); + EXPECT_EQ(true, ion.read.empty()); + EXPECT_EQ(false, ion.write.empty()); + EXPECT_EQ(true, ion.uses_valence()); + EXPECT_EQ(false, ion.verifies_valence()); + EXPECT_EQ("zx", ion.valence_var.spelling); + } +} + TEST(Parser, parse_kinetic) { char str[] = "KINETIC kin {\n" diff --git a/test/unit/CMakeLists.txt b/test/unit/CMakeLists.txt index 696683f1f61400f18d9f0bcbaaca5b54003d7f56..0acd83afab66b24301fd54e336c0c27e81ca8c50 100644 --- a/test/unit/CMakeLists.txt +++ b/test/unit/CMakeLists.txt @@ -2,6 +2,8 @@ set(test_mechanisms celsius_test + test_cl_valence + test_ca_read_valence ) include(${PROJECT_SOURCE_DIR}/mechanisms/BuildModules.cmake) diff --git a/test/unit/mod/test_ca_read_valence.mod b/test/unit/mod/test_ca_read_valence.mod new file mode 100644 index 0000000000000000000000000000000000000000..97f40707a6c61aeb700ae2c3e1d05033f02ce18b --- /dev/null +++ b/test/unit/mod/test_ca_read_valence.mod @@ -0,0 +1,22 @@ +: Test mechanism for checking ionic valence read + +NEURON { + SUFFIX test_ca_read_valence + USEION ca READ ica VALENCE zca +} + +PARAMETER {} + +ASSIGNED {} + +STATE { + record_zca +} + +INITIAL { + record_zca = zca +} + +BREAKPOINT { +} + diff --git a/test/unit/mod/test_cl_valence.mod b/test/unit/mod/test_cl_valence.mod new file mode 100644 index 0000000000000000000000000000000000000000..0613f5ad5329ee3eac4f99278afeff78701a29e9 --- /dev/null +++ b/test/unit/mod/test_cl_valence.mod @@ -0,0 +1,19 @@ +: Test mechanism for verifying ionic valence + +NEURON { + SUFFIX test_cl_valence + USEION cl WRITE icl VALENCE -1 +} + +PARAMETER {} + +ASSIGNED {} + +INITIAL {} + +STATE {} + +BREAKPOINT { + icl = 1.23*(v-4.56) +} + diff --git a/test/unit/test_fvm_layout.cpp b/test/unit/test_fvm_layout.cpp index 9731ff5603bb485ceac10b1ee82a14c2aba630cc..a7f471c8638de0326767cd532d4ea8c6a301d87b 100644 --- a/test/unit/test_fvm_layout.cpp +++ b/test/unit/test_fvm_layout.cpp @@ -12,6 +12,7 @@ #include "util/span.hpp" #include "common.hpp" +#include "unit_test_catalogue.hpp" #include "../common_cells.hpp" using namespace std::string_literals; @@ -292,8 +293,9 @@ TEST(fvm_layout, mech_index) { cells[1].add_synapse({2, 0.4}, "exp2syn"); cells[1].add_synapse({3, 0.4}, "expsyn"); + cable_cell_global_properties gprop; fvm_discretization D = fvm_discretize(cells); - fvm_mechanism_data M = fvm_build_mechanism_data(global_default_catalogue(), cells, D); + fvm_mechanism_data M = fvm_build_mechanism_data(gprop, cells, D); auto& hh_config = M.mechanisms.at("hh"); auto& expsyn_config = M.mechanisms.at("expsyn"); @@ -350,9 +352,14 @@ TEST(fvm_layout, coalescing_synapses) { return m; }; + cable_cell_global_properties gprop_no_coalesce; + gprop_no_coalesce.coalesce_synapses = false; + + cable_cell_global_properties gprop_coalesce; + gprop_coalesce.coalesce_synapses = true; + { cable_cell cell = make_cell_ball_and_stick(); - bool coalesce_synapses = true; // Add synapses of two varieties. cell.add_synapse({1, 0.3}, "expsyn"); @@ -361,7 +368,7 @@ TEST(fvm_layout, coalescing_synapses) { cell.add_synapse({1, 0.9}, "expsyn"); fvm_discretization D = fvm_discretize({cell}); - fvm_mechanism_data M = fvm_build_mechanism_data(global_default_catalogue(), {cell}, D, coalesce_synapses); + fvm_mechanism_data M = fvm_build_mechanism_data(gprop_coalesce, {cell}, D); auto &expsyn_config = M.mechanisms.at("expsyn"); EXPECT_EQ(ivec({1, 2, 3, 4}), expsyn_config.cv); @@ -369,7 +376,6 @@ TEST(fvm_layout, coalescing_synapses) { } { cable_cell cell = make_cell_ball_and_stick(); - bool coalesce_synapses = true; // Add synapses of two varieties. cell.add_synapse({1, 0.3}, "expsyn"); @@ -378,7 +384,7 @@ TEST(fvm_layout, coalescing_synapses) { cell.add_synapse({1, 0.9}, "exp2syn"); fvm_discretization D = fvm_discretize({cell}); - fvm_mechanism_data M = fvm_build_mechanism_data(global_default_catalogue(), {cell}, D, coalesce_synapses); + fvm_mechanism_data M = fvm_build_mechanism_data(gprop_coalesce, {cell}, D); auto &expsyn_config = M.mechanisms.at("expsyn"); EXPECT_EQ(ivec({1, 3}), expsyn_config.cv); @@ -390,7 +396,6 @@ TEST(fvm_layout, coalescing_synapses) { } { cable_cell cell = make_cell_ball_and_stick(); - bool coalesce_synapses = false; // Add synapses of two varieties. cell.add_synapse({1, 0.3}, "expsyn"); @@ -399,7 +404,7 @@ TEST(fvm_layout, coalescing_synapses) { cell.add_synapse({1, 0.9}, "expsyn"); fvm_discretization D = fvm_discretize({cell}); - fvm_mechanism_data M = fvm_build_mechanism_data(global_default_catalogue(), {cell}, D, coalesce_synapses); + fvm_mechanism_data M = fvm_build_mechanism_data(gprop_no_coalesce, {cell}, D); auto &expsyn_config = M.mechanisms.at("expsyn"); EXPECT_EQ(ivec({1, 2, 3, 4}), expsyn_config.cv); @@ -407,7 +412,6 @@ TEST(fvm_layout, coalescing_synapses) { } { cable_cell cell = make_cell_ball_and_stick(); - bool coalesce_synapses = false; // Add synapses of two varieties. cell.add_synapse({1, 0.3}, "expsyn"); @@ -416,7 +420,7 @@ TEST(fvm_layout, coalescing_synapses) { cell.add_synapse({1, 0.9}, "exp2syn"); fvm_discretization D = fvm_discretize({cell}); - fvm_mechanism_data M = fvm_build_mechanism_data(global_default_catalogue(), {cell}, D, coalesce_synapses); + fvm_mechanism_data M = fvm_build_mechanism_data(gprop_no_coalesce, {cell}, D); auto &expsyn_config = M.mechanisms.at("expsyn"); EXPECT_EQ(ivec({1, 3}), expsyn_config.cv); @@ -436,7 +440,7 @@ TEST(fvm_layout, coalescing_synapses) { cell.add_synapse({1, 0.7}, "expsyn"); fvm_discretization D = fvm_discretize({cell}); - fvm_mechanism_data M = fvm_build_mechanism_data(global_default_catalogue(), {cell}, D); + fvm_mechanism_data M = fvm_build_mechanism_data(gprop_coalesce, {cell}, D); auto &expsyn_config = M.mechanisms.at("expsyn"); EXPECT_EQ(ivec({1, 3}), expsyn_config.cv); @@ -452,7 +456,7 @@ TEST(fvm_layout, coalescing_synapses) { cell.add_synapse({1, 0.7}, syn_desc("expsyn", 0.1, 0.2)); fvm_discretization D = fvm_discretize({cell}); - fvm_mechanism_data M = fvm_build_mechanism_data(global_default_catalogue(), {cell}, D); + fvm_mechanism_data M = fvm_build_mechanism_data(gprop_coalesce, {cell}, D); auto &expsyn_config = M.mechanisms.at("expsyn"); EXPECT_EQ(ivec({1, 1, 3}), expsyn_config.cv); @@ -474,7 +478,7 @@ TEST(fvm_layout, coalescing_synapses) { cell.add_synapse({1, 0.3}, syn_desc("expsyn", 1, 2)); fvm_discretization D = fvm_discretize({cell}); - fvm_mechanism_data M = fvm_build_mechanism_data(global_default_catalogue(), {cell}, D); + fvm_mechanism_data M = fvm_build_mechanism_data(gprop_coalesce, {cell}, D); auto &expsyn_config = M.mechanisms.at("expsyn"); EXPECT_EQ(ivec({1, 1, 3, 3}), expsyn_config.cv); @@ -499,7 +503,7 @@ TEST(fvm_layout, coalescing_synapses) { cell.add_synapse({1, 0.7}, syn_desc_2("exp2syn", 2, 2)); fvm_discretization D = fvm_discretize({cell}); - fvm_mechanism_data M = fvm_build_mechanism_data(global_default_catalogue(), {cell}, D); + fvm_mechanism_data M = fvm_build_mechanism_data(gprop_coalesce, {cell}, D); auto &expsyn_config = M.mechanisms.at("expsyn"); EXPECT_EQ(ivec({1, 1}), expsyn_config.cv); @@ -544,8 +548,9 @@ TEST(fvm_layout, synapse_targets) { cells[1].add_synapse({3, 0.4}, syn_desc("expsyn", 5)); cells[1].add_synapse({3, 0.7}, syn_desc("exp2syn", 6)); + cable_cell_global_properties gprop; fvm_discretization D = fvm_discretize(cells); - fvm_mechanism_data M = fvm_build_mechanism_data(global_default_catalogue(), cells, D); + fvm_mechanism_data M = fvm_build_mechanism_data(gprop, cells, D); ASSERT_EQ(1u, M.mechanisms.count("expsyn")); ASSERT_EQ(1u, M.mechanisms.count("exp2syn")); @@ -701,8 +706,9 @@ TEST(fvm_layout, density_norm_area) { expected_gl[8] = seg3_gl; expected_gl[9] = seg3_gl; + cable_cell_global_properties gprop; fvm_discretization D = fvm_discretize(cells); - fvm_mechanism_data M = fvm_build_mechanism_data(global_default_catalogue(), cells, D); + fvm_mechanism_data M = fvm_build_mechanism_data(gprop, cells, D); // Check CV area assumptions. // Note: area integrator used here and in `fvm_multicell` may differ, and so areas computed may @@ -732,6 +738,27 @@ TEST(fvm_layout, density_norm_area) { EXPECT_TRUE(testing::seq_almost_eq<double>(expected_gl, gl)); } +TEST(fvm_layout, valence_verify) { + std::vector<cable_cell> cells(1); + cable_cell& c = cells[0]; + auto soma = c.add_soma(6); + + soma->add_mechanism("test_cl_valence"); + fvm_discretization D = fvm_discretize(cells); + + mechanism_catalogue testcat = make_unit_test_catalogue(); + cable_cell_global_properties gprop; + gprop.catalogue = &testcat; + + EXPECT_THROW(fvm_build_mechanism_data(gprop, cells, D), cable_cell_error); + + gprop.ion_default["cl"] = { -2, 1., 1. }; + EXPECT_THROW(fvm_build_mechanism_data(gprop, cells, D), cable_cell_error); + + gprop.ion_default["cl"] = { -1, 1., 1. }; + EXPECT_NO_THROW(fvm_build_mechanism_data(gprop, cells, D)); +} + TEST(fvm_layout, ion_weights) { // Create a cell with 4 segments: // - Soma (segment 0) plus three dendrites (1, 2, 3) meeting at a branch point. @@ -783,6 +810,8 @@ TEST(fvm_layout, ion_weights) { {1./3}, {1./3, 1./2, 0.}, {1./4, 0., 0.}, {0., 0., 0., 0.}, {3./4, 0.} }; + cable_cell_global_properties gprop; + for (auto run: count_along(mech_segs)) { std::vector<cable_cell> cells(1); cable_cell& c = cells[0]; @@ -793,7 +822,7 @@ TEST(fvm_layout, ion_weights) { } fvm_discretization D = fvm_discretize(cells); - fvm_mechanism_data M = fvm_build_mechanism_data(global_default_catalogue(), cells, D); + fvm_mechanism_data M = fvm_build_mechanism_data(gprop, cells, D); ASSERT_EQ(1u, M.ions.count("ca"s)); auto& ca = M.ions.at("ca"s); diff --git a/test/unit/test_fvm_lowered.cpp b/test/unit/test_fvm_lowered.cpp index 34adc17e0b2c4f75bd3f288bf13adc4dfa25d658..011939b08ad0894c31ba11205b19dfa0278f7e4b 100644 --- a/test/unit/test_fvm_lowered.cpp +++ b/test/unit/test_fvm_lowered.cpp @@ -27,6 +27,7 @@ #include "util/rangeutil.hpp" #include "common.hpp" +#include "unit_test_catalogue.hpp" #include "../common_cells.hpp" #include "../simple_recipes.hpp" @@ -486,6 +487,35 @@ TEST(fvm_lowered, derived_mechs) { } } +// Test that ion charge is propagated into mechanism variable. + +TEST(fvm_lowered, read_valence) { + std::vector<cable_cell> cells(1); + + cable_cell& c = cells[0]; + c.add_soma(6.0)->add_mechanism("test_ca_read_valence"); + + cable1d_recipe rec(cells); + rec.catalogue() = make_unit_test_catalogue(); + + std::vector<target_handle> targets; + std::vector<fvm_index_type> cell_to_intdom; + probe_association_map<probe_handle> probe_map; + + execution_context context; + fvm_cell fvcell(context); + fvcell.initialize({0}, rec, cell_to_intdom, targets, probe_map); + + // test_ca_read_valence initialization should write ca ion valence + // to state variable 'record_zca': + + auto mech_ptr = dynamic_cast<multicore::mechanism*>(find_mechanism(fvcell, "test_ca_read_valence")); + auto opt_record_zca_ptr = util::value_by_key((mech_ptr->*private_field_table_ptr)(), "record_zca"s); + + ASSERT_TRUE(opt_record_zca_ptr); + auto& record_zca = *opt_record_zca_ptr.value(); + ASSERT_EQ(2.0, record_zca[0]); +} // Test area-weighted linear combination of ion species concentrations diff --git a/test/unit/unit_test_catalogue.cpp b/test/unit/unit_test_catalogue.cpp index 52a2f72c50adea56d7a97e9d7086b2d0ecfa3259..fd3af1db45d0b9e09398a8534f93e21e7c5292ac 100644 --- a/test/unit/unit_test_catalogue.cpp +++ b/test/unit/unit_test_catalogue.cpp @@ -8,20 +8,30 @@ #include "unit_test_catalogue.hpp" #include "mechanisms/celsius_test.hpp" +#include "mechanisms/test_cl_valence.hpp" +#include "mechanisms/test_ca_read_valence.hpp" #include "../gtest.h" +#ifndef ARB_GPU_ENABLED +#define ADD_MECH(c, x)\ +c.add(#x, mechanism_##x##_info());\ +c.register_implementation(#x, make_mechanism_##x<multicore::backend>()); +#else +#define ADD_MECH(c, x)\ +c.add(#x, mechanism_##x##_info());\ +c.register_implementation(#x, make_mechanism_##x<multicore::backend>());\ +c.register_implementation(#x, make_mechanism_##x<gpu::backend>()); +#endif + using namespace arb; mechanism_catalogue make_unit_test_catalogue() { mechanism_catalogue cat; - cat.add("celsius_test", mechanism_celsius_test_info()); - - cat.register_implementation("celsius_test", make_mechanism_celsius_test<multicore::backend>()); -#ifdef ARB_GPU_ENABLED - cat.register_implementation("celsius_test", make_mechanism_celsius_test<gpu::backend>()); -#endif + ADD_MECH(cat, celsius_test) + ADD_MECH(cat, test_cl_valence) + ADD_MECH(cat, test_ca_read_valence) return cat; }