diff --git a/arbor/cable_cell.cpp b/arbor/cable_cell.cpp index 2d643732fd54e0f1c6bcf72a723627dcea4f4504..ca7cc23b15e14bd734271f00a283076925a9abc9 100644 --- a/arbor/cable_cell.cpp +++ b/arbor/cable_cell.cpp @@ -88,8 +88,16 @@ struct cable_cell_impl { return region_map.get<mechanism_desc>()[desc.name()]; } - mcable_map<initial_ion_data>& get_region_map(const initial_ion_data& init) { - return region_map.get<initial_ion_data>()[init.ion]; + mcable_map<init_int_concentration>& get_region_map(const init_int_concentration& init) { + return region_map.get<init_int_concentration>()[init.ion]; + } + + mcable_map<init_ext_concentration>& get_region_map(const init_ext_concentration& init) { + return region_map.get<init_ext_concentration>()[init.ion]; + } + + mcable_map<init_reversal_potential>& get_region_map(const init_reversal_potential& init) { + return region_map.get<init_reversal_potential>()[init.ion]; } template <typename Property> @@ -168,7 +176,8 @@ void cable_cell::paint(const region& target, proptype prop) {\ } ARB_PP_FOREACH(FWD_PAINT,\ mechanism_desc, init_membrane_potential, axial_resistivity,\ - temperature_K, membrane_capacitance, initial_ion_data) + temperature_K, membrane_capacitance, init_int_concentration, + init_ext_concentration, init_reversal_potential) // Forward place methods to implementation class. diff --git a/arbor/cable_cell_param.cpp b/arbor/cable_cell_param.cpp index b20a93a1b9067e3d0ef36106603380d8ff05f9a3..a83d5d31db8dbfbc7dad24de8fd21089de4a4777 100644 --- a/arbor/cable_cell_param.cpp +++ b/arbor/cable_cell_param.cpp @@ -38,13 +38,13 @@ void check_global_properties(const cable_cell_global_properties& G) { for (const auto& kv: param.ion_data) { auto& ion = kv.first; const cable_cell_ion_data& data = kv.second; - if (std::isnan(data.init_int_concentration)) { + if (!data.init_int_concentration) { throw cable_cell_error("missing init_int_concentration for ion "+ion); } - if (std::isnan(data.init_ext_concentration)) { + if (!data.init_ext_concentration) { throw cable_cell_error("missing init_ext_concentration for ion "+ion); } - if (std::isnan(data.init_reversal_potential) && !param.reversal_potential_method.count(ion)) { + if (!data.init_reversal_potential && !param.reversal_potential_method.count(ion)) { throw cable_cell_error("missing init_reversal_potential or reversal_potential_method for ion "+ion); } } diff --git a/arbor/fvm_layout.cpp b/arbor/fvm_layout.cpp index e370372a7cf4ce90337ac60f6fcb01b985a6c465..451002a97c27341fe7279c6ab3fd8d81031ec0bf 100644 --- a/arbor/fvm_layout.cpp +++ b/arbor/fvm_layout.cpp @@ -1096,7 +1096,9 @@ fvm_mechanism_data fvm_build_mechanism_data(const cable_cell_global_properties& // Ions: - auto initial_ion_data_map = cell.region_assignments().get<initial_ion_data>(); + auto initial_iconc_map = cell.region_assignments().get<init_int_concentration>(); + auto initial_econc_map = cell.region_assignments().get<init_ext_concentration>(); + auto initial_rvpot_map = cell.region_assignments().get<init_reversal_potential>(); for (const auto& ion_cvs: ion_support) { const std::string& ion = ion_cvs.first; @@ -1111,12 +1113,20 @@ fvm_mechanism_data fvm_build_mechanism_data(const cable_cell_global_properties& config.reset_econc.resize(n_cv); config.init_revpot.resize(n_cv); - cable_cell_ion_data ion_data = *(value_by_key(dflt.ion_data, ion) | value_by_key(global_dflt.ion_data, ion)); - double dflt_iconc = ion_data.init_int_concentration; - double dflt_econc = ion_data.init_ext_concentration; - double dflt_rvpot = ion_data.init_reversal_potential; + auto global_ion_data = value_by_key(global_dflt.ion_data, ion).value(); + auto dflt_iconc = global_ion_data.init_int_concentration.value(); + auto dflt_econc = global_ion_data.init_ext_concentration.value(); + auto dflt_rvpot = global_ion_data.init_reversal_potential.value(); - const mcable_map<initial_ion_data>& ion_on_cable = initial_ion_data_map[ion]; + if (auto ion_data = value_by_key(dflt.ion_data, ion)) { + dflt_iconc = ion_data.value().init_int_concentration.value_or(dflt_iconc); + dflt_econc = ion_data.value().init_ext_concentration.value_or(dflt_econc); + dflt_rvpot = ion_data.value().init_reversal_potential.value_or(dflt_rvpot); + } + + const mcable_map<init_int_concentration>& iconc_on_cable = initial_iconc_map[ion]; + const mcable_map<init_ext_concentration>& econc_on_cable = initial_econc_map[ion]; + const mcable_map<init_reversal_potential>& rvpot_on_cable = initial_rvpot_map[ion]; auto pw_times = [](const pw_elements<double>& a, const pw_elements<double>& b) { return zip(a, b, [](double left, double right, pw_element<double> a, pw_element<double> b) { return a.second*b.second; }); @@ -1127,9 +1137,9 @@ fvm_mechanism_data fvm_build_mechanism_data(const cable_cell_global_properties& if (D.cv_area[cv]==0) continue; for (mcable c: D.geometry.cables(cv)) { - auto iconc = pw_over_cable(ion_on_cable, c, dflt_iconc, [](auto x) { return x.initial.init_int_concentration; }); - auto econc = pw_over_cable(ion_on_cable, c, dflt_econc, [](auto x) { return x.initial.init_ext_concentration; }); - auto rvpot = pw_over_cable(ion_on_cable, c, dflt_rvpot, [](auto x) { return x.initial.init_reversal_potential; }); + auto iconc = pw_over_cable(iconc_on_cable, c, dflt_iconc); + auto econc = pw_over_cable(econc_on_cable, c, dflt_econc); + auto rvpot = pw_over_cable(rvpot_on_cable, c, dflt_rvpot); config.reset_iconc[i] += embedding.integrate_area(c.branch, iconc); config.reset_econc[i] += embedding.integrate_area(c.branch, econc); diff --git a/arbor/include/arbor/cable_cell.hpp b/arbor/include/arbor/cable_cell.hpp index 574750a60a20f14feae01411793b20e904b9a912..70b42e0e8bf16c526267e3d18de0d3ddfb184b10 100644 --- a/arbor/include/arbor/cable_cell.hpp +++ b/arbor/include/arbor/cable_cell.hpp @@ -191,7 +191,8 @@ struct cable_cell_impl; template <typename T> using region_assignment = std::conditional_t< - std::is_same<T, mechanism_desc>::value || std::is_same<T, initial_ion_data>::value, + std::is_same<T, mechanism_desc>::value || std::is_same<T, init_int_concentration>::value || + std::is_same<T, init_ext_concentration>::value || std::is_same<T, init_reversal_potential>::value, std::unordered_map<std::string, mcable_map<T>>, mcable_map<T>>; @@ -215,7 +216,8 @@ using location_assignment = using cable_cell_region_map = static_typed_map<region_assignment, mechanism_desc, init_membrane_potential, axial_resistivity, - temperature_K, membrane_capacitance, initial_ion_data>; + temperature_K, membrane_capacitance, init_int_concentration, + init_ext_concentration, init_reversal_potential>; using cable_cell_location_map = static_typed_map<location_assignment, mechanism_desc, i_clamp, gap_junction_site, threshold_detector>; @@ -275,6 +277,18 @@ public: default_parameters.ion_data[prop.ion] = prop.initial; } + void set_default(init_int_concentration prop) { + default_parameters.ion_data[prop.ion].init_int_concentration = prop.value; + } + + void set_default(init_ext_concentration prop) { + default_parameters.ion_data[prop.ion].init_ext_concentration = prop.value; + } + + void set_default(init_reversal_potential prop) { + default_parameters.ion_data[prop.ion].init_reversal_potential = prop.value; + } + void set_default(ion_reversal_potential_method prop) { default_parameters.reversal_potential_method[prop.ion] = prop.method; } @@ -292,7 +306,9 @@ public: void paint(const region&, axial_resistivity); void paint(const region&, temperature_K); void paint(const region&, membrane_capacitance); - void paint(const region&, initial_ion_data); + void paint(const region&, init_int_concentration); + void paint(const region&, init_ext_concentration); + void paint(const region&, init_reversal_potential); // Synapses. lid_range place(const locset&, mechanism_desc); diff --git a/arbor/include/arbor/cable_cell_param.hpp b/arbor/include/arbor/cable_cell_param.hpp index 804ab131cfe3a2820ca28ddef47394f388955c3f..6f733ad2a41cec0738ef8a65e215dbdbb0a72e00 100644 --- a/arbor/include/arbor/cable_cell_param.hpp +++ b/arbor/include/arbor/cable_cell_param.hpp @@ -21,13 +21,14 @@ struct cable_cell_error: arbor_exception { // Ion inital concentration and reversal potential // parameters, as used in cable_cell_parameter_set, -// and set locally via painting initial_ion_data -// (see below). +// and set locally via painting init_int_concentration, +// init_ext_concentration and init_reversal_potential +// separately (see below). struct cable_cell_ion_data { - double init_int_concentration = NAN; - double init_ext_concentration = NAN; - double init_reversal_potential = NAN; + util::optional<double> init_int_concentration; + util::optional<double> init_ext_concentration; + util::optional<double> init_reversal_potential; }; // Current clamp description for stimulus specification. @@ -72,6 +73,21 @@ struct membrane_capacitance { double value = NAN; // [F/m²] }; +struct init_int_concentration { + std::string ion; + double value = NAN; +}; + +struct init_ext_concentration { + std::string ion; + double value = NAN; +}; + +struct init_reversal_potential { + std::string ion; + double value = NAN; +}; + // Mechanism description, viz. mechanism name and // (non-global) parameter settings. Used to assign // density and point mechanisms to segments and diff --git a/arbor/morph/region.cpp b/arbor/morph/region.cpp index 68ff60d380413f2ee087b28d0276fa6c4e2fb488..e78a59b4c8da8ce974a5e0323876d1927bf7c3b8 100644 --- a/arbor/morph/region.cpp +++ b/arbor/morph/region.cpp @@ -547,7 +547,7 @@ mextent thingify_(const named_& n, const mprovider& p) { } std::ostream& operator<<(std::ostream& o, const named_& x) { - return o << "(region \"" << x.name << "\")"; + return o << "(region " << x.name << ")"; } // Adds all cover points to a region. diff --git a/python/cells.cpp b/python/cells.cpp index 30d72bf3a78ef7178d12b8c148782be40655f07f..cc867f6a5a368f0ec7efe9ebebea88843b9503e3 100644 --- a/python/cells.cpp +++ b/python/cells.cpp @@ -165,7 +165,6 @@ global_props_shim::global_props_shim(): std::string to_string(const global_props_shim& G) { std::string s = "{arbor.cable_global_properties"; - auto nan_is_none = [](double x) {return x==x? std::to_string(x): "None";}; const auto& P = G.props; const auto& D = P.default_parameters; const auto& I = D.ion_data; @@ -182,9 +181,9 @@ std::string to_string(const global_props_shim& G) { "'"+D.reversal_potential_method.at(ion.first).name()+"'": "None"; s += util::pprintf("\n {name: '{}', valence: {}, int_con: {}, ext_con: {}, rev_pot: {}, rev_pot_method: {}}", ion.first, ion.second, - nan_is_none(props.init_int_concentration), - nan_is_none(props.init_ext_concentration), - nan_is_none(props.init_reversal_potential), + props.init_int_concentration, + props.init_ext_concentration, + props.init_reversal_potential, method); } } @@ -554,12 +553,20 @@ void register_cells(pybind11::module& m) { " rL: axial resistivity [Ω·cm].\n" " tempK: temperature [Kelvin].") + + // Paint ion species initial conditions on a region. .def("paint", - [](arb::cable_cell& c, const char* region, const arb::initial_ion_data& d) { - c.paint(region, d); + [](arb::cable_cell& c, const char* region, const char* name, + optional<double> int_con, optional<double> ext_con, optional<double> rev_pot) { + if (int_con) c.paint(region, arb::init_int_concentration{name, *int_con}); + if (ext_con) c.paint(region, arb::init_int_concentration{name, *ext_con}); + if (rev_pot) c.paint(region, arb::init_int_concentration{name, *rev_pot}); }, - "region"_a, "ion_data"_a, + "region"_a, "ion_name"_a, + pybind11::arg_v("int_con", pybind11::none(), "Intial internal concentration [mM]"), + pybind11::arg_v("ext_con", pybind11::none(), "Intial external concentration [mM]"), + pybind11::arg_v("rev_pot", pybind11::none(), "Intial reversal potential [mV]"), "Set ion species properties conditions on a region.") // Place synapses .def("place", diff --git a/test/unit/test_fvm_layout.cpp b/test/unit/test_fvm_layout.cpp index 0f55d42991f9c80279df36fb895fbce4e5cd4a4c..8403a30d4f829926c5d04812db44a930240bace4 100644 --- a/test/unit/test_fvm_layout.cpp +++ b/test/unit/test_fvm_layout.cpp @@ -769,8 +769,8 @@ TEST(fvm_layout, ion_weights) { gprop.catalogue = &testcat; gprop.default_parameters = neuron_parameter_defaults; - fvm_value_type cai = gprop.default_parameters.ion_data["ca"].init_int_concentration; - fvm_value_type cao = gprop.default_parameters.ion_data["ca"].init_ext_concentration; + fvm_value_type cai = gprop.default_parameters.ion_data["ca"].init_int_concentration.value(); + fvm_value_type cao = gprop.default_parameters.ion_data["ca"].init_ext_concentration.value(); for (auto& v: expected_init_iconc) { for (auto& iconc: v) {