diff --git a/arbor/fvm_layout.cpp b/arbor/fvm_layout.cpp index 27d59688c706ee05a87431e0fbda478228c90eb2..94a8e1071f2e8c4f3eacd5089975aea4a3d2fb35 100644 --- a/arbor/fvm_layout.cpp +++ b/arbor/fvm_layout.cpp @@ -16,8 +16,6 @@ #include <arbor/morph/mprovider.hpp> #include <arbor/morph/morphology.hpp> -#include <iostream> - #include "fvm_layout.hpp" #include "threading/threading.hpp" #include "util/maputil.hpp" @@ -315,11 +313,18 @@ ARB_ARBOR_API fvm_cv_discretization fvm_cv_discretize(const cable_cell& cell, co id_map.insert(k, 1.0/v.value); } } - // Fetch defaults, either global or per cell - auto gd = *(value_by_key(dflt.ion_data, ion).value().diffusivity - | value_by_key(global_dflt.ion_data, ion).value().diffusivity); - if (gd <= 0.0) { - throw cable_cell_error{util::pprintf("Illegal global diffusivity '{}' for ion '{}'.", gd, ion)}; + arb_value_type def = 0.0; + if (auto data = value_by_key(global_dflt.ion_data, ion); + data && data->diffusivity) { + def = data->diffusivity.value(); + } + if (auto data = value_by_key(dflt.ion_data, ion); + data && data->diffusivity) { + def = data->diffusivity.value(); + } + if (def <= 0.0) { + throw cable_cell_error{util::pprintf("Illegal global diffusivity '{}' for ion '{}'; possibly unset." + " Please define a positive global or cell default.", def, ion)}; } // Write inverse diffusivity / diffuse resistivity map @@ -328,7 +333,7 @@ ARB_ARBOR_API fvm_cv_discretization fvm_cv_discretize(const cable_cell& cell, co msize_t n_branch = D.geometry.n_branch(0); id.reserve(n_branch); for (msize_t i = 0; i<n_branch; ++i) { - auto pw = pw_over_cable(id_map, mcable{i, 0., 1.}, 1.0/gd); + auto pw = pw_over_cable(id_map, mcable{i, 0., 1.}, 1.0/def); id[0].push_back(pw); } // Prepare conductivity map diff --git a/example/diffusion/diffusion.cpp b/example/diffusion/diffusion.cpp index 6b0bb382ec0e4da5f5311542c0c10459c987e930..4a1454baba7b2c22a4b11528a0ac0b4fbf7d42fa 100644 --- a/example/diffusion/diffusion.cpp +++ b/example/diffusion/diffusion.cpp @@ -25,6 +25,7 @@ struct linear: public recipe { linear(double ext, double dx, double Xi, double beta): l{ext}, d{dx}, i{Xi}, b{beta} { gprop.default_parameters = neuron_parameter_defaults; gprop.default_parameters.discretization = cv_policy_max_extent{d}; + gprop.add_ion("bla", 1, 23, 42, 0, 0); } cell_size_type num_cells() const override { return 1; } @@ -37,12 +38,13 @@ struct linear: public recipe { // -----|----- segment_tree tree; tree.append(mnpos, { -l, 0, 0, 3}, {l, 0, 0, 3}, 1); + tree.append(0, { -l, 0, 0, 3}, {l, 0, 0, 3}, 2); // Setup decor decor; decor.set_default(init_int_concentration{"na", i}); - decor.set_default(ion_diffusivity{"na", b}); - decor.place("(location 0 0.5)"_ls, synapse("inject/x=na", {{"alpha", 200.0*l}}), "Zap"); - decor.paint("(all)"_reg, density("decay/x=na")); + decor.paint("(tag 1)"_reg, ion_diffusivity{"na", b}); + decor.place("(location 0 0.5)"_ls, synapse("inject/x=bla", {{"alpha", 200.0*l}}), "Zap"); + decor.paint("(all)"_reg, density("decay/x=bla")); return cable_cell({tree}, decor); } diff --git a/test/unit/test_diffusion.cpp b/test/unit/test_diffusion.cpp index de56e11aa0c8235117e94b8cbdb7877b46656a22..aa0e07b48334c78044b8bef7be3f2928667d90e7 100644 --- a/test/unit/test_diffusion.cpp +++ b/test/unit/test_diffusion.cpp @@ -288,3 +288,64 @@ TEST(diffusion, decay_by_event) { { 0.100000, 0.900000, 0.000000}}; EXPECT_TRUE(run(rec, exp)); } + +TEST(diffusion, setting_diffusivity) { + // Skeleton recipe + struct R: public recipe { + R() { + gprop.default_parameters = neuron_parameter_defaults; + // make a two region tree + tree.append(mnpos, { -1, 0, 0, 3}, {1, 0, 0, 3}, 1); + tree.append(0, { -1, 0, 0, 3}, {1, 0, 0, 3}, 2); + // Utilise diffusive ions + dec.place("(location 0 0.5)"_ls, synapse("inject/x=bla", {{"alpha", 200.0}}), "Zap"); + dec.paint("(all)"_reg, density("decay/x=bla")); + } + + cell_size_type num_cells() const override { return 1; } + cell_kind get_cell_kind(cell_gid_type) const override { return cell_kind::cable; } + std::any get_global_properties(cell_kind) const override { return gprop; } + util::unique_any get_cell_description(cell_gid_type) const override { return cable_cell({tree}, dec); } + + cable_cell_global_properties gprop; + segment_tree tree; + decor dec; + }; + + // BAD: Trying to use a diffusive ion, but b=0. + { + R r; + r.gprop.add_ion("bla", 1, 23, 42, 0, 0); + EXPECT_THROW(simulation(r).run(1, 1), illegal_diffusive_mechanism); + } + // BAD: Trying to use a partially diffusive ion + { + R r; + r.gprop.add_ion("bla", 1, 23, 42, 0, 0); + r.dec.paint("(tag 1)"_reg, ion_diffusivity{"bla", 13}); + EXPECT_THROW(simulation(r).run(1, 1), cable_cell_error); + } + // OK: Using the global default + { + R r; + r.gprop.add_ion("bla", 1, 23, 42, 0, 8); + r.dec.paint("(tag 1)"_reg, ion_diffusivity{"bla", 13}); + EXPECT_NO_THROW(simulation(r).run(1, 1)); + } + // OK: Using the cell default + { + R r; + r.gprop.add_ion("bla", 1, 23, 42, 0, 0); + r.dec.set_default(ion_diffusivity{"bla", 8}); + r.dec.paint("(tag 1)"_reg, ion_diffusivity{"bla", 13}); + EXPECT_NO_THROW(simulation(r).run(1, 1)); + } + // BAD: Using an unknown species + { + R r; + r.dec.set_default(ion_diffusivity{"bla", 8}); + r.dec.paint("(tag 1)"_reg, ion_diffusivity{"bla", 13}); + EXPECT_THROW(simulation(r).run(1, 1), cable_cell_error); + } + +}