diff --git a/doc/fig/performance_ratio_nonVec_vs_vec_compartmental.png b/doc/fig/performance_ratio_nonVec_vs_vec_compartmental.png new file mode 100644 index 000000000..7d53db935 Binary files /dev/null and b/doc/fig/performance_ratio_nonVec_vs_vec_compartmental.png differ diff --git a/doc/running/running_nest_compartmental.rst b/doc/running/running_nest_compartmental.rst index 712fe263e..da3a4583d 100644 --- a/doc/running/running_nest_compartmental.rst +++ b/doc/running/running_nest_compartmental.rst @@ -79,7 +79,7 @@ As an example for a HH-type channel: All of the currents within a compartment (marked by ``@mechanism::channel``) are added up within a compartment. -For a complete example, please see `cm_default.nestml `_ and its associated unit test, `compartmental_model_test.py `_. +For a complete example, please see `cm_default.nestml `_ and its associated unit test, `test__compartmental_model.py `_. Concentration description @@ -111,12 +111,14 @@ As an example a description of a calcium concentration model where we pretend th The only difference here is that the equation that is marked with the ``@mechanism::concentration`` descriptor is not an inline equation but an ODE. This is because in case of the ion-channel what we want to simulate is the current which relies on the evolution of some state variables like gating variables in case of the HH-models, and the compartment voltage. The concentration though can be more simply described by an evolving state directly. -For a complete example, please see `concmech.nestml `_ and its associated unit test, `compartmental_model_test.py `_. +For a complete example, please see `concmech.nestml `_ and its associated unit test, `test__concmech_model.py `_. Synapse description ------------------- -Here synapse models are based on convolutions over a buffer of incoming spikes. This means that the equation for the current-contribution must contain a convolve() call and a description of the kernel used for that convolution is needed. The descriptor for synapses is ``@mechanism::receptor``. +Here synapse models are based on convolutions over a buffer of incoming spikes. This means that the equation for the +current-contribution must contain a convolve() call and a description of the kernel used for that convolution is needed. +The descriptor for synapses is ``@mechanism::receptor``. .. code-block:: nestml @@ -133,13 +135,49 @@ Here synapse models are based on convolutions over a buffer of incoming spikes. input: <- spike -For a complete example, please see `concmech.nestml `_ and its associated unit test, `compartmental_model_test.py `_. +For a complete example, please see `concmech.nestml `_ and its associated unit test, `test__concmech_model.py `_. + +Continuous input description +---------------------------- + +The continuous inputs are defined by an inline with the descriptor @mechanism::continuous_input. This inline needs to +include one input of type continuous and may include any states, parameters and functions. + +.. code-block:: nestml + + model : + equations: + inline real = \ + > \ + @mechanism::continuous_input + + input: + real <- continuous + +For a complete example, please see `continuous_test.nestml `_ and its associated unit test, `test__continuous_input.py `_. Mechanism interdependence ------------------------- Above examples of explicit interdependence inbetween concentration and channel models where already described. Note that it is not necessary to describe the basic interaction inherent through the contribution to the overall current of the compartment. During a simulation step all currents of channels and synapses are added up and contribute to the change of the membrane potential (v_comp) in the next timestep. Thereby one must only express a dependence explicitly if the mechanism depends on the activity of a specific channel- or synapse-type amongst multiple in a given compartment or some concentration. +Technical Notes +--------------- + +We have put an emphasis on delivering good performance for neurons with high spatial complexity. We utilize vectorization, therefore, you should compile NEST with the OpenMP flag enabled. This, of course, can only be utilized if your hardware supports SIMD instructions. In that case, you can expect a performance improvement of about 3/4th of the theoretical maximum. + +Let's say you have an AVX2 SIMD instruction set available, which can fit 4 doubles (4*64-bit) into its vector register. In this case you can expect about a 3x performance improvement as long as your neuron has enough compartments. We vectorize the simulation steps of all instances of the same mechanism you have defined in your NESTML model, meaning that you will get a better complexity/performance ratio the more instances of the same mechanism are used. + +Here is a small benchmark example that shows the performance ratio (y-axis) as the number of compartments per neuron (x-axis) increases. + +.. figure:: https://raw.githubusercontent.com/nest/nestml/master/doc/fig/performance_ratio_nonVec_vs_vec_compartmental.png + :width: 326px + :height: 203px + :align: left + :target: # + +Be aware that we are using the -ffast-math flag when compiling the model by default. This can potentially lead to precision problems and inconsistencies across different systems. If you encounter unexpected results or want to be on the safe side, you can disable this by removing the flag from the CMakeLists.txt, which is part of the generated code. Note, however, that this may inhibit the compiler's ability to vectorize parts of the code in some cases. See also -------- diff --git a/extras/convert_cm_default_to_template.py b/extras/convert_cm_default_to_template.py index 92873d628..821955c0d 100644 --- a/extras/convert_cm_default_to_template.py +++ b/extras/convert_cm_default_to_template.py @@ -37,7 +37,7 @@ def get_replacement_patterns(): # file names 'cm_default' : '{{neuronSpecificFileNamesCmSyns[\"main\"]}}', 'cm_tree' : '{{neuronSpecificFileNamesCmSyns[\"tree\"]}}', - 'cm_compartmentcurrents': '{{neuronSpecificFileNamesCmSyns[\"compartmentcurrents\"]}}', + 'cm_neuroncurrents': '{{neuronSpecificFileNamesCmSyns[\"neuroncurrents\"]}}', # class names 'CompTree' : 'CompTree{{cm_unique_suffix}}', 'Compartment' : 'Compartment{{cm_unique_suffix}}', diff --git a/pynestml/cocos/co_co_cm_continuous_input_model.py b/pynestml/cocos/co_co_cm_continuous_input_model.py new file mode 100644 index 000000000..c3e6eb2fa --- /dev/null +++ b/pynestml/cocos/co_co_cm_continuous_input_model.py @@ -0,0 +1,36 @@ +# -*- coding: utf-8 -*- +# +# co_co_cm_continuous_input_model.py +# +# This file is part of NEST. +# +# Copyright (C) 2004 The NEST Initiative +# +# NEST is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 2 of the License, or +# (at your option) any later version. +# +# NEST is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with NEST. If not, see . + +from pynestml.cocos.co_co import CoCo +from pynestml.meta_model.ast_model import ASTModel +from pynestml.utils.continuous_input_processing import ContinuousInputProcessing + + +class CoCoCmContinuousInputModel(CoCo): + @classmethod + def check_co_co(cls, neuron: ASTModel): + """ + Checks if this compartmental condition applies to the handed over neuron. + If yes, it checks the presence of expected functions and declarations. + :param neuron: a single neuron instance. + :type neuron: ast_neuron + """ + return ContinuousInputProcessing.check_co_co(neuron) diff --git a/pynestml/cocos/co_cos_manager.py b/pynestml/cocos/co_cos_manager.py index 83900539b..01d008890 100644 --- a/pynestml/cocos/co_cos_manager.py +++ b/pynestml/cocos/co_cos_manager.py @@ -25,6 +25,7 @@ from pynestml.cocos.co_co_inline_expression_not_assigned_to import CoCoInlineExpressionNotAssignedTo from pynestml.cocos.co_co_input_port_not_assigned_to import CoCoInputPortNotAssignedTo from pynestml.cocos.co_co_cm_channel_model import CoCoCmChannelModel +from pynestml.cocos.co_co_cm_continuous_input_model import CoCoCmContinuousInputModel from pynestml.cocos.co_co_convolve_cond_correctly_built import CoCoConvolveCondCorrectlyBuilt from pynestml.cocos.co_co_convolve_has_correct_parameter import CoCoConvolveHasCorrectParameter from pynestml.cocos.co_co_input_port_not_assigned_to import CoCoInputPortNotAssignedTo @@ -148,6 +149,7 @@ def check_compartmental_model(cls, neuron: ASTModel) -> None: CoCoCmChannelModel.check_co_co(neuron) CoCoCmConcentrationModel.check_co_co(neuron) CoCoCmSynapseModel.check_co_co(neuron) + CoCoCmContinuousInputModel.check_co_co(neuron) @classmethod def check_inline_expressions_have_rhs(cls, model: ASTModel): diff --git a/pynestml/codegeneration/nest_compartmental_code_generator.py b/pynestml/codegeneration/nest_compartmental_code_generator.py index 2e0fc37b6..69a32bd9d 100644 --- a/pynestml/codegeneration/nest_compartmental_code_generator.py +++ b/pynestml/codegeneration/nest_compartmental_code_generator.py @@ -18,7 +18,7 @@ # # You should have received a copy of the GNU General Public License # along with NEST. If not, see . - +import shutil from typing import Any, Dict, List, Mapping, Optional import datetime @@ -53,10 +53,14 @@ from pynestml.meta_model.ast_variable import ASTVariable from pynestml.symbol_table.symbol_table import SymbolTable from pynestml.symbols.symbol import SymbolKind +from pynestml.utils.ast_vector_parameter_setter_and_printer import ASTVectorParameterSetterAndPrinter +from pynestml.utils.ast_vector_parameter_setter_and_printer_factory import ASTVectorParameterSetterAndPrinterFactory from pynestml.transformers.inline_expression_expansion_transformer import InlineExpressionExpansionTransformer from pynestml.utils.mechanism_processing import MechanismProcessing from pynestml.utils.channel_processing import ChannelProcessing from pynestml.utils.concentration_processing import ConcentrationProcessing +from pynestml.utils.continuous_input_processing import ContinuousInputProcessing +from pynestml.utils.con_in_info_enricher import ConInInfoEnricher from pynestml.utils.conc_info_enricher import ConcInfoEnricher from pynestml.utils.ast_utils import ASTUtils from pynestml.utils.chan_info_enricher import ChanInfoEnricher @@ -96,8 +100,8 @@ class NESTCompartmentalCodeGenerator(CodeGenerator): "path": "resources_nest_compartmental/cm_neuron", "model_templates": { "neuron": [ - "cm_compartmentcurrents_@NEURON_NAME@.cpp.jinja2", - "cm_compartmentcurrents_@NEURON_NAME@.h.jinja2", + "cm_neuroncurrents_@NEURON_NAME@.cpp.jinja2", + "cm_neuroncurrents_@NEURON_NAME@.h.jinja2", "@NEURON_NAME@.cpp.jinja2", "@NEURON_NAME@.h.jinja2", "cm_tree_@NEURON_NAME@.cpp.jinja2", @@ -151,7 +155,7 @@ def setup_printers(self): self._nest_printer = CppPrinter(expression_printer=self._printer) self._nest_variable_printer_no_origin = NESTVariablePrinter(None, with_origin=False, - with_vector_parameter=False, + with_vector_parameter=True, enforce_getter=False) self._printer_no_origin = CppExpressionPrinter( simple_expression_printer=CppSimpleExpressionPrinter(variable_printer=self._nest_variable_printer_no_origin, @@ -245,7 +249,7 @@ def _get_module_namespace(self, neurons: List[ASTModel]) -> Dict: neuron_name_to_filename = dict() for neuron in neurons: neuron_name_to_filename[neuron.get_name()] = { - "compartmentcurrents": self.get_cm_syns_compartmentcurrents_file_prefix(neuron), + "neuroncurrents": self.get_cm_syns_neuroncurrents_file_prefix(neuron), "main": self.get_cm_syns_main_file_prefix(neuron), "tree": self.get_cm_syns_tree_file_prefix(neuron) } @@ -261,6 +265,9 @@ def _get_module_namespace(self, neurons: List[ASTModel]) -> Dict: def get_cm_syns_compartmentcurrents_file_prefix(self, neuron): return "cm_compartmentcurrents_" + neuron.get_name() + def get_cm_syns_neuroncurrents_file_prefix(self, neuron): + return "cm_neuroncurrents_" + neuron.get_name() + def get_cm_syns_main_file_prefix(self, neuron): return neuron.get_name() @@ -525,7 +532,7 @@ def compute_name_of_generated_file(self, jinja_file_name, neuron): jinja_file_name).split(".")[0] file_name_calculators = { - "CompartmentCurrents": self.get_cm_syns_compartmentcurrents_file_prefix, + "NeuronCurrents": self.get_cm_syns_neuroncurrents_file_prefix, "Tree": self.get_cm_syns_tree_file_prefix, "Main": self.get_cm_syns_main_file_prefix, } @@ -588,6 +595,7 @@ def _get_neuron_model_namespace(self, neuron: ASTModel) -> Dict: namespace["nest_printer"] = self._nest_printer namespace["nestml_printer"] = NESTMLPrinter() namespace["type_symbol_printer"] = self._type_symbol_printer + namespace["vector_printer_factory"] = ASTVectorParameterSetterAndPrinterFactory(neuron, self._printer_no_origin) # NESTML syntax keywords namespace["PyNestMLLexer"] = {} @@ -700,14 +708,19 @@ def _get_neuron_model_namespace(self, neuron: ASTModel) -> Dict: namespace["conc_info"] = ConcentrationProcessing.get_mechs_info(neuron) namespace["conc_info"] = ConcInfoEnricher.enrich_with_additional_info(neuron, namespace["conc_info"]) + namespace["con_in_info"] = ContinuousInputProcessing.get_mechs_info(neuron) + namespace["con_in_info"] = ConInInfoEnricher.enrich_with_additional_info(neuron, namespace["con_in_info"]) + chan_info_string = MechanismProcessing.print_dictionary(namespace["chan_info"], 0) syns_info_string = MechanismProcessing.print_dictionary(namespace["syns_info"], 0) conc_info_string = MechanismProcessing.print_dictionary(namespace["conc_info"], 0) - code, message = Messages.get_mechs_dictionary_info(chan_info_string, syns_info_string, conc_info_string) + con_in_info_string = MechanismProcessing.print_dictionary(namespace["con_in_info"], 0) + + code, message = Messages.get_mechs_dictionary_info(chan_info_string, syns_info_string, conc_info_string, con_in_info_string) Logger.log_message(None, code, message, None, LoggingLevel.DEBUG) neuron_specific_filenames = { - "compartmentcurrents": self.get_cm_syns_compartmentcurrents_file_prefix(neuron), + "neuroncurrents": self.get_cm_syns_neuroncurrents_file_prefix(neuron), "main": self.get_cm_syns_main_file_prefix(neuron), "tree": self.get_cm_syns_tree_file_prefix(neuron)} diff --git a/pynestml/codegeneration/printers/nest_variable_printer.py b/pynestml/codegeneration/printers/nest_variable_printer.py index 48bbef59e..df2b6438b 100644 --- a/pynestml/codegeneration/printers/nest_variable_printer.py +++ b/pynestml/codegeneration/printers/nest_variable_printer.py @@ -20,6 +20,7 @@ # along with NEST. If not, see . from __future__ import annotations +from typing import Dict, Optional from pynestml.utils.ast_utils import ASTUtils @@ -49,6 +50,7 @@ def __init__(self, expression_printer: ExpressionPrinter, with_origin: bool = Tr self.with_vector_parameter = with_vector_parameter self.enforce_getter = enforce_getter self.variables_special_cases = variables_special_cases + self.cpp_variable_suffix = "" def print_variable(self, variable: ASTVariable) -> str: """ @@ -104,7 +106,9 @@ def print_variable(self, variable: ASTVariable) -> str: s = "" if not units_conversion_factor == 1: s += "(" + str(units_conversion_factor) + " * " - s += "B_." + self._print_buffer_value(variable) + if self.cpp_variable_suffix == "": + s += "B_." + s += self._print_buffer_value(variable) if not units_conversion_factor == 1: s += ")" return s @@ -112,17 +116,17 @@ def print_variable(self, variable: ASTVariable) -> str: if symbol.is_inline_expression: # there might not be a corresponding defined state variable; insist on calling the getter function if self.enforce_getter: - return "get_" + self._print(variable, symbol, with_origin=False) + vector_param + "()" + return "get_" + self._print(variable, symbol, with_origin=False) + vector_param + "()" + self.cpp_variable_suffix # modification to not enforce getter function: else: - return self._print(variable, symbol, with_origin=False) + return self._print(variable, symbol, with_origin=False) + self.cpp_variable_suffix assert not symbol.is_kernel(), "Cannot print kernel; kernel should have been converted during code generation" if symbol.is_state() or symbol.is_inline_expression: - return self._print(variable, symbol, with_origin=self.with_origin) + vector_param + return self._print(variable, symbol, with_origin=self.with_origin) + vector_param + self.cpp_variable_suffix - return self._print(variable, symbol, with_origin=self.with_origin) + vector_param + return self._print(variable, symbol, with_origin=self.with_origin) + vector_param + self.cpp_variable_suffix def _print_delay_variable(self, variable: ASTVariable) -> str: """ @@ -153,6 +157,9 @@ def _print_buffer_value(self, variable: ASTVariable) -> str: var_name += "_" + str(variable.get_vector_parameter()) return "spike_inputs_grid_sum_[" + var_name + " - MIN_SPIKE_RECEPTOR]" + if self.cpp_variable_suffix: + return variable_symbol.get_symbol_name() + self.cpp_variable_suffix + return variable_symbol.get_symbol_name() + '_grid_sum_' def _print(self, variable: ASTVariable, symbol, with_origin: bool = True) -> str: diff --git a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/@NEURON_NAME@.cpp.jinja2 b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/@NEURON_NAME@.cpp.jinja2 index 83e82f2d1..9fc7bf66f 100644 --- a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/@NEURON_NAME@.cpp.jinja2 +++ b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/@NEURON_NAME@.cpp.jinja2 @@ -76,7 +76,7 @@ nest::{{neuronSpecificFileNamesCmSyns["main"]}}::{{neuronSpecificFileNamesCmSyns * ---------------------------------------------------------------- */ void -{{neuronSpecificFileNamesCmSyns["main"]}}::get_status( DictionaryDatum& statusdict ) const +{{neuronSpecificFileNamesCmSyns["main"]}}::get_status( DictionaryDatum& statusdict ) { def< double >( statusdict, names::V_th, V_th_ ); ArchivingNode::get_status( statusdict ); @@ -99,7 +99,7 @@ void compartment_ad.push_back( dd ); // add receptor info - compartment->compartment_currents.add_receptor_info( receptor_ad, compartment->comp_index ); + c_tree_.neuron_currents.add_receptor_info( receptor_ad, compartment->comp_index ); } // add compartment info and receptor info to the status dictionary def< ArrayDatum >( statusdict, names::compartments, compartment_ad ); @@ -233,15 +233,14 @@ nest::{{neuronSpecificFileNamesCmSyns["main"]}}::add_receptor_( DictionaryDatum& syn_buffers_.push_back( buffer ); // add the receptor to the compartment - Compartment{{cm_unique_suffix}}* compartment = c_tree_.get_compartment( compartment_idx ); if ( dd->known( names::params ) ) { - compartment->compartment_currents.add_synapse( - receptor_type, syn_idx, getValue< DictionaryDatum >( dd, names::params ) ); + c_tree_.neuron_currents.add_mechanism( + receptor_type, compartment_idx, getValue< DictionaryDatum >( dd, names::params ), syn_idx ); } else { - compartment->compartment_currents.add_synapse( receptor_type, syn_idx ); + c_tree_.neuron_currents.add_mechanism( receptor_type, compartment_idx, syn_idx ); } } @@ -312,13 +311,13 @@ nest::{{neuronSpecificFileNamesCmSyns["main"]}}::update( Time const& origin, con for ( long lag = from; lag < to; ++lag ) { - const double v_0_prev = c_tree_.get_root()->v_comp; + const double v_0_prev = *(c_tree_.get_root()->v_comp); c_tree_.construct_matrix( lag ); c_tree_.solve_matrix(); // threshold crossing - if ( c_tree_.get_root()->v_comp >= V_th_ && v_0_prev < V_th_ ) + if ( *(c_tree_.get_root()->v_comp) >= V_th_ && v_0_prev < V_th_ ) { set_spiketime( Time::step( origin.get_steps() + lag + 1 ) ); @@ -353,8 +352,11 @@ nest::{{neuronSpecificFileNamesCmSyns["main"]}}::handle( CurrentEvent& e ) const double c = e.get_current(); const double w = e.get_weight(); - Compartment{{cm_unique_suffix}}* compartment = c_tree_.get_compartment_opt( e.get_rport() ); - compartment->currents.add_value( e.get_rel_delivery_steps( kernel().simulation_manager.get_slice_origin() ), w * c ); + assert( e.get_delay_steps() > 0 ); + assert( ( e.get_rport() >= 0 ) && ( ( size_t ) e.get_rport() < syn_buffers_.size() ) ); + + syn_buffers_[ e.get_rport() ].add_value( + e.get_rel_delivery_steps( kernel().simulation_manager.get_slice_origin() ), c*w ); } void diff --git a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/@NEURON_NAME@.h.jinja2 b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/@NEURON_NAME@.h.jinja2 index 642a724a7..a3c9f3665 100644 --- a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/@NEURON_NAME@.h.jinja2 +++ b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/@NEURON_NAME@.h.jinja2 @@ -29,7 +29,7 @@ #include "nest_types.h" #include "universal_data_logger.h" -#include "{{neuronSpecificFileNamesCmSyns["compartmentcurrents"]}}.h" +#include "{{neuronSpecificFileNamesCmSyns["neuroncurrents"]}}.h" #include "{{neuronSpecificFileNamesCmSyns["tree"]}}.h" namespace nest @@ -251,7 +251,7 @@ public: size_t handles_test_event( CurrentEvent&, size_t ); size_t handles_test_event( DataLoggingRequest&, size_t ); - void get_status( DictionaryDatum& ) const; + void get_status( DictionaryDatum& ); void set_status( const DictionaryDatum& ); private: diff --git a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_compartmentcurrents_@NEURON_NAME@.cpp.jinja2 b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_compartmentcurrents_@NEURON_NAME@.cpp.jinja2 deleted file mode 100644 index 98626fc96..000000000 --- a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_compartmentcurrents_@NEURON_NAME@.cpp.jinja2 +++ /dev/null @@ -1,417 +0,0 @@ -{#- -cm_compartmentcurrents_@NEURON_NAME@.cpp.jinja2 - -This file is part of NEST. - -Copyright (C) 2004 The NEST Initiative - -NEST is free software: you can redistribute it and/or modify -it under the terms of the GNU General Public License as published by -the Free Software Foundation, either version 2 of the License, or -(at your option) any later version. - -NEST is distributed in the hope that it will be useful, -but WITHOUT ANY WARRANTY; without even the implied warranty of -MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -GNU General Public License for more details. - -You should have received a copy of the GNU General Public License -along with NEST. If not, see . -#} -{%- if tracing %}/* generated by {{self._TemplateReference__context.name}} */ {% endif %} -{%- import 'directives_cpp/FunctionDeclaration.jinja2' as function_declaration with context %} -#include "{{neuronSpecificFileNamesCmSyns["compartmentcurrents"]}}.h" - -{%- set current_conductance_name_prefix = "g" %} -{%- set current_equilibrium_name_prefix = "e" %} -{% macro render_dynamic_channel_variable_name(variable_type, ion_channel_name) %} - {%- if variable_type == "gbar" %} - {{ current_conductance_name_prefix~"_"~ion_channel_name }} - {%- elif variable_type == "e" %} - {{ current_equilibrium_name_prefix~"_"~ion_channel_name }} - {%- endif %} -{%- endmacro %} - -{%- macro render_state_variable_name(pure_variable_name, ion_channel_name) %} - {{ pure_variable_name~"_"~ion_channel_name }} -{%- endmacro %} - -{% macro render_time_resolution_variable(synapse_info) %} -{# we assume here that there is only one such variable ! #} -{%- for analytic_helper_name, analytic_helper_info in synapse_info["analytic_helpers"].items() %} -{%- if analytic_helper_info["is_time_resolution"] %} - {{ analytic_helper_name }} -{%- endif %} -{%- endfor %} -{%- endmacro %} - -{% macro render_function_return_type(function) %} -{%- with %} - {%- set symbol = function.get_scope().resolve_to_symbol(function.get_name(), SymbolKind.FUNCTION) %} - {{ types_printer.print(symbol.get_return_type()) }} -{%- endwith %} -{%- endmacro %} - -{% macro render_inline_expression_type(inline_expression) %} -{%- with %} - {%- set symbol = inline_expression.get_scope().resolve_to_symbol(inline_expression.variable_name, SymbolKind.VARIABLE) %} - {{ types_printer.print(symbol.get_type_symbol()) }} -{%- endwith %} -{%- endmacro %} - -{% macro render_static_channel_variable_name(variable_type, ion_channel_name) %} - -{%- for ion_channel_nm, channel_info in chan_info.items() %} - {%- if ion_channel_nm == ion_channel_name %} - {%- for variable_tp, variable_info in channel_info["channel_parameters"].items() %} - {%- if variable_tp == variable_type %} - {%- set variable = variable_info["parameter_block_variable"] %} - {{ variable.name }} - {%- endif %} - {%- endfor %} - {%- endif %} -{%- endfor %} - -{%- endmacro %} - -{% macro render_channel_function(function, ion_channel_name) %} -{{ function_declaration.FunctionDeclaration(function, "nest::"~ion_channel_name~cm_unique_suffix~"::") }} -{ -{%- filter indent(2,True) %} -{%- with ast = function.get_block() %} -{%- include "directives_cpp/Block.jinja2" %} -{%- endwith %} -{%- endfilter %} -} -{%- endmacro %} - - -{%- for ion_channel_name, channel_info in chan_info.items() %} - -// {{ion_channel_name}} channel ////////////////////////////////////////////////////////////////// -nest::{{ion_channel_name}}{{cm_unique_suffix}}::{{ion_channel_name}}{{cm_unique_suffix}}() - -{%- for pure_variable_name, variable_info in channel_info["States"].items() %} -// state variable {{pure_variable_name -}} -{%- set variable = variable_info["ASTVariable"] %} -{%- set rhs_expression = variable_info["rhs_expression"] %} -{% if loop.first %}: {% else %}, {% endif %} -{{- variable.name}}({{ printer_no_origin.print(rhs_expression) -}}) -{%- endfor %} - -{% for variable_type, variable_info in channel_info["Parameters"].items() %} -// channel parameter {{variable_type -}} -{%- set variable = variable_info["ASTVariable"] %} -{%- set rhs_expression = variable_info["rhs_expression"] %} -,{{- variable.name }}({{ printer_no_origin.print(rhs_expression) -}}) -{%- endfor %} -{} - -nest::{{ion_channel_name}}{{cm_unique_suffix}}::{{ion_channel_name}}{{cm_unique_suffix}}(const DictionaryDatum& channel_params) - -{%- for pure_variable_name, variable_info in channel_info["States"].items() %} -// state variable {{pure_variable_name -}} -{%- set variable = variable_info["ASTVariable"] %} -{%- set rhs_expression = variable_info["rhs_expression"] %} -{% if loop.first %}: {% else %}, {% endif %} -{{- variable.name}}({{ printer_no_origin.print(rhs_expression) -}}) -{%- endfor %} - -{% for variable_type, variable_info in channel_info["Parameters"].items() %} -// channel parameter {{variable_type -}} -{%- set variable = variable_info["ASTVariable"] %} -{%- set rhs_expression = variable_info["rhs_expression"] %} -,{{- variable.name }}({{ printer_no_origin.print(rhs_expression) -}}) -{%- endfor %} -// update {{ion_channel_name}} channel parameters -{ - {%- for variable_type, variable_info in channel_info["Parameters"].items() %} - {%- set variable = variable_info["ASTVariable"] %} - {%- set dynamic_variable = render_dynamic_channel_variable_name(variable_type, ion_channel_name) %} //have to remove??????????? - // {{ion_channel_name}} channel parameter {{dynamic_variable }} - if( channel_params->known( "{{variable.name}}" ) ) - {{variable.name}} = getValue< double >( channel_params, "{{variable.name}}" ); - {%- endfor %} -} - -void -nest::{{ion_channel_name}}{{cm_unique_suffix}}::append_recordables(std::map< Name, double* >* recordables, - const long compartment_idx) -{ - // add state variables to recordables map - {%- for pure_variable_name, variable_info in channel_info["States"].items() %} - {%- set variable = variable_info["ASTVariable"] %} - ( *recordables )[ Name( "{{variable.name}}" + std::to_string(compartment_idx) )] = &{{variable.name}}; - {%- endfor %} - ( *recordables )[ Name( "i_tot_{{ion_channel_name}}" + std::to_string(compartment_idx) )] = &i_tot_{{ion_channel_name}}; -} - -std::pair< double, double > nest::{{ion_channel_name}}{{cm_unique_suffix}}::f_numstep(const double v_comp{% for ode in channel_info["Dependencies"]["concentrations"] %}, double {{ode.lhs.name}}{% endfor %} - {% for inline in channel_info["Dependencies"]["receptors"] %}, double {{inline.variable_name}}{% endfor %} - {% for inline in channel_info["Dependencies"]["channels"] %}, double {{inline.variable_name}}{% endfor %}) -{ - double g_val = 0., i_val = 0.; - - if({%- for key_zero_param in channel_info["RootInlineKeyZeros"] %} {{ key_zero_param }} > 1e-9 && {%- endfor %} true ){ - {% if channel_info["ODEs"].items()|length %} double {{ printer_no_origin.print(channel_info["time_resolution_var"]) }} = Time::get_resolution().get_ms(); {% endif %} - - {%- for ode_variable, ode_info in channel_info["ODEs"].items() %} - {%- for propagator, propagator_info in ode_info["transformed_solutions"][0]["propagators"].items() %} - double {{ propagator }} = {{ printer_no_origin.print(propagator_info["init_expression"]) }}; - {%- endfor %} - {%- for state, state_solution_info in ode_info["transformed_solutions"][0]["states"].items() %} - {{state}} = {{ printer_no_origin.print(state_solution_info["update_expression"]) }}; - {%- endfor %} - {%- endfor %} - - {%- set inline_expression = channel_info["root_expression"] %} - {%- set inline_expression_d = channel_info["inline_derivative"] %} - // compute the conductance of the {{ion_channel_name}} channel - this->i_tot_{{ion_channel_name}} = {{ printer_no_origin.print(inline_expression.get_expression()) }}; - // derivative - double d_i_tot_dv = {{ printer_no_origin.print(inline_expression_d) }}; - - g_val = - d_i_tot_dv / 2.; - i_val = this->i_tot_{{ion_channel_name}} - d_i_tot_dv * v_comp / 2.; - } - return std::make_pair(g_val, i_val); - -} - -{%- for function in channel_info["Functions"] %} -{{render_channel_function(function, ion_channel_name)}} -{%- endfor %} - -double nest::{{ion_channel_name}}{{cm_unique_suffix}}::get_current_{{ion_channel_name}}(){ - return this->i_tot_{{ion_channel_name}}; -} - -// {{ion_channel_name}} channel end /////////////////////////////////////////////////////////// -{% endfor %} -//////////////////////////////////////////////////////////////////////////////// - -{%- for synapse_name, synapse_info in syns_info.items() %} -// {{synapse_name}} synapse //////////////////////////////////////////////////////////////// -nest::{{synapse_name}}{{cm_unique_suffix}}::{{synapse_name}}{{cm_unique_suffix}}( const long syn_index ) - {%- for param_name, param_declaration in synapse_info["Parameters"].items() %} - {% if loop.first %}: {% else %}, {% endif %} - {{ param_name }}({{ printer_no_origin.print(param_declaration["rhs_expression"]) }}) - {%- endfor %} -{ - syn_idx = syn_index; -} - -nest::{{synapse_name}}{{cm_unique_suffix}}::{{synapse_name}}{{cm_unique_suffix}}( const long syn_index, const DictionaryDatum& receptor_params ) - {%- for param_name, param_declaration in synapse_info["Parameters"].items() %} - {% if loop.first %}: {% else %}, {% endif %} - {{ param_name }}({{ printer_no_origin.print(param_declaration["rhs_expression"]) }}) - {%- endfor %} -{ - syn_idx = syn_index; - - // update parameters - {%- for param_name, param_declaration in synapse_info["Parameters"].items() %} - if( receptor_params->known( "{{param_name}}" ) ) - {{param_name}} = getValue< double >( receptor_params, "{{param_name}}" ); - {%- endfor %} -} - -void -nest::{{synapse_name}}{{cm_unique_suffix}}::append_recordables(std::map< Name, double* >* recordables) -{ - {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} - ( *recordables )[ Name( "{{convolution_info["kernel"]["name"]}}" + std::to_string(syn_idx) )] = &{{convolution}}; - {%- endfor %} - ( *recordables )[ Name( "i_tot_{{synapse_name}}" + std::to_string(syn_idx) )] = &i_tot_{{synapse_name}}; -} - -{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} -void nest::{{synapse_name}}{{cm_unique_suffix}}::calibrate() -{%- else %} -void nest::{{synapse_name}}{{cm_unique_suffix}}::pre_run_hook() -{%- endif %} -{ - - const double {{render_time_resolution_variable(synapse_info)}} = Time::get_resolution().get_ms(); - - // set propagators to ode toolbox returned value - {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} - {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["propagators"].items()%} - {{state_variable_name}} = {{ printer_no_origin.print(state_variable_info["init_expression"]) }}; - {%- endfor %} - {%- endfor %} - - // initial values for user defined states - // warning: this shadows class variables - {%- for state_name, state_declaration in synapse_info["States"].items() %} - double {{state_name}} = {{ printer_no_origin.print(state_declaration["rhs_expression"]) }}; - {%- endfor %} - - // initial values for kernel state variables, set to zero - {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} - {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["kernel_states"].items()%} - {{state_variable_name}} = 0; - {%- endfor %} - {%- endfor %} - - // user declared internals in order they were declared - {%- for internal_name, internal_declaration in synapse_info["internals_used_declared"] %} - {{internal_name}} = {{ printer_no_origin.print(internal_declaration.get_expression()) }}; - {%- endfor %} - - {{synapse_info["buffer_name"]}}_->clear(); -} - -std::pair< double, double > nest::{{synapse_name}}{{cm_unique_suffix}}::f_numstep( const double v_comp, const long lag {% for ode in synapse_info["Dependencies"]["concentrations"] %}, double {{ode.lhs.name}}{% endfor %} - {% for inline in synapse_info["Dependencies"]["receptors"] %}, double {{inline.variable_name}}{% endfor %} - {% for inline in synapse_info["Dependencies"]["channels"] %}, double {{inline.variable_name}}{% endfor %}) -{ - // get spikes - double s_val = {{synapse_info["buffer_name"]}}_->get_value( lag ); // * g_norm_; - - //update ODE state variable - {% if synapse_info["ODEs"].items()|length %} double {{ printer_no_origin.print(synapse_info["time_resolution_var"]) }} = Time::get_resolution().get_ms(); {% endif %} - {%- for ode_variable, ode_info in synapse_info["ODEs"].items() %} - {%- for propagator, propagator_info in ode_info["transformed_solutions"][0]["propagators"].items() %} - double {{ propagator }} = {{ printer_no_origin.print(propagator_info["init_expression"]) }}; - {%- endfor %} - {%- for state, state_solution_info in ode_info["transformed_solutions"][0]["states"].items() %} - {{state}} = {{ printer_no_origin.print(state_solution_info["update_expression"]) }}; - {%- endfor %} - {%- endfor %} - - // update kernel state variable / compute synaptic conductance - {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} - {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["kernel_states"].items() %} - {{state_variable_name}} = {{ printer_no_origin.print(state_variable_info["update_expression"]) }}; - {{state_variable_name}} += s_val * {{ printer_no_origin.print(state_variable_info["init_expression"]) }}; - - {%- endfor %} - {%- endfor %} - - // total current - // this expression should be the transformed inline expression - this->i_tot_{{synapse_name}} = {{ printer_no_origin.print(synapse_info["root_expression"].get_expression()) }}; - - // derivative of that expression - // voltage derivative of total current - // compute derivative with respect to current with sympy - double d_i_tot_dv = {{ printer_no_origin.print(synapse_info["inline_expression_d"]) }}; - - // for numerical integration - double g_val = - d_i_tot_dv / 2.; - double i_val = this->i_tot_{{synapse_name}} - d_i_tot_dv * v_comp / 2.; - - return std::make_pair(g_val, i_val); - -} - -{%- for function in synapse_info["functions_used"] %} -{{ function_declaration.FunctionDeclaration(function, "nest::"~synapse_name~cm_unique_suffix~"::") }} -{ -{%- filter indent(2,True) %} -{%- with ast = function.get_block() %} -{%- include "directives_cpp/Block.jinja2" %} -{%- endwith %} -{%- endfilter %} -} -{%- endfor %} - - double nest::{{synapse_name}}{{cm_unique_suffix}}::get_current_{{synapse_name}}(){ - return this->i_tot_{{synapse_name}}; - } - -// {{synapse_name}} synapse end /////////////////////////////////////////////////////////// -{%- endfor %} - -//////////////////////////////// concentrations -{%- for concentration_name, concentration_info in conc_info.items() %} - -// {{ concentration_name }} concentration ////////////////////////////////////////////////////////////////// -nest::{{ concentration_name }}{{cm_unique_suffix}}::{{ concentration_name }}{{cm_unique_suffix}}(): -{%- set states_written = False %} -{%- for pure_variable_name, variable_info in concentration_info["States"].items() %} -// state variable {{pure_variable_name -}} -{%- set variable = variable_info["ASTVariable"] %} -{%- set rhs_expression = variable_info["rhs_expression"] %} -{% if loop.first %} {%- set states_written = True %} {% else %}, {% endif %} -{{- variable.name}}({{ printer_no_origin.print(rhs_expression) -}}) -{%- endfor %} - -{% for variable_type, variable_info in concentration_info["Parameters"].items() %} -// channel parameter {{variable_type -}} -{%- set variable = variable_info["ASTVariable"] %} -{%- set rhs_expression = variable_info["rhs_expression"] %} -{% if loop.first %} {% if states_written %}, {% endif %} {% else %}, {% endif %} -{{- variable.name }}({{ printer_no_origin.print(rhs_expression) -}}) -{%- endfor %} -{} - -nest::{{ concentration_name }}{{cm_unique_suffix}}::{{ concentration_name }}{{cm_unique_suffix}}(const DictionaryDatum& concentration_params): -{%- set states_written = False %} -{%- for pure_variable_name, variable_info in concentration_info["States"].items() %} -// state variable {{pure_variable_name -}} -{%- set variable = variable_info["ASTVariable"] %} -{%- set rhs_expression = variable_info["rhs_expression"] %} -{% if loop.first %} {%- set states_written = True %} {% else %}, {% endif %} -{{- variable.name}}({{ printer_no_origin.print(rhs_expression) -}}) -{%- endfor %} - -{% for variable_type, variable_info in concentration_info["Parameters"].items() %} -// channel parameter {{variable_type -}} -{%- set variable = variable_info["ASTVariable"] %} -{%- set rhs_expression = variable_info["rhs_expression"] %} -{% if loop.first %} {% if states_written %}, {% endif %} {% else %}, {% endif %} -{{- variable.name }}({{ printer_no_origin.print(rhs_expression) -}}) -{%- endfor %} -// update {{ concentration_name }} concentration parameters -{ - {%- for variable_type, variable_info in concentration_info["Parameters"].items() %} - {%- set variable = variable_info["ASTVariable"] %} - {%- set dynamic_variable = render_dynamic_channel_variable_name(variable_type, concentration_name) %} //have to remove??????????? - // {{ concentration_name }} concentration parameter {{dynamic_variable }} - if( concentration_params->known( "{{variable.name}}" ) ) - {{variable.name}} = getValue< double >( concentration_params, "{{variable.name}}" ); - {%- endfor %} -} - -void -nest::{{ concentration_name }}{{cm_unique_suffix}}::append_recordables(std::map< Name, double* >* recordables, - const long compartment_idx) -{ - // add state variables to recordables map - {%- for pure_variable_name, variable_info in concentration_info["States"].items() %} - {%- set variable = variable_info["ASTVariable"] %} - ( *recordables )[ Name( "{{variable.name}}" + std::to_string(compartment_idx) )] = &{{variable.name}}; - {%- endfor %} - ( *recordables )[ Name( "{{concentration_name}}" + std::to_string(compartment_idx) )] = &{{concentration_name}}; -} - -void nest::{{ concentration_name }}{{cm_unique_suffix}}::f_numstep(const double v_comp{% for ode in concentration_info["Dependencies"]["concentrations"] %}, double {{ode.lhs.name}}{% endfor %} - {% for inline in concentration_info["Dependencies"]["receptors"] %}, double {{inline.variable_name}}{% endfor %} - {% for inline in concentration_info["Dependencies"]["channels"] %}, double {{inline.variable_name}}{% endfor %}) -{ - if({%- for key_zero_param in concentration_info["RootInlineKeyZeros"] %} {{ key_zero_param }} > 1e-9 && {%- endfor %} true ){ - double {{ printer_no_origin.print(concentration_info["time_resolution_var"]) }} = Time::get_resolution().get_ms(); - - {%- for ode_variable, ode_info in concentration_info["ODEs"].items() %} - {%- for propagator, propagator_info in ode_info["transformed_solutions"][0]["propagators"].items() %} - double {{ propagator }} = {{ printer_no_origin.print(propagator_info["init_expression"]) }}; - {%- endfor %} - {%- for state, state_solution_info in ode_info["transformed_solutions"][0]["states"].items() %} - {{state}} = {{ printer_no_origin.print(state_solution_info["update_expression"]) }}; - {%- endfor %} - {%- endfor %} - } -} - -{%- for function in concentration_info["Functions"] %} -{{render_channel_function(function, concentration_name)}} -{%- endfor %} - -double nest::{{concentration_name}}{{cm_unique_suffix}}::get_concentration_{{concentration_name}}(){ - return this->{{concentration_name}}; -} - -// {{concentration_name}} concentration end /////////////////////////////////////////////////////////// -{% endfor %} diff --git a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_compartmentcurrents_@NEURON_NAME@.h.jinja2 b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_compartmentcurrents_@NEURON_NAME@.h.jinja2 deleted file mode 100644 index 508d3331b..000000000 --- a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_compartmentcurrents_@NEURON_NAME@.h.jinja2 +++ /dev/null @@ -1,470 +0,0 @@ -{#- -cm_compartmentcurrents_@NEURON_NAME@.h.jinja2 - -This file is part of NEST. - -Copyright (C) 2004 The NEST Initiative - -NEST is free software: you can redistribute it and/or modify -it under the terms of the GNU General Public License as published by -the Free Software Foundation, either version 2 of the License, or -(at your option) any later version. - -NEST is distributed in the hope that it will be useful, -but WITHOUT ANY WARRANTY; without even the implied warranty of -MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -GNU General Public License for more details. - -You should have received a copy of the GNU General Public License -along with NEST. If not, see . -#} -{%- if tracing %}/* generated by {{self._TemplateReference__context.name}} */ {% endif %} -{%- import 'directives_cpp/FunctionDeclaration.jinja2' as function_declaration with context %} -#ifndef SYNAPSES_NEAT_H_{{cm_unique_suffix | upper }} -#define SYNAPSES_NEAT_H_{{cm_unique_suffix | upper }} - -#include - -#include "ring_buffer.h" - -{% macro render_variable_type(variable) %} -{%- with %} - {%- set symbol = variable.get_scope().resolve_to_symbol(variable.name, SymbolKind.VARIABLE) %} - {{ types_printer.print(symbol.type_symbol) }} -{%- endwith %} -{%- endmacro %} - -namespace nest -{ - -{%- for ion_channel_name, channel_info in chan_info.items() %} - -class {{ion_channel_name}}{{cm_unique_suffix}}{ -private: - // states - {%- for pure_variable_name, variable_info in channel_info["States"].items() %} - {%- set variable = variable_info["ASTVariable"] %} - {%- set rhs_expression = variable_info["rhs_expression"] %} - {{ render_variable_type(variable) }} {{ variable.name }} = {{ printer_no_origin.print(rhs_expression) }}; - {%- endfor %} - - // parameters - {%- for pure_variable_name, variable_info in channel_info["Parameters"].items() %} - {%- set variable = variable_info["ASTVariable"] %} - {%- set rhs_expression = variable_info["rhs_expression"] %} - {{ render_variable_type(variable) }} {{ variable.name }} = {{ printer_no_origin.print(rhs_expression) }}; - {%- endfor %} - - // ion-channel root-inline value - double i_tot_{{ion_channel_name}} = 0; - -public: - // constructor, destructor - {{ion_channel_name}}{{cm_unique_suffix}}(); - {{ion_channel_name}}{{cm_unique_suffix}}(const DictionaryDatum& channel_params); - ~{{ion_channel_name}}{{cm_unique_suffix}}(){}; - - // initialization channel -{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} - void calibrate() { -{%- else %} - void pre_run_hook() { -{%- endif %} - // states - {%- for pure_variable_name, variable_info in channel_info["States"].items() %} - {%- set variable = variable_info["ASTVariable"] %} - {%- set rhs_expression = variable_info["rhs_expression"] %} - {{ variable.name }} = {{ printer_no_origin.print(rhs_expression) }}; - {%- endfor %} - }; - void append_recordables(std::map< Name, double* >* recordables, - const long compartment_idx); - - // numerical integration step - std::pair< double, double > f_numstep( const double v_comp{% for ode in channel_info["Dependencies"]["concentrations"] %}, double {{ode.lhs.name}}{% endfor %}{% if channel_info["Dependencies"]["receptors"]|length %} - {% endif %}{% for inline in channel_info["Dependencies"]["receptors"] %}, double {{inline.variable_name}}{% endfor %}{% if channel_info["Dependencies"]["channels"]|length %} - {% endif %}{% for inline in channel_info["Dependencies"]["channels"] %}, double {{inline.variable_name}}{% endfor %}); - - // function declarations - -{%- for function in channel_info["Functions"] %} - {{ function_declaration.FunctionDeclaration(function) }}; -{%- endfor %} - - // root_inline getter - double get_current_{{ion_channel_name}}(); - -}; -{% endfor %} - - -////////////////////////////////////////////////// synapses - -{% macro render_time_resolution_variable(synapse_info) %} -{# we assume here that there is only one such variable ! #} -{%- for analytic_helper_name, analytic_helper_info in synapse_info["analytic_helpers"].items() %} -{%- if analytic_helper_info["is_time_resolution"] %} - {{ analytic_helper_name }} -{%- endif %} -{%- endfor %} -{%- endmacro %} - -{%- for synapse_name, synapse_info in syns_info.items() %} - -class {{synapse_name}}{{cm_unique_suffix}}{ -private: - // global synapse index - long syn_idx = 0; - - // propagators, initialized via pre_run_hook() or calibrate() - {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} - {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["propagators"].items()%} - double {{state_variable_name}}; - {%- endfor %} - {%- endfor %} - - // kernel state variables, initialized via pre_run_hook() or calibrate() - {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} - {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["kernel_states"].items()%} - double {{state_variable_name}}; - {%- endfor %} - {%- endfor %} - - // user defined parameters, initialized via pre_run_hook() or calibrate() - {%- for param_name, param_declaration in synapse_info["Parameters"].items() %} - double {{param_name}}; - {%- endfor %} - - // states - {%- for pure_variable_name, variable_info in synapse_info["States"].items() %} - {%- set variable = variable_info["ASTVariable"] %} - {%- set rhs_expression = variable_info["rhs_expression"] %} - {{ render_variable_type(variable) }} {{ variable.name }} = {{ printer_no_origin.print(rhs_expression) }}; - {%- endfor %} - double i_tot_{{synapse_name}} = 0; - - // user declared internals in order they were declared, initialized via pre_run_hook() or calibrate() - {%- for internal_name, internal_declaration in synapse_info["internals_used_declared"] %} - double {{internal_name}}; - {%- endfor %} - - - - // spike buffer - RingBuffer* {{synapse_info["buffer_name"]}}_; - -public: - // constructor, destructor - {{synapse_name}}{{cm_unique_suffix}}( const long syn_index); - {{synapse_name}}{{cm_unique_suffix}}( const long syn_index, const DictionaryDatum& receptor_params); - ~{{synapse_name}}{{cm_unique_suffix}}(){}; - - long - get_syn_idx() - { - return syn_idx; - }; - - // numerical integration step - std::pair< double, double > f_numstep( const double v_comp, const long lag {% for ode in synapse_info["Dependencies"]["concentrations"] %}, double {{ode.lhs.name}}{% endfor %}{% if synapse_info["Dependencies"]["receptors"]|length %} - {% endif %}{% for inline in synapse_info["Dependencies"]["receptors"] %}, double {{inline.variable_name}}{% endfor %}{% if synapse_info["Dependencies"]["channels"]|length %} - {% endif %}{% for inline in synapse_info["Dependencies"]["channels"] %}, double {{inline.variable_name}}{% endfor %}); - - // calibration -{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} - void calibrate(); -{%- else %} - void pre_run_hook(); -{%- endif %} - void append_recordables(std::map< Name, double* >* recordables); - void set_buffer_ptr( std::vector< RingBuffer >& syn_buffers ) - { - {{synapse_info["buffer_name"]}}_ = &syn_buffers[ syn_idx ]; - }; - - // function declarations - {%- for function in synapse_info["Functions"] %} - {{ function_declaration.FunctionDeclaration(function, "") -}}; - - {% endfor %} - - // root_inline getter - double get_current_{{synapse_name}}(); -}; - - -{% endfor %} - -///////////////////////////////////////////// concentrations - -{%- for concentration_name, concentration_info in conc_info.items() %} - -class {{ concentration_name }}{{cm_unique_suffix}}{ -private: - // parameters - {%- for pure_variable_name, variable_info in concentration_info["Parameters"].items() %} - {%- set variable = variable_info["ASTVariable"] %} - {%- set rhs_expression = variable_info["rhs_expression"] %} - {{ render_variable_type(variable) }} {{ variable.name }} = {{ printer_no_origin.print(rhs_expression) }}; - {%- endfor %} - - // states - {%- for pure_variable_name, variable_info in concentration_info["States"].items() %} - {%- set variable = variable_info["ASTVariable"] %} - {%- set rhs_expression = variable_info["rhs_expression"] %} - {{ render_variable_type(variable) }} {{ variable.name }} = {{ printer_no_origin.print(rhs_expression) }}; - {%- endfor %} - - // concentration value (root-ode state) - double {{concentration_name}} = 0; - -public: - // constructor, destructor - {{ concentration_name }}{{cm_unique_suffix}}(); - {{ concentration_name }}{{cm_unique_suffix}}(const DictionaryDatum& concentration_params); - ~{{ concentration_name }}{{cm_unique_suffix}}(){}; - - // initialization channel -{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} - void calibrate() { -{%- else %} - void pre_run_hook() { -{%- endif %} - // states - {%- for pure_variable_name, variable_info in concentration_info["States"].items() %} - {%- set variable = variable_info["ASTVariable"] %} - {%- set rhs_expression = variable_info["rhs_expression"] %} - {{ variable.name }} = {{ printer_no_origin.print(rhs_expression) }}; - {%- endfor %} - }; - void append_recordables(std::map< Name, double* >* recordables, - const long compartment_idx); - - // numerical integration step - void f_numstep( const double v_comp{% for ode in concentration_info["Dependencies"]["concentrations"] %}, double {{ode.lhs.name}}{% endfor %}{% if concentration_info["Dependencies"]["receptors"]|length %} - {% endif %}{% for inline in concentration_info["Dependencies"]["receptors"] %}, double {{inline.variable_name}}{% endfor %}{% if concentration_info["Dependencies"]["channels"]|length %} - {% endif %}{% for inline in concentration_info["Dependencies"]["channels"] %}, double {{inline.variable_name}}{% endfor %}); - - // function declarations -{%- for function in concentration_info["Functions"] %} - {{ function_declaration.FunctionDeclaration(function) }}; -{%- endfor %} - - // root_ode getter - double get_concentration_{{concentration_name}}(); - -}; -{% endfor %} - -///////////////////////////////////////////// currents - -{%- set channel_suffix = "_chan_" %} -{%- set concentration_suffix = "_conc_" %} - -class CompartmentCurrents{{cm_unique_suffix}} { -private: - // ion channels -{% with %} - {%- for ion_channel_name, channel_info in chan_info.items() %} - {{ion_channel_name}}{{cm_unique_suffix}} {{ion_channel_name}}{{channel_suffix}}; - {% endfor %} -{% endwith %} - - // synapses - {%- for synapse_name, synapse_info in syns_info.items() %} - std::vector < {{synapse_name}}{{cm_unique_suffix}} > {{synapse_name}}_syns_; - {% endfor %} - - //concentrations -{% with %} - {%- for concentration_name, concentration_info in conc_info.items() %} - {{concentration_name}}{{cm_unique_suffix}} {{concentration_name}}{{concentration_suffix}}; - {% endfor %} -{% endwith %} - -public: - CompartmentCurrents{{cm_unique_suffix}}(){}; - explicit CompartmentCurrents{{cm_unique_suffix}}(const DictionaryDatum& compartment_params) - { - {%- for ion_channel_name, channel_info in chan_info.items() %} - {{ion_channel_name}}{{channel_suffix}} = {{ion_channel_name}}{{cm_unique_suffix}}( compartment_params ); - {% endfor %} - - {%- for concentration_name, concentration_info in conc_info.items() %} - {{ concentration_name }}{{concentration_suffix}} = {{ concentration_name }}{{cm_unique_suffix}}( compartment_params ); - {% endfor %} - }; - ~CompartmentCurrents{{cm_unique_suffix}}(){}; - -{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} - void calibrate() { -{%- else %} - void pre_run_hook() { -{%- endif %} - // initialization of ion channels - {%- for ion_channel_name, channel_info in chan_info.items() %} -{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} - {{ion_channel_name}}{{channel_suffix}}.calibrate(); -{%- else %} - {{ion_channel_name}}{{channel_suffix}}.pre_run_hook(); -{%- endif %} - {% endfor %} - - // initialization of concentrations - {%- for concentration_name, concentration_info in conc_info.items() %} -{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} - {{ concentration_name }}{{concentration_suffix}}.calibrate(); -{%- else %} - {{ concentration_name }}{{concentration_suffix}}.pre_run_hook(); -{%- endif %} - {% endfor %} - - // initialization of synapses - {%- for synapse_name, synapse_info in syns_info.items() %} - // initialization of {{synapse_name}} synapses - for( auto syn_it = {{synapse_name}}_syns_.begin(); - syn_it != {{synapse_name}}_syns_.end(); - ++syn_it ) - { -{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} - syn_it->calibrate(); -{%- else %} - syn_it->pre_run_hook(); -{%- endif %} - } - {% endfor %} - }; - - void add_synapse( const std::string& type, const long syn_idx ) - { - {%- for synapse_name, synapse_info in syns_info.items() %} - {% if not loop.first %}else{% endif %} if ( type == "{{synapse_name}}" ) - { - {{synapse_name}}_syns_.push_back( {{synapse_name}}{{cm_unique_suffix}}( syn_idx ) ); - } - {% endfor %} - else - { - assert( false ); - } - }; - void add_synapse( const std::string& type, const long syn_idx, const DictionaryDatum& receptor_params ) - { - {%- for synapse_name, synapse_info in syns_info.items() %} - {% if not loop.first %}else{% endif %} if ( type == "{{synapse_name}}" ) - { - {{synapse_name}}_syns_.push_back( {{synapse_name}}{{cm_unique_suffix}}( syn_idx, receptor_params ) ); - } - {% endfor %} - else - { - assert( false ); - } - }; - - void - add_receptor_info( ArrayDatum& ad, const long compartment_index ) - { - {%- for synapse_name, synapse_info in syns_info.items() %} - for( auto syn_it = {{synapse_name}}_syns_.begin(); syn_it != {{synapse_name}}_syns_.end(); syn_it++) - { - DictionaryDatum dd = DictionaryDatum( new Dictionary ); - def< long >( dd, names::receptor_idx, syn_it->get_syn_idx() ); - def< long >( dd, names::comp_idx, compartment_index ); - def< std::string >( dd, names::receptor_type, "{{synapse_name}}" ); - ad.push_back( dd ); - } - {% endfor %} - }; - - void - set_syn_buffers( std::vector< RingBuffer >& syn_buffers ) - { - // spike buffers for synapses - {%- for synapse_name, synapse_info in syns_info.items() %} - for( auto syn_it = {{synapse_name}}_syns_.begin(); syn_it != {{synapse_name}}_syns_.end(); syn_it++) - syn_it->set_buffer_ptr( syn_buffers ); - {% endfor %} - }; - - std::map< Name, double* > - get_recordables( const long compartment_idx ) - { - std::map< Name, double* > recordables; - - // append ion channel state variables to recordables - {%- for ion_channel_name, channel_info in chan_info.items() %} - {{ion_channel_name}}{{channel_suffix}}.append_recordables( &recordables, compartment_idx ); - {% endfor %} - - // append concentration state variables to recordables - {%- for concentration_name, concentration_info in conc_info.items() %} - {{concentration_name}}{{concentration_suffix}}.append_recordables( &recordables, compartment_idx ); - {% endfor %} - - // append synapse state variables to recordables - {%- for synapse_name, synapse_info in syns_info.items() %} - for( auto syn_it = {{synapse_name}}_syns_.begin(); syn_it != {{synapse_name}}_syns_.end(); syn_it++) - syn_it->append_recordables( &recordables ); - {% endfor %} - - return recordables; - }; - - std::pair< double, double > - f_numstep( const double v_comp, const long lag ) - { - std::pair< double, double > gi(0., 0.); - double g_val = 0.; - double i_val = 0.; -{%- for synapse_name, synapse_info in syns_info.items() %} - double {{synapse_name}}{{channel_suffix}}current_sum = 0; - for( auto syn_it = {{synapse_name}}_syns_.begin(); - syn_it != {{synapse_name}}_syns_.end(); - ++syn_it ) - { - {{synapse_name}}{{channel_suffix}}current_sum += syn_it->get_current_{{synapse_name}}(); - } -{% endfor %} - - {%- for concentration_name, concentration_info in conc_info.items() %} - // computation of {{ concentration_name }} concentration - {{ concentration_name }}{{concentration_suffix}}.f_numstep( v_comp{% for ode in concentration_info["Dependencies"]["concentrations"] %}, {{ode.lhs.name}}{{concentration_suffix}}.get_concentration_{{ode.lhs.name}}(){% endfor %}{% if concentration_info["Dependencies"]["receptors"]|length %} - {% endif %}{% for inline in concentration_info["Dependencies"]["receptors"] %}, {{inline.variable_name}}{{channel_suffix}}_current_sum{% endfor %}{% if concentration_info["Dependencies"]["channels"]|length %} - {% endif %}{% for inline in concentration_info["Dependencies"]["channels"] %}, {{inline.variable_name}}{{channel_suffix}}.get_current_{{inline.variable_name}}(){% endfor %}); - - {% endfor %} - - {%- for ion_channel_name, channel_info in chan_info.items() %} - // contribution of {{ion_channel_name}} channel - gi = {{ion_channel_name}}{{channel_suffix}}.f_numstep( v_comp{% for ode in channel_info["Dependencies"]["concentrations"] %}, {{ode.lhs.name}}{{concentration_suffix}}.get_concentration_{{ode.lhs.name}}(){% endfor %}{% if channel_info["Dependencies"]["receptors"]|length %} - {% endif %}{% for inline in channel_info["Dependencies"]["receptors"] %}, {{inline.variable_name}}{{channel_suffix}}_current_sum{% endfor %}{% if channel_info["Dependencies"]["channels"]|length %} - {% endif %}{% for inline in channel_info["Dependencies"]["channels"] %}, {{inline.variable_name}}{{channel_suffix}}.get_current_{{inline.variable_name}}(){% endfor %}); - - g_val += gi.first; - i_val += gi.second; - - {% endfor %} - - {%- for synapse_name, synapse_info in syns_info.items() %} - // contribution of {{synapse_name}} synapses - for( auto syn_it = {{synapse_name}}_syns_.begin(); - syn_it != {{synapse_name}}_syns_.end(); - ++syn_it ) - { - gi = syn_it->f_numstep( v_comp, lag {% for ode in synapse_info["Dependencies"]["concentrations"] %}, {{ode.lhs.name}}{{concentration_suffix}}.get_concentration_{{ode.lhs.name}}(){% endfor %}{% if synapse_info["Dependencies"]["receptors"]|length %} - {% endif %}{% for inline in synapse_info["Dependencies"]["receptors"] %}, {{inline.variable_name}}{{channel_suffix}}_current_sum{% endfor %}{% if synapse_info["Dependencies"]["channels"]|length %} - {% endif %}{% for inline in synapse_info["Dependencies"]["channels"] %}, {{inline.variable_name}}{{channel_suffix}}.get_current{{inline.variable_name}}(){% endfor %}); - - g_val += gi.first; - i_val += gi.second; - } - {% endfor %} - - return std::make_pair(g_val, i_val); - }; -}; - -} // namespace - -#endif /* #ifndef SYNAPSES_NEAT_H_{{cm_unique_suffix | upper }} */ diff --git a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_neuroncurrents_@NEURON_NAME@.cpp.jinja2 b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_neuroncurrents_@NEURON_NAME@.cpp.jinja2 new file mode 100644 index 000000000..fdf83b974 --- /dev/null +++ b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_neuroncurrents_@NEURON_NAME@.cpp.jinja2 @@ -0,0 +1,1030 @@ +{#- +cm_neuroncurrents_@NEURON_NAME@.cpp.jinja2 + +This file is part of NEST. + +Copyright (C) 2004 The NEST Initiative + +NEST is free software: you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation, either version 2 of the License, or +(at your option) any later version. + +NEST is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License +along with NEST. If not, see . +#} +{%- if tracing %}/* generated by {{self._TemplateReference__context.name}} */ {% endif -%} +{%- import 'directives_cpp/FunctionDeclaration.jinja2' as function_declaration with context %} +#include "{{neuronSpecificFileNamesCmSyns["neuroncurrents"]}}.h" + +#define NEAR_ZERO 1e-9 + +{%- set current_conductance_name_prefix = "g" %} +{%- set current_equilibrium_name_prefix = "e" %} +{% macro render_dynamic_channel_variable_name(variable_type, ion_channel_name) -%} + {%- if variable_type == "gbar" -%} + {{ current_conductance_name_prefix~"_"~ion_channel_name }} + {%- elif variable_type == "e" -%} + {{ current_equilibrium_name_prefix~"_"~ion_channel_name }} + {%- endif -%} +{%- endmacro -%} + +{%- macro render_state_variable_name(pure_variable_name, ion_channel_name) -%} + {{ pure_variable_name~"_"~ion_channel_name }} +{%- endmacro -%} + +{% macro render_time_resolution_variable(synapse_info) -%} +{# we assume here that there is only one such variable ! #} +{%- with %} +{%- for analytic_helper_name, analytic_helper_info in synapse_info["analytic_helpers"].items() -%} +{%- if analytic_helper_info["is_time_resolution"] -%} + {{ analytic_helper_name }} +{%- endif -%} +{%- endfor -%} +{% endwith %} +{%- endmacro %} + +{% macro render_function_return_type(function) -%} +{%- with -%} + {%- set symbol = function.get_scope().resolve_to_symbol(function.get_name(), SymbolKind.FUNCTION) -%} + {{ types_printer.print(symbol.get_return_type()) }} +{%- endwith -%} +{%- endmacro -%} + +{% macro render_inline_expression_type(inline_expression) -%} +{%- with -%} + {%- set symbol = inline_expression.get_scope().resolve_to_symbol(inline_expression.variable_name, SymbolKind.VARIABLE) -%} + {{ types_printer.print(symbol.get_type_symbol()) }} +{%- endwith -%} +{%- endmacro -%} + +{% macro render_static_channel_variable_name(variable_type, ion_channel_name) -%} + +{%- with %} +{%- for ion_channel_nm, channel_info in chan_info.items() -%} + {%- if ion_channel_nm == ion_channel_name -%} + {%- for variable_tp, variable_info in channel_info["channel_parameters"].items() -%} + {%- if variable_tp == variable_type -%} + {%- set variable = variable_info["parameter_block_variable"] -%} + {{ variable.name }} + {%- endif -%} + {%- endfor -%} + {%- endif -%} +{%- endfor -%} +{% endwith %} + +{%- endmacro %} + +{% macro render_channel_function(function, ion_channel_name) -%} +{%- with %} +inline {{ function_declaration.FunctionDeclaration(function, "nest::"~ion_channel_name~cm_unique_suffix~"::") }} +{ +{%- filter indent(2,True) %} +{%- with ast = function.get_block() %} +{%- include "directives_cpp/Block.jinja2" %} +{%- endwith %} +{%- endfilter %} +} +{% endwith %} +{%- endmacro %} + +{% macro render_vectorized_channel_function(function, ion_channel_name) -%} +{%- with %} +{{ vectorized_function_declaration.FunctionDeclaration(function, "nest::"~ion_channel_name~cm_unique_suffix~"::") }} +{ +{%- filter indent(2,True) %} +{%- with ast = function.get_block() %} +{%- include "directives_cpp/VectorizedBlock.jinja2" %} +{%- endwith %} +{%- endfilter %} +} +{% endwith %} +{%- endmacro %} + +{%- macro vectorized_function_call(ast_function, ion_channel_name) -%} +{%- with function_symbol = ast_function.get_scope().resolve_to_symbol(ast_function.get_name(), SymbolKind.FUNCTION) -%} +{%- if function_symbol is none -%} +{{ raise('Cannot resolve the method ' + ast_function.get_name()) }} +{%- endif %} +{{ "std::vector< " + type_symbol_printer.print(function_symbol.get_return_type()) + " >" | replace('.', '::') }} {{ ast_function.get_name() }}_v(neuron_{{ ion_channel_name }}_channel_count); +{{ ast_function.get_name() }}( +{%- for param in ast_function.get_parameters() %} +{%- with typeSym = param.get_data_type().get_type_symbol() -%} +{%- filter indent(1, True) -%} +{{ param.get_name() }} +{%- if not loop.last -%} +, +{%- endif -%} +{%- endfilter -%} +{%- endwith -%} +{%- endfor -%} +, {{ ast_function.get_name() }}_v ); +{%- endwith -%} +{%- endmacro -%} + + +{%- with %} +{%- for ion_channel_name, channel_info in chan_info.items() %} + +// {{ion_channel_name}} channel ////////////////////////////////////////////////////////////////// +void nest::{{ion_channel_name}}{{cm_unique_suffix}}::new_channel(std::size_t comp_ass) +{ + //Check whether the channel will contribute at all based on initial key-parameters. If not then don't add the channel. + bool channel_contributing = true; + {%- for key_zero_param in channel_info["RootInlineKeyZeros"] %} + {% for variable_type, variable_info in channel_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {% if key_zero_param == variable.name %} + if(std::abs({{ printer_no_origin.print(rhs_expression) }}) <= NEAR_ZERO){ + channel_contributing = false; + } + {% endif %} + {%- endfor %} + {%- endfor %} + + if(channel_contributing){ + neuron_{{ ion_channel_name }}_channel_count++; + i_tot_{{ion_channel_name}}.push_back(0); + compartment_association.push_back(comp_ass); + + {%- for pure_variable_name, variable_info in channel_info["States"].items() %} + // state variable {{pure_variable_name}} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+ion_channel_name+"_channel_count").print(rhs_expression) -}}); + {%- endfor %} + + {% for variable_type, variable_info in channel_info["Parameters"].items() %} + // channel parameter {{variable_type }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+ion_channel_name+"_channel_count").print(rhs_expression) -}}); + {%- endfor %} + + {% for variable_type, variable_info in channel_info["Internals"].items() %} + // channel parameter {{variable_type }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+ion_channel_name+"_channel_count").print(rhs_expression) -}}); + {%- endfor %} + } +} + +void nest::{{ion_channel_name}}{{cm_unique_suffix}}::new_channel(std::size_t comp_ass, const DictionaryDatum& channel_params) +/* update {{ion_channel_name}} channel parameters and states */ +{ + //Check whether the channel will contribute at all based on initial key-parameters. If not then don't add the channel. + bool channel_contributing = true; + {%- for key_zero_param in channel_info["RootInlineKeyZeros"] %} + if( channel_params->known( "{{key_zero_param}}" ) ){ + if(std::abs(getValue< double >( channel_params, "{{key_zero_param}}" )) <= NEAR_ZERO){ + channel_contributing = false; + } + }else{ + {% for variable_type, variable_info in channel_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {% if key_zero_param == variable.name %} + if(std::abs({{ printer_no_origin.print(rhs_expression) }}) <= NEAR_ZERO){ + channel_contributing = false; + } + {% endif %} + {%- endfor %} + } + {%- endfor %} + + if(channel_contributing){ + neuron_{{ ion_channel_name }}_channel_count++; + compartment_association.push_back(comp_ass); + i_tot_{{ion_channel_name}}.push_back(0); + + {%- for pure_variable_name, variable_info in channel_info["States"].items() %} + // state variable {{pure_variable_name }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+ion_channel_name+"_channel_count").print(rhs_expression) -}}); + {%- endfor %} + + {%- with %} + {%- for variable_type, variable_info in channel_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set dynamic_variable = render_dynamic_channel_variable_name(variable_type, ion_channel_name) %} + // {{ion_channel_name}} channel parameter {{dynamic_variable }} + if( channel_params->known( "{{variable.name}}" ) ) + {{variable.name}}[neuron_{{ ion_channel_name }}_channel_count-1] = getValue< double >( channel_params, "{{variable.name}}" ); + {%- endfor %} + {% endwith %} + + {%- with %} + {%- for variable_type, variable_info in channel_info["ODEs"].items() %} + {%- set variable_name = variable_type %} + {%- set dynamic_variable = render_dynamic_channel_variable_name(variable_type, ion_channel_name) %} + // {{ion_channel_name}} channel ODE state {{dynamic_variable }} + if( channel_params->known( "{{variable_name}}" ) ) + {{variable_name}}[neuron_{{ ion_channel_name }}_channel_count-1] = getValue< double >( channel_params, "{{variable_name}}" ); + {%- endfor %} + {% endwith %} + + {% for variable_type, variable_info in channel_info["Parameters"].items() %} + // channel parameter {{variable_type }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+ion_channel_name+"_channel_count").print(rhs_expression) -}}); + {%- endfor %} + + {%- with %} + {%- for variable_type, variable_info in channel_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set dynamic_variable = render_dynamic_channel_variable_name(variable_type, ion_channel_name) %} + // {{ion_channel_name}} channel parameter {{dynamic_variable }} + if( channel_params->known( "{{variable.name}}" ) ) + {{variable.name}}[neuron_{{ ion_channel_name }}_channel_count-1] = getValue< double >( channel_params, "{{variable.name}}" ); + {%- endfor %} + {% endwith %} + + {%- for pure_variable_name, variable_info in channel_info["Internals"].items() %} + // state variable {{pure_variable_name }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+ion_channel_name+"_channel_count").print(rhs_expression) -}}); + {%- endfor %} + } +} + +void +nest::{{ion_channel_name}}{{cm_unique_suffix}}::append_recordables(std::map< Name, double* >* recordables, + const long compartment_idx) +{ + // add state variables to recordables map + bool found_rec = false; + {%- with %} + {%- for pure_variable_name, variable_info in channel_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + found_rec = false; + for(size_t chan_id = 0; chan_id < neuron_{{ ion_channel_name }}_channel_count; chan_id++){ + if(compartment_association[chan_id] == compartment_idx){ + ( *recordables )[ Name( std::string("{{variable.name}}") + std::to_string(compartment_idx))] = &{{variable.name}}[chan_id]; + found_rec = true; + } + } + if(!found_rec) ( *recordables )[ Name( std::string("{{variable.name}}") + std::to_string(compartment_idx))] = &zero_recordable; + {%- endfor %} + {% endwith %} + found_rec = false; + for(size_t chan_id = 0; chan_id < neuron_{{ ion_channel_name }}_channel_count; chan_id++){ + if(compartment_association[chan_id] == compartment_idx){ + ( *recordables )[ Name( std::string("i_tot_{{ion_channel_name}}") + std::to_string(compartment_idx))] = &i_tot_{{ion_channel_name}}[chan_id]; + found_rec = true; + } + } + if(!found_rec) ( *recordables )[ Name( std::string("i_tot_{{ion_channel_name}}") + std::to_string(compartment_idx))] = &zero_recordable; +} + +std::pair< std::vector< double >, std::vector< double > > nest::{{ion_channel_name}}{{cm_unique_suffix}}::f_numstep(std::vector< double > v_comp{% for ode in channel_info["Dependencies"]["concentrations"] %}, std::vector< double > {{ode.lhs.name}}{% endfor %}{% if channel_info["Dependencies"]["receptors"]|length %} + {% endif %}{% for inline in channel_info["Dependencies"]["receptors"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if channel_info["Dependencies"]["channels"]|length %} + {% endif %}{% for inline in channel_info["Dependencies"]["channels"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if channel_info["Dependencies"]["continuous"]|length %} + {% endif %}{% for inline in channel_info["Dependencies"]["continuous"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}) +{ + std::vector< double > g_val(neuron_{{ ion_channel_name }}_channel_count, 0.); + std::vector< double > i_val(neuron_{{ ion_channel_name }}_channel_count, 0.); + + std::vector< double > d_i_tot_dv(neuron_{{ ion_channel_name }}_channel_count, 0.); + + {% if channel_info["ODEs"].items()|length %} std::vector< double > {{ printer_no_origin.print(channel_info["time_resolution_var"]) }}(neuron_{{ ion_channel_name }}_channel_count, Time::get_resolution().get_ms()); {% endif %} + + {%- for ode_variable, ode_info in channel_info["ODEs"].items() %} + {%- for propagator, propagator_info in ode_info["transformed_solutions"][0]["propagators"].items() %} + std::vector< double > {{ propagator }}(neuron_{{ ion_channel_name }}_channel_count, 0); + {%- endfor %} + {%- endfor %} + #pragma omp simd + for(std::size_t i = 0; i < neuron_{{ ion_channel_name }}_channel_count; i++){ + {%- for ode_variable, ode_info in channel_info["ODEs"].items() %} + {%- for propagator, propagator_info in ode_info["transformed_solutions"][0]["propagators"].items() %} + {{ propagator }}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(propagator_info["init_expression"]) }}; + {%- endfor %} + {%- for state, state_solution_info in ode_info["transformed_solutions"][0]["states"].items() %} + {{state}}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(state_solution_info["update_expression"]) }}; + {%- endfor %} + {%- endfor %} + + {%- set inline_expression = channel_info["root_expression"] %} + {%- set inline_expression_d = channel_info["inline_derivative"] %} + + // compute the conductance of the {{ion_channel_name}} channel + this->i_tot_{{ion_channel_name}}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(inline_expression.get_expression()) }}; + + // derivative + d_i_tot_dv[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(inline_expression_d) }}; + g_val[i] = - d_i_tot_dv[i]; + i_val[i] = this->i_tot_{{ion_channel_name}}[i] - d_i_tot_dv[i] * v_comp[i]; + } + return std::make_pair(g_val, i_val); + +} + +{%- for function in channel_info["Functions"] %} +{{render_channel_function(function, ion_channel_name)}} +{%- endfor %} +void nest::{{ion_channel_name}}{{cm_unique_suffix}}::get_currents_per_compartment(std::vector< double >& compartment_to_current){ + for(std::size_t comp_id = 0; comp_id < compartment_to_current.size(); comp_id++){ + compartment_to_current[comp_id] = 0; + } + for(std::size_t chan_id = 0; chan_id < neuron_{{ ion_channel_name }}_channel_count; chan_id++){ + compartment_to_current[this->compartment_association[chan_id]] += this->i_tot_{{ion_channel_name}}[chan_id]; + } +} + +std::vector< double > nest::{{ion_channel_name}}{{cm_unique_suffix}}::distribute_shared_vector(std::vector< double > shared_vector){ + std::vector< double > distributed_vector(this->neuron_{{ ion_channel_name }}_channel_count, 0.0); + for(std::size_t chan_id = 0; chan_id < this->neuron_{{ ion_channel_name }}_channel_count; chan_id++){ + distributed_vector[chan_id] = shared_vector[compartment_association[chan_id]]; + } + return distributed_vector; +} + +// {{ion_channel_name}} channel end /////////////////////////////////////////////////////////// +{% endfor -%} +{% endwith %} +//////////////////////////////////////////////////////////////////////////////// + +//////////////////////////////// concentrations +{%- with %} +{%- for concentration_name, concentration_info in conc_info.items() %} + +// {{ concentration_name }} concentration ///////////////////////////////////////////////////// + +void nest::{{concentration_name}}{{cm_unique_suffix}}::new_concentration(std::size_t comp_ass) +{ + //Check whether the concentration will contribute at all based on initial key-parameters. If not then don't add the concentration. + bool concentration_contributing = true; + {%- for key_zero_param in concentration_info["RootInlineKeyZeros"] %} + {% for variable_type, variable_info in concentration_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {% if key_zero_param == variable.name %} + if(std::abs({{ printer_no_origin.print(rhs_expression) }}) <= NEAR_ZERO){ + concentration_contributing = false; + } + {% endif %} + {%- endfor %} + {%- endfor %} + + if(concentration_contributing){ + neuron_{{ concentration_name }}_concentration_count++; + {{concentration_name}}.push_back(0); + compartment_association.push_back(comp_ass); + + {%- for pure_variable_name, variable_info in concentration_info["States"].items() %} + // state variable {{pure_variable_name }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+concentration_name+"_concentration_count").print(rhs_expression) -}}); + {%- endfor %} + + {% for variable_type, variable_info in concentration_info["Parameters"].items() %} + // channel parameter {{variable_type }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+concentration_name+"_concentration_count").print(rhs_expression) -}}); + {%- endfor %} + + {% for variable_type, variable_info in concentration_info["Internals"].items() %} + // channel parameter {{variable_type }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+concentration_name+"_concentration_count").print(rhs_expression) -}}); + {%- endfor %} + } +} + +void nest::{{concentration_name}}{{cm_unique_suffix}}::new_concentration(std::size_t comp_ass, const DictionaryDatum& concentration_params) +{ + //Check whether the concentration will contribute at all based on initial key-parameters. If not then don't add the concentration. + bool concentration_contributing = true; + {%- for key_zero_param in concentration_info["RootInlineKeyZeros"] %} + if( concentration_params->known( "{{key_zero_param}}" ) ){ + if(std::abs(getValue< double >( concentration_params, "{{key_zero_param}}" )) <= NEAR_ZERO){ + concentration_contributing = false; + } + }else{ + {% for variable_type, variable_info in concentration_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {% if key_zero_param == variable.name %} + if(std::abs({{ printer_no_origin.print(rhs_expression) }}) <= NEAR_ZERO){ + concentration_contributing = false; + } + {% endif %} + {%- endfor %} + } + {%- endfor %} + + if(concentration_contributing){ + neuron_{{ concentration_name }}_concentration_count++; + {{concentration_name}}.push_back(0); + compartment_association.push_back(comp_ass); + + {%- for pure_variable_name, variable_info in concentration_info["States"].items() %} + // state variable {{pure_variable_name }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+concentration_name+"_concentration_count").print(rhs_expression) -}}); + {%- endfor %} + + {%- with %} + {%- for variable_type, variable_info in concentration_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set dynamic_variable = render_dynamic_channel_variable_name(variable_type, ion_channel_name) %} + // {{ion_channel_name}} channel parameter {{dynamic_variable }} + if( concentration_params->known( "{{variable.name}}" ) ) + {{variable.name}}[neuron_{{ concentration_name }}_concentration_count-1] = getValue< double >( concentration_params, "{{variable.name}}" ); + {%- endfor %} + {% endwith %} + + {%- with %} + {%- for variable_type, variable_info in concentration_info["ODEs"].items() %} + {%- set variable_name = variable_type %} + {%- set dynamic_variable = render_dynamic_channel_variable_name(variable_type, ion_channel_name) %} + // {{concentration_name}} concentration ODE state {{dynamic_variable }} + if( concentration_params->known( "{{variable_name}}" ) ) + {{variable_name}}[neuron_{{ concentration_name }}_concentration_count-1] = getValue< double >( concentration_params, "{{variable_name}}" ); + {%- endfor %} + {% endwith %} + + {% for variable_type, variable_info in concentration_info["Parameters"].items() %} + // channel parameter {{variable_type }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+concentration_name+"_concentration_count").print(rhs_expression) -}}); + {%- endfor %} + + {%- with %} + {%- for variable_type, variable_info in concentration_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set dynamic_variable = render_dynamic_channel_variable_name(variable_type, concentration_name) %} + // {{ concentration_name }} concentration parameter {{dynamic_variable }} + if( concentration_params->known( "{{variable.name}}" ) ) + {{variable.name}}[neuron_{{ concentration_name }}_concentration_count-1] = getValue< double >( concentration_params, "{{variable.name}}" ); + {%- endfor %} + {% endwith %} + + {%- for pure_variable_name, variable_info in concentration_info["Internals"].items() %} + // state variable {{pure_variable_name }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+concentration_name+"_concentration_count").print(rhs_expression) -}}); + {%- endfor %} + } +} + +void +nest::{{ concentration_name }}{{cm_unique_suffix}}::append_recordables(std::map< Name, double* >* recordables, + const long compartment_idx) +{ + // add state variables to recordables map + bool found_rec = false; + {%- with %} + {%- for pure_variable_name, variable_info in concentration_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + found_rec = false; + for(size_t conc_id = 0; conc_id < neuron_{{ concentration_name }}_concentration_count; conc_id++){ + if(compartment_association[conc_id] == compartment_idx){ + ( *recordables )[ Name( std::string("{{variable.name}}") + std::to_string(compartment_idx))] = &{{variable.name}}[conc_id]; + found_rec = true; + } + } + if(!found_rec) ( *recordables )[ Name( std::string("{{variable.name}}") + std::to_string(compartment_idx))] = &zero_recordable; + {%- endfor %} + {% endwith %} + found_rec = false; + for(size_t conc_id = 0; conc_id < neuron_{{ concentration_name }}_concentration_count; conc_id++){ + if(compartment_association[conc_id] == compartment_idx){ + ( *recordables )[ Name( std::string("{{concentration_name}}") + std::to_string(compartment_idx))] = &{{concentration_name}}[conc_id]; + found_rec = true; + } + } + if(!found_rec) ( *recordables )[ Name( std::string("{{concentration_name}}") + std::to_string(compartment_idx))] = &zero_recordable; +} + +void nest::{{ concentration_name }}{{cm_unique_suffix}}::f_numstep(std::vector< double > v_comp{% for ode in concentration_info["Dependencies"]["concentrations"] %}, std::vector< double > {{ode.lhs.name}}{% endfor %}{% if concentration_info["Dependencies"]["receptors"]|length %} + {% endif %}{% for inline in concentration_info["Dependencies"]["receptors"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if concentration_info["Dependencies"]["channels"]|length %} + {% endif %}{% for inline in concentration_info["Dependencies"]["channels"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if concentration_info["Dependencies"]["continuous"]|length %} + {% endif %}{% for inline in concentration_info["Dependencies"]["continuous"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}) +{ + std::vector< double > {{ printer_no_origin.print(concentration_info["time_resolution_var"]) }}(neuron_{{ concentration_name }}_concentration_count, Time::get_resolution().get_ms()); + + {%- for ode_variable, ode_info in concentration_info["ODEs"].items() %} + {%- for propagator, propagator_info in ode_info["transformed_solutions"][0]["propagators"].items() %} + std::vector< double > {{ propagator }}(neuron_{{ concentration_name }}_concentration_count, 0); + {%- endfor %} + {%- endfor %} + + #pragma omp simd + for(std::size_t i = 0; i < neuron_{{ concentration_name }}_concentration_count; i++){ + {%- for ode_variable, ode_info in concentration_info["ODEs"].items() %} + {%- for propagator, propagator_info in ode_info["transformed_solutions"][0]["propagators"].items() %} + {{ propagator }}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(propagator_info["init_expression"]) }}; + {%- endfor %} + {%- for state, state_solution_info in ode_info["transformed_solutions"][0]["states"].items() %} + {{state}}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(state_solution_info["update_expression"]) }}; + {%- endfor %} + {%- endfor %} + } +} + +{%- for function in concentration_info["Functions"] %} +{{render_channel_function(function, concentration_name)}} +{%- endfor %} + +void nest::{{concentration_name}}{{cm_unique_suffix}}::get_concentrations_per_compartment(std::vector< double >& compartment_to_concentration){ + for(std::size_t comp_id = 0; comp_id < compartment_to_concentration.size(); comp_id++){ + compartment_to_concentration[comp_id] = 0; + } + for(std::size_t conc_id = 0; conc_id < neuron_{{ concentration_name }}_concentration_count; conc_id++){ + compartment_to_concentration[this->compartment_association[conc_id]] += this->{{concentration_name}}[conc_id]; + } +} + +std::vector< double > nest::{{concentration_name}}{{cm_unique_suffix}}::distribute_shared_vector(std::vector< double > shared_vector){ + std::vector< double > distributed_vector(this->neuron_{{ concentration_name }}_concentration_count, 0.0); + for(std::size_t conc_id = 0; conc_id < this->neuron_{{ concentration_name }}_concentration_count; conc_id++){ + distributed_vector[conc_id] = shared_vector[compartment_association[conc_id]]; + } + return distributed_vector; +} + +// {{concentration_name}} concentration end /////////////////////////////////////////////////////////// +{% endfor -%} +{% endwith %} + + +////////////////////////////////////// synapses + +{%- for synapse_name, synapse_info in syns_info.items() %} +// {{synapse_name}} synapse //////////////////////////////////////////////////////////////// + +void nest::{{synapse_name}}{{cm_unique_suffix}}::new_synapse(std::size_t comp_ass, const long syn_index) +{ + neuron_{{ synapse_name }}_synapse_count++; + i_tot_{{synapse_name}}.push_back(0); + compartment_association.push_back(comp_ass); + syn_idx.push_back(syn_index); + + {%- for pure_variable_name, variable_info in synapse_info["States"].items() %} + // state variable {{pure_variable_name }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+synapse_name+"_synapse_count").print(rhs_expression) -}}); + {%- endfor %} + + {% for variable_type, variable_info in synapse_info["Parameters"].items() %} + // synapse parameter {{variable_type }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+synapse_name+"_synapse_count").print(rhs_expression) -}}); + {%- endfor %} + + // set propagators to ode toolbox returned value + {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} + {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["propagators"].items()%} + {{state_variable_name}}.push_back(0); + {%- endfor %} + {%- endfor %} + + // initial values for kernel state variables, set to zero + {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} + {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["kernel_states"].items()%} + {{state_variable_name}}.push_back(0); + {%- endfor %} + {%- endfor %} + + // user declared internals in order they were declared + {%- for internal_name, internal_declaration in synapse_info["internals_used_declared"] %} + {{internal_name}}.push_back(0); + {%- endfor %} +} + +void nest::{{synapse_name}}{{cm_unique_suffix}}::new_synapse(std::size_t comp_ass, const long syn_index, const DictionaryDatum& synapse_params) +/* update {{synapse}} synapse parameters and states */ +{ + neuron_{{ synapse_name }}_synapse_count++; + compartment_association.push_back(comp_ass); + i_tot_{{synapse_name}}.push_back(0); + syn_idx.push_back(syn_index); + + {%- for pure_variable_name, variable_info in synapse_info["States"].items() %} + // state variable {{pure_variable_name }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+synapse_name+"_synapse_count").print(rhs_expression) -}}); + {%- endfor %} + {%- with %} + {%- for variable_type, variable_info in synapse_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + if( synapse_params->known( "{{variable.name}}" ) ) + {{variable.name}}[neuron_{{ synapse_name }}_synapse_count-1] = getValue< double >( synapse_params, "{{variable.name}}" ); + {%- endfor %} + {% endwith %} + + {%- with %} + {%- for variable_type, variable_info in synapse_info["ODEs"].items() %} + {%- set variable_name = variable_type %} + {%- set dynamic_variable = render_dynamic_channel_variable_name(variable_type, ion_channel_name) %} + // {{concentration_name}} concentration ODE state {{dynamic_variable }} + if( synapse_params->known( "{{variable_name}}" ) ) + {{variable_name}}[neuron_{{ synapse_name }}_synapse_count-1] = getValue< double >( synapse_params, "{{variable_name}}" ); + {%- endfor %} + {% endwith %} + + {% for variable_type, variable_info in synapse_info["Parameters"].items() %} + // synapse parameter {{variable_type }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+synapse_name+"_synapse_count").print(rhs_expression) -}}); + {%- endfor %} + + {%- with %} + {%- for variable_type, variable_info in synapse_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + if( synapse_params->known( "{{variable.name}}" ) ) + {{variable.name}}[neuron_{{ synapse_name }}_synapse_count-1] = getValue< double >( synapse_params, "{{variable.name}}" ); + {%- endfor %} + {% endwith %} + + // set propagators to ode toolbox returned value + {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} + {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["propagators"].items()%} + {{state_variable_name}}.push_back(0); + {%- endfor %} + {%- endfor %} + + // initial values for kernel state variables, set to zero + {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} + {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["kernel_states"].items()%} + {{state_variable_name}}.push_back(0); + {%- endfor %} + {%- endfor %} + + // user declared internals in order they were declared + {%- for internal_name, internal_declaration in synapse_info["internals_used_declared"] %} + {{internal_name}}.push_back(0); + {%- endfor %} +} + +void +nest::{{synapse_name}}{{cm_unique_suffix}}::append_recordables(std::map< Name, double* >* recordables, const long compartment_idx) +{ + {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} + for(size_t syns_id = 0; syns_id < neuron_{{ synapse_name }}_synapse_count; syns_id++){ + if(compartment_association[syns_id] == compartment_idx){ + ( *recordables )[ Name( "{{convolution_info["kernel"]["name"]}}" + std::to_string(syns_id) )] = &{{convolution}}[syns_id]; + } + } + {%- endfor %} + for(size_t syns_id = 0; syns_id < neuron_{{ synapse_name }}_synapse_count; syns_id++){ + if(compartment_association[syns_id] == compartment_idx){ + ( *recordables )[ Name( "i_tot_{{synapse_name}}" + std::to_string(syns_id) )] = &i_tot_{{synapse_name}}[syns_id]; + } + } +} + +{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} +void nest::{{synapse_name}}{{cm_unique_suffix}}::calibrate() +{%- else %} +void nest::{{synapse_name}}{{cm_unique_suffix}}::pre_run_hook() +{%- endif %} +{ + + std::vector< double > {{ printer_no_origin.print(synapse_info["analytic_helpers"]["__h"]["ASTVariable"]) }}(neuron_{{ synapse_name }}_synapse_count, Time::get_resolution().get_ms()); + + {%- for state_name, state_declaration in synapse_info["States"].items() %} + std::vector< double > {{state_name}} = (neuron_{{ synapse_name }}_synapse_count, {{ printer_no_origin.print(state_declaration["rhs_expression"])}}); + {%- endfor %} + + for(std::size_t i = 0; i < neuron_{{ synapse_name }}_synapse_count; i++){ + // set propagators to ode toolbox returned value + {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} + {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["propagators"].items()%} + {{state_variable_name}}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(state_variable_info["init_expression"]) }}; + {%- endfor %} + {%- endfor %} + + // initial values for kernel state variables, set to zero + {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} + {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["kernel_states"].items()%} + {{state_variable_name}}[i] = 0; + {%- endfor %} + {%- endfor %} + + // user declared internals in order they were declared + {%- for internal_name, internal_declaration in synapse_info["internals_used_declared"] %} + {{internal_name}}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(internal_declaration.get_expression()) }}; + {%- endfor %} + + {{synapse_info["buffer_name"]}}_[i]->clear(); + } +} + +std::pair< std::vector< double >, std::vector< double > > nest::{{synapse_name}}{{cm_unique_suffix}}::f_numstep( std::vector< double > v_comp, const long lag {% for ode in synapse_info["Dependencies"]["concentrations"] %}, std::vector< double > {{ode.lhs.name}}{% endfor %}{% if synapse_info["Dependencies"]["receptors"]|length %} + {% endif %}{% for inline in synapse_info["Dependencies"]["receptors"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if synapse_info["Dependencies"]["channels"]|length %} + {% endif %}{% for inline in synapse_info["Dependencies"]["channels"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if synapse_info["Dependencies"]["continuous"]|length %} + {% endif %}{% for inline in synapse_info["Dependencies"]["continuous"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}) +{ + std::vector< double > g_val(neuron_{{ synapse_name }}_synapse_count, 0.); + std::vector< double > i_val(neuron_{{ synapse_name }}_synapse_count, 0.); + std::vector< double > d_i_tot_dv(neuron_{{ synapse_name }}_synapse_count, 0.); + + {%- for ode_variable, ode_info in synapse_info["ODEs"].items() %} + {%- for propagator, propagator_info in ode_info["transformed_solutions"][0]["propagators"].items() %} + std::vector< double > {{ propagator }}(neuron_{{ synapse_name }}_synapse_count, 0); + {%- endfor %} + {%- endfor %} + + {% if synapse_info["ODEs"].items()|length %} std::vector< double > {{ printer_no_origin.print(synapse_info["analytic_helpers"]["__h"]["ASTVariable"]) }}(neuron_{{ synapse_name }}_synapse_count, Time::get_resolution().get_ms()); {% endif %} + + std::vector < double > s_val(neuron_{{ synapse_name }}_synapse_count, 0); + + for(std::size_t i = 0; i < neuron_{{ synapse_name }}_synapse_count; i++){ + // get spikes + s_val[i] = {{synapse_info["buffer_name"]}}_[i]->get_value( lag ); // * g_norm_; + } + + //update ODE state variable + #pragma omp simd + for(std::size_t i = 0; i < neuron_{{ synapse_name }}_synapse_count; i++){ + {%- for ode_variable, ode_info in synapse_info["ODEs"].items() %} + {%- for propagator, propagator_info in ode_info["transformed_solutions"][0]["propagators"].items() %} + {{ propagator }}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(propagator_info["init_expression"]) }}; + {%- endfor %} + {%- for state, state_solution_info in ode_info["transformed_solutions"][0]["states"].items() %} + {{state}}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(state_solution_info["update_expression"]) }}; + {%- endfor %} + {%- endfor %} + + // update kernel state variable / compute synaptic conductance + {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} + {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["kernel_states"].items() %} + {{state_variable_name}}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(state_variable_info["update_expression"]) }}; + {{state_variable_name}}[i] += s_val[i] * {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(state_variable_info["init_expression"]) }}; + {%- endfor %} + {%- endfor %} + + // total current + // this expression should be the transformed inline expression + + this->i_tot_{{synapse_name}}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(synapse_info["root_expression"].get_expression()) }}; + + // derivative of that expression + // voltage derivative of total current + // compute derivative with respect to current with sympy + d_i_tot_dv[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(synapse_info["inline_expression_d"]) }}; + + // for numerical integration + g_val[i] = - d_i_tot_dv[i]; + i_val[i] = this->i_tot_{{synapse_name}}[i] - d_i_tot_dv[i] * v_comp[i]; + } + + return std::make_pair(g_val, i_val); + +} + +{%- for function in synapse_info["functions_used"] %} +inline {{ function_declaration.FunctionDeclaration(function, "nest::"~synapse_name~cm_unique_suffix~"::") }} +{ +{%- filter indent(2,True) %} +{%- with ast = function.get_block() %} +{%- include "directives/Block.jinja2" %} +{%- endwith %} +{%- endfilter %} +} +{%- endfor %} + +void nest::{{synapse_name}}{{cm_unique_suffix}}::get_currents_per_compartment(std::vector< double >& compartment_to_current){ + for(std::size_t comp_id = 0; comp_id < compartment_to_current.size(); comp_id++){ + compartment_to_current[comp_id] = 0; + } + for(std::size_t syn_id = 0; syn_id < neuron_{{ synapse_name }}_synapse_count; syn_id++){ + compartment_to_current[this->compartment_association[syn_id]] += this->i_tot_{{synapse_name}}[syn_id]; + } +} + +std::vector< double > nest::{{synapse_name}}{{cm_unique_suffix}}::distribute_shared_vector(std::vector< double > shared_vector){ + std::vector< double > distributed_vector(this->neuron_{{ synapse_name }}_synapse_count, 0.0); + for(std::size_t syn_id = 0; syn_id < this->neuron_{{ synapse_name }}_synapse_count; syn_id++){ + distributed_vector[syn_id] = shared_vector[compartment_association[syn_id]]; + } + return distributed_vector; +} + +// {{synapse_name}} synapse end /////////////////////////////////////////////////////////// +{%- endfor %} + + +////////////////////////////////////// continuous inputs + +{%- for continuous_name, continuous_info in con_in_info.items() %} +// {{continuous_name}} continuous input /////////////////////////////////////////////////// + +void nest::{{continuous_name}}{{cm_unique_suffix}}::new_continuous_input(std::size_t comp_ass, const long con_in_index) +{ + neuron_{{ continuous_name }}_continuous_input_count++; + i_tot_{{continuous_name}}.push_back(0); + compartment_association.push_back(comp_ass); + continuous_idx.push_back(con_in_index); + + {%- for pure_variable_name, variable_info in continuous_info["States"].items() %} + // state variable {{pure_variable_name }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+continuous_name+"_continuous_input_count").print(rhs_expression) -}}); + {%- endfor %} + + {% for variable_type, variable_info in continuous_info["Parameters"].items() %} + // parameter {{variable_type }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+continuous_name+"_continuous_input_count").print(rhs_expression) -}}); + {%- endfor %} + + // user declared internals in order they were declared + {%- for internal_name, internal_declaration in continuous_info["internals_used_declared"] %} + {{internal_name}}.push_back(0); + {%- endfor %} +} + +void nest::{{continuous_name}}{{cm_unique_suffix}}::new_continuous_input(std::size_t comp_ass, const long con_in_index, const DictionaryDatum& con_in_params) +/* update {{continuous_name}} continuous input parameters and states */ +{ + neuron_{{ continuous_name }}_continuous_input_count++; + compartment_association.push_back(comp_ass); + i_tot_{{continuous_name}}.push_back(0); + continuous_idx.push_back(con_in_index); + + {%- for pure_variable_name, variable_info in continuous_info["States"].items() %} + // state variable {{pure_variable_name }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+continuous_name+"_continuous_input_count").print(rhs_expression) -}}); + {%- endfor %} + {%- for variable_type, variable_info in continuous_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + if( con_in_params->known( "{{variable.name}}" ) ) + {{variable.name}}[neuron_{{ continuous_name }}_continuous_input_count-1] = getValue< double >( con_in_params, "{{variable.name}}" ); + {%- endfor %} + + {%- with %} + {%- for variable_type, variable_info in continuous_info["ODEs"].items() %} + {%- set variable_name = variable_type %} + {%- set dynamic_variable = render_dynamic_channel_variable_name(variable_type, ion_channel_name) %} + // {{continuous_name}} concentration ODE state {{dynamic_variable }} + if( con_in_params->known( "{{variable_name}}" ) ) + {{variable_name}}[neuron_{{ continuous_name }}_continuous_input_count-1] = getValue< double >( con_in_params, "{{variable_name}}" ); + {%- endfor %} + {% endwith %} + + {% for variable_type, variable_info in continuous_info["Parameters"].items() %} + // continuous parameter {{variable_type }} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name}}.push_back({{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("neuron_"+continuous_name+"_continuous_input_count").print(rhs_expression) -}}); + {%- endfor %} + + {%- with %} + {%- for variable_type, variable_info in continuous_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + if( con_in_params->known( "{{variable.name}}" ) ) + {{variable.name}}[neuron_{{ continuous_name }}_continuous_input_count-1] = getValue< double >( con_in_params, "{{variable.name}}" ); + {%- endfor %} + {% endwith %} + + // user declared internals in order they were declared + {%- for internal_name, internal_declaration in continuous_info["internals_used_declared"] %} + {{internal_name}}.push_back(0); + {%- endfor %} +} + +void +nest::{{continuous_name}}{{cm_unique_suffix}}::append_recordables(std::map< Name, double* >* recordables, const long compartment_idx) +{ + for(size_t con_in_id = 0; con_in_id < neuron_{{ continuous_name }}_continuous_input_count; con_in_id++){ + if(compartment_association[con_in_id] == compartment_idx){ + ( *recordables )[ Name( "i_tot_{{continuous_name}}" + std::to_string(con_in_id) )] = &i_tot_{{continuous_name}}[con_in_id]; + } + } +} + +{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} +void nest::{{continuous_name}}{{cm_unique_suffix}}::calibrate() +{%- else %} +void nest::{{continuous_name}}{{cm_unique_suffix}}::pre_run_hook() +{%- endif %} +{ + + // user declared internals in order they were declared + {%- for internal_name, internal_declaration in continuous_info["Internals"] %} + {{internal_name}}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(internal_declaration.get_expression()) }}; + {%- endfor %} + + for(std::size_t i = 0; i < neuron_{{ continuous_name }}_continuous_input_count; i++){ + {% for port_name, port_info in continuous_info["Continuous"].items() %} + {{port_name}}_[i]->clear(); + {% endfor %} + } +} + +std::pair< std::vector< double >, std::vector< double > > nest::{{continuous_name}}{{cm_unique_suffix}}::f_numstep( std::vector< double > v_comp, const long lag {% for ode in continuous_info["Dependencies"]["concentrations"] %}, std::vector< double > {{ode.lhs.name}}{% endfor %}{% if continuous_info["Dependencies"]["receptors"]|length %} + {% endif %}{% for inline in continuous_info["Dependencies"]["receptors"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if continuous_info["Dependencies"]["channels"]|length %} + {% endif %}{% for inline in continuous_info["Dependencies"]["channels"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if continuous_info["Dependencies"]["continuous"]|length %} + {% endif %}{% for inline in continuous_info["Dependencies"]["continuous"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}) +{ + std::vector< double > g_val(neuron_{{ continuous_name }}_continuous_input_count, 0.); + std::vector< double > i_val(neuron_{{ continuous_name }}_continuous_input_count, 0.); + std::vector< double > d_i_tot_dv(neuron_{{ continuous_name }}_continuous_input_count, 0.); + + {% if continuous_info["ODEs"].items()|length %} + std::vector< double > {{ printer_no_origin.print(continuous_info["time_resolution_var"]) }}(neuron_{{ continuous_name }}_continuous_input_count, Time::get_resolution().get_ms()); + {% endif %} + + {%- for ode_variable, ode_info in continuous_info["ODEs"].items() %} + {%- for propagator, propagator_info in ode_info["transformed_solutions"][0]["propagators"].items() %} + std::vector< double > {{ propagator }}(neuron_{{ continuous_name }}_continuous_input_count, 0); + {%- endfor %} + {%- endfor %} + + {% for port_name, port_info in continuous_info["Continuous"].items() %} + std::vector< double > {{ port_name }}(neuron_{{ continuous_name }}_continuous_input_count, 0.); + {% endfor %} + #pragma omp simd + for(std::size_t i = 0; i < neuron_{{ continuous_name }}_continuous_input_count; i++){ + {% for port_name, port_info in continuous_info["Continuous"].items() %} + {{ port_name }}[i] = {{ port_name }}_[i]->get_value( lag ); + {% endfor %} + } + + #pragma omp simd + for(std::size_t i = 0; i < neuron_{{ continuous_name }}_continuous_input_count; i++){ + //update ODE state variable + {%- for ode_variable, ode_info in continuous_info["ODEs"].items() %} + {%- for propagator, propagator_info in ode_info["transformed_solutions"][0]["propagators"].items() %} + {{ propagator }}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(propagator_info["init_expression"]) }}; + {%- endfor %} + {%- for state, state_solution_info in ode_info["transformed_solutions"][0]["states"].items() %} + {{state}}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(state_solution_info["update_expression"]) }}; + {%- endfor %} + {%- endfor %} + + // total current + // this expression should be the transformed inline expression + this->i_tot_{{continuous_name}}[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(continuous_info["root_expression"].get_expression()) }}; + + // derivative of that expression + // voltage derivative of total current + // compute derivative with respect to current with sympy + d_i_tot_dv[i] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("i").print(continuous_info["inline_derivative"]) }}; + + // for numerical integration + g_val[i] = - d_i_tot_dv[i]; + i_val[i] = this->i_tot_{{continuous_name}}[i] - d_i_tot_dv[i] * v_comp[i]; + } + + return std::make_pair(g_val, i_val); + +} + +{%- for function in continuous_info["Functions"] %} +inline {{ function_declaration.FunctionDeclaration(function, "nest::"~continuous_name~cm_unique_suffix~"::") }} +{ +{%- filter indent(2,True) %} +{%- with ast = function.get_block() %} +{%- include "directives/Block.jinja2" %} +{%- endwith %} +{%- endfilter %} +} +{%- endfor %} + +void nest::{{continuous_name}}{{cm_unique_suffix}}::get_currents_per_compartment(std::vector< double >& compartment_to_current){ + for(std::size_t comp_id = 0; comp_id < compartment_to_current.size(); comp_id++){ + compartment_to_current[comp_id] = 0; + } + for(std::size_t con_in_id = 0; con_in_id < neuron_{{ continuous_name }}_continuous_input_count; con_in_id++){ + compartment_to_current[this->compartment_association[con_in_id]] += this->i_tot_{{continuous_name}}[con_in_id]; + } +} + +std::vector< double > nest::{{continuous_name}}{{cm_unique_suffix}}::distribute_shared_vector(std::vector< double > shared_vector){ + std::vector< double > distributed_vector(this->neuron_{{ continuous_name }}_continuous_input_count, 0.0); + for(std::size_t con_in_id = 0; con_in_id < this->neuron_{{ continuous_name }}_continuous_input_count; con_in_id++){ + distributed_vector[con_in_id] = shared_vector[compartment_association[con_in_id]]; + } + return distributed_vector; +} + +// {{continuous_name}} continuous input end /////////////////////////////////////////////////////////// +{%- endfor %} diff --git a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_neuroncurrents_@NEURON_NAME@.h.jinja2 b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_neuroncurrents_@NEURON_NAME@.h.jinja2 new file mode 100644 index 000000000..576abcda0 --- /dev/null +++ b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_neuroncurrents_@NEURON_NAME@.h.jinja2 @@ -0,0 +1,929 @@ +{#- +cm_compartmentcurrents_@NEURON_NAME@.h.jinja2 + +This file is part of NEST. + +Copyright (C) 2004 The NEST Initiative + +NEST is free software: you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation, either version 2 of the License, or +(at your option) any later version. + +NEST is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License +along with NEST. If not, see . +#} + +{%- if tracing %}/* generated by {{self._TemplateReference__context.name}} */ {% endif -%} +{%- import 'directives_cpp/FunctionDeclaration.jinja2' as function_declaration with context %} +#ifndef SYNAPSES_NEAT_H_{{cm_unique_suffix | upper }} +#define SYNAPSES_NEAT_H_{{cm_unique_suffix | upper }} + +#include +#include +#include + +#include "ring_buffer.h" + +{% macro render_variable_type(variable) -%} +{%- with -%} + {%- set symbol = variable.get_scope().resolve_to_symbol(variable.name, SymbolKind.VARIABLE) -%} + {{ types_printer.print(symbol.type_symbol) }} +{%- endwith -%} +{%- endmacro %} + +namespace nest +{ + +///////////////////////////////////// channels + +{%- with %} +{%- for ion_channel_name, channel_info in chan_info.items() %} + +class {{ion_channel_name}}{{cm_unique_suffix}}{ +private: + // states + {%- for pure_variable_name, variable_info in channel_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + std::vector< {{ render_variable_type(variable) }} > {{ variable.name }} = {}; + {%- endfor %} + + // parameters + {%- for pure_variable_name, variable_info in channel_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + std::vector< {{ render_variable_type(variable) }} > {{ variable.name }} = {}; + {%- endfor %} + + // internals + {%- for pure_variable_name, variable_info in channel_info["Internals"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + std::vector< {{ render_variable_type(variable) }} > {{ variable.name }} = {}; + {%- endfor %} + + // ion-channel root-inline value + std::vector< double > i_tot_{{ion_channel_name}} = {}; + + //zero recordable variable in case of zero contribution channel + double zero_recordable = 0; + +public: + // constructor, destructor + {{ion_channel_name}}{{cm_unique_suffix}}(){}; + ~{{ion_channel_name}}{{cm_unique_suffix}}(){}; + + void new_channel(std::size_t comp_ass); + void new_channel(std::size_t comp_ass, const DictionaryDatum& channel_params); + + //number of channels + std::size_t neuron_{{ ion_channel_name }}_channel_count = 0; + + std::vector< size_t > compartment_association = {}; + + // initialization channel +{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} + void calibrate() { +{%- else %} + void pre_run_hook() { +{%- endif %} + }; + + void append_recordables(std::map< Name, double* >* recordables, const long compartment_idx); + + // numerical integration step + std::pair< std::vector< double >, std::vector< double > > f_numstep( std::vector< double > v_comp{% for ode in channel_info["Dependencies"]["concentrations"] %}, std::vector< double > {{ode.lhs.name}}{% endfor %}{% if channel_info["Dependencies"]["receptors"]|length %} + {% endif %}{% for inline in channel_info["Dependencies"]["receptors"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if channel_info["Dependencies"]["channels"]|length %} + {% endif %}{% for inline in channel_info["Dependencies"]["channels"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if channel_info["Dependencies"]["continuous"]|length %} + {% endif %}{% for inline in channel_info["Dependencies"]["continuous"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}); + + // function declarations + +{%- for function in channel_info["Functions"] %} + #pragma omp declare simd + __attribute__((always_inline)) inline {{ function_declaration.FunctionDeclaration(function) }}; +{%- endfor %} + + // root_inline getter + void get_currents_per_compartment(std::vector< double >& compartment_to_current); + + std::vector< double > distribute_shared_vector(std::vector< double > shared_vector); + +}; +{% endfor -%} +{% endwith -%} + + +///////////////////////////////////////////// concentrations + +{%- with %} +{%- for concentration_name, concentration_info in conc_info.items() %} + +class {{ concentration_name }}{{cm_unique_suffix}}{ +private: + // parameters + {%- for pure_variable_name, variable_info in concentration_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + std::vector< {{ render_variable_type(variable) }} > {{ variable.name }} = {}; + {%- endfor %} + + // states + {%- for pure_variable_name, variable_info in concentration_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + std::vector< {{ render_variable_type(variable) }} > {{ variable.name }} = {}; + {%- endfor %} + + // internals + {%- for pure_variable_name, variable_info in concentration_info["Internals"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + std::vector< {{ render_variable_type(variable) }} > {{ variable.name }} = {}; + {%- endfor %} + + // concentration value (root-ode state) + std::vector< double > {{concentration_name}} = {}; + + //zero recordable variable in case of zero contribution concentration + double zero_recordable = 0; + +public: + // constructor, destructor + {{ concentration_name }}{{cm_unique_suffix}}(){}; + ~{{ concentration_name }}{{cm_unique_suffix}}(){}; + + void new_concentration(std::size_t comp_ass); + void new_concentration(std::size_t comp_ass, const DictionaryDatum& concentration_params); + + //number of channels + std::size_t neuron_{{ concentration_name }}_concentration_count = 0; + + std::vector< size_t > compartment_association = {}; + + // initialization concentration +{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} + void calibrate() { +{%- else %} + void pre_run_hook() { +{%- endif %} + for(std::size_t concentration_id = 0; concentration_id < neuron_{{ concentration_name }}_concentration_count; concentration_id++){ + // states + {%- for pure_variable_name, variable_info in concentration_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + {{ variable.name }}[concentration_id] = {{ vector_printer_factory.create_ast_vector_parameter_setter_and_printer("concentration_id").print(rhs_expression) }}; + {%- endfor %} + } + }; + void append_recordables(std::map< Name, double* >* recordables, const long compartment_idx); + + // numerical integration step + void f_numstep( std::vector< double > v_comp{% for ode in concentration_info["Dependencies"]["concentrations"] %}, std::vector< double > {{ode.lhs.name}}{% endfor %}{% if concentration_info["Dependencies"]["receptors"]|length %} + {% endif %}{% for inline in concentration_info["Dependencies"]["receptors"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if concentration_info["Dependencies"]["channels"]|length %} + {% endif %}{% for inline in concentration_info["Dependencies"]["channels"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if concentration_info["Dependencies"]["continuous"]|length %} + {% endif %}{% for inline in concentration_info["Dependencies"]["continuous"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}); + + // function declarations +{%- for function in concentration_info["Functions"] %} + #pragma omp declare simd + __attribute__((always_inline)) inline {{ function_declaration.FunctionDeclaration(function) }}; +{%- endfor %} + + // root_ode getter + void get_concentrations_per_compartment(std::vector< double >& compartment_to_concentration); + + std::vector< double > distribute_shared_vector(std::vector< double > shared_vector); + +}; +{% endfor -%} +{% endwith -%} + + +////////////////////////////////////////////////// synapses + +{% macro render_time_resolution_variable(synapse_info) -%} +{# we assume here that there is only one such variable ! #} +{%- with %} +{%- for analytic_helper_name, analytic_helper_info in synapse_info["analytic_helpers"].items() -%} +{%- if analytic_helper_info["is_time_resolution"] -%} + {{ analytic_helper_name }} +{%- endif -%} +{%- endfor -%} +{% endwith %} +{%- endmacro %} + +{%- with %} +{%- for synapse_name, synapse_info in syns_info.items() %} + +class {{synapse_name}}{{cm_unique_suffix}}{ +private: + // global synapse index + std::vector< long > syn_idx = {}; + + // propagators, initialized via pre_run_hook() or calibrate() + {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} + {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["propagators"].items()%} + std::vector< double > {{state_variable_name}}; + {%- endfor %} + {%- endfor %} + + // kernel state variables, initialized via pre_run_hook() or calibrate() + {%- for convolution, convolution_info in synapse_info["convolutions"].items() %} + {%- for state_variable_name, state_variable_info in convolution_info["analytic_solution"]["kernel_states"].items()%} + std::vector< double > {{state_variable_name}}; + {%- endfor %} + {%- endfor %} + + // user defined parameters, initialized via pre_run_hook() or calibrate() + {%- for param_name, param_declaration in synapse_info["Parameters"].items() %} + std::vector< double > {{param_name}}; + {%- endfor %} + + // states + {%- for pure_variable_name, variable_info in synapse_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + std::vector<{{ render_variable_type(variable) }}> {{ variable.name }} = {} + }; + {%- endfor %} + + std::vector< double > i_tot_{{synapse_name}} = {}; + + // user declared internals in order they were declared, initialized via pre_run_hook() or calibrate() + {%- for internal_name, internal_declaration in synapse_info["internals_used_declared"] %} + std::vector< double > {{internal_name}}; + {%- endfor %} + + + + // spike buffer + std::vector< RingBuffer* > {{synapse_info["buffer_name"]}}_; + +public: + // constructor, destructor + {{synapse_name}}{{cm_unique_suffix}}(){}; + ~{{synapse_name}}{{cm_unique_suffix}}(){}; + + void new_synapse(std::size_t comp_ass, const long syn_index); + void new_synapse(std::size_t comp_ass, const long syn_index, const DictionaryDatum& synapse_params); + + //number of synapses + std::size_t neuron_{{ synapse_name }}_synapse_count = 0; + + std::vector< size_t > compartment_association = {}; + + // numerical integration step + std::pair< std::vector< double >, std::vector< double > > f_numstep( std::vector< double > v_comp, const long lag {% for ode in synapse_info["Dependencies"]["concentrations"] %}, std::vector< double > {{ode.lhs.name}}{% endfor %}{% if synapse_info["Dependencies"]["receptors"]|length %} + {% endif %}{% for inline in synapse_info["Dependencies"]["receptors"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if synapse_info["Dependencies"]["channels"]|length %} + {% endif %}{% for inline in synapse_info["Dependencies"]["channels"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if synapse_info["Dependencies"]["continuous"]|length %} + {% endif %}{% for inline in synapse_info["Dependencies"]["continuous"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}); + + // calibration +{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} + void calibrate(); +{%- else %} + void pre_run_hook(); +{%- endif %} + void append_recordables(std::map< Name, double* >* recordables, const long compartment_idx); + void set_buffer_ptr( std::vector< RingBuffer >& syn_buffers ) + { + for(std::size_t i = 0; i < syn_idx.size(); i++){ + {{synapse_info["buffer_name"]}}_.push_back(&(syn_buffers[syn_idx[i]])); + } + }; + + // function declarations + {%- for function in synapse_info["Functions"] %} + #pragma omp declare simd + __attribute__((always_inline)) inline {{ function_declaration.FunctionDeclaration(function) -}}; + + {% endfor %} + + // root_inline getter + void get_currents_per_compartment(std::vector< double >& compartment_to_current); + + std::vector< double > distribute_shared_vector(std::vector< double > shared_vector); + +}; + +{% endfor -%} +{% endwith -%} + + +////////////////////////////////////////////////// continuous inputs + +{%- with %} +{%- for continuous_name, continuous_info in con_in_info.items() %} + +class {{continuous_name}}{{cm_unique_suffix}}{ +private: + // global continuous input index + std::vector< long > continuous_idx = {}; + + // user defined parameters, initialized via pre_run_hook() or calibrate() + {%- for param_name, param_declaration in continuous_info["Parameters"].items() %} + std::vector< double > {{param_name}}; + {%- endfor %} + + // states + {%- for pure_variable_name, variable_info in continuous_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + {%- set rhs_expression = variable_info["rhs_expression"] %} + std::vector<{{ render_variable_type(variable) }}> {{ variable.name }} = {} + }; + {%- endfor %} + + std::vector< double > i_tot_{{continuous_name}} = {}; + + // user declared internals in order they were declared, initialized via pre_run_hook() or calibrate() + {%- for internal_name, internal_declaration in continuous_info["Internals"] %} + std::vector< double > {{internal_name}}; + {%- endfor %} + + + + // continuous buffer + {% for port_name, port_info in continuous_info["Continuous"].items() %} + std::vector< RingBuffer* > {{ port_name }}_; + {% endfor %} + +public: + // constructor, destructor + {{continuous_name}}{{cm_unique_suffix}}(){}; + ~{{continuous_name}}{{cm_unique_suffix}}(){}; + + void new_continuous_input(std::size_t comp_ass, const long con_in_index); + void new_continuous_input(std::size_t comp_ass, const long con_in_index, const DictionaryDatum& con_in_params); + + //number of continuous inputs + std::size_t neuron_{{ continuous_name }}_continuous_input_count = 0; + + std::vector< size_t > compartment_association = {}; + + // numerical integration step + std::pair< std::vector< double >, std::vector< double > > f_numstep( std::vector< double > v_comp, const long lag {% for ode in continuous_info["Dependencies"]["concentrations"] %}, std::vector< double > {{ode.lhs.name}}{% endfor %}{% if continuous_info["Dependencies"]["receptors"]|length %} + {% endif %}{% for inline in continuous_info["Dependencies"]["receptors"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if continuous_info["Dependencies"]["channels"]|length %} + {% endif %}{% for inline in continuous_info["Dependencies"]["channels"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}{% if continuous_info["Dependencies"]["continuous"]|length %} + {% endif %}{% for inline in continuous_info["Dependencies"]["continuous"] %}, std::vector< double > {{inline.variable_name}}{% endfor %}); + + // calibration +{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} + void calibrate(); +{%- else %} + void pre_run_hook(); +{%- endif %} + void append_recordables(std::map< Name, double* >* recordables, const long compartment_idx); + void set_buffer_ptr( std::vector< RingBuffer >& continuous_buffers ) + { + for(std::size_t i = 0; i < continuous_idx.size(); i++){ + {% for port_name, port_info in continuous_info["Continuous"].items() %} + {{port_name}}_.push_back(&(continuous_buffers[continuous_idx[i]])); + {% endfor %} + } + }; + + // function declarations + {%- for function in continuous_info["Functions"] %} + #pragma omp declare simd + __attribute__((always_inline)) inline {{ function_declaration.FunctionDeclaration(function) -}}; + + {% endfor %} + + // root_inline getter + void get_currents_per_compartment(std::vector< double >& compartment_to_current); + + std::vector< double > distribute_shared_vector(std::vector< double > shared_vector); + +}; + +{% endfor -%} +{% endwith -%} + + +///////////////////////////////////////////// currents + +{%- set channel_suffix = "_chan_" %} +{%- set concentration_suffix = "_conc_" %} +{%- set synapse_suffix = "_syn_" %} +{%- set continuous_suffix = "_con_in_" %} + +class NeuronCurrents{{cm_unique_suffix}} { +private: + //mechanisms + // ion channels +{% with %} + {%- for ion_channel_name, channel_info in chan_info.items() %} + {{ion_channel_name}}{{cm_unique_suffix}} {{ion_channel_name}}{{channel_suffix}}; + {% endfor -%} +{% endwith %} + // concentrations +{% with %} + {%- for concentration_name, concentration_info in conc_info.items() %} + {{concentration_name}}{{cm_unique_suffix}} {{concentration_name}}{{concentration_suffix}}; + {% endfor -%} +{% endwith %} + // synapses +{% with %} + {%- for synapse_name, synapse_info in syns_info.items() %} + {{synapse_name}}{{cm_unique_suffix}} {{synapse_name}}{{synapse_suffix}}; + {% endfor -%} +{% endwith %} + // continuous inputs +{% with %} + {%- for continuous_name, continuous_info in con_in_info.items() %} + {{continuous_name}}{{cm_unique_suffix}} {{continuous_name}}{{continuous_suffix}}; + {% endfor -%} +{% endwith %} + + //number of compartments + std::size_t compartment_number = 0; + + //interdependency shared reference vectors and consecutive area vectors + // ion channels +{% with %} + {%- for ion_channel_name, channel_info in chan_info.items() %} + std::vector < double > {{ion_channel_name}}{{channel_suffix}}_shared_current; + std::vector < std::pair< std::size_t, int > > {{ion_channel_name}}{{channel_suffix}}_con_area; + {% endfor -%} +{% endwith %} + // concentrations +{% with %} + {%- for concentration_name, concentration_info in conc_info.items() %} + std::vector < double > {{concentration_name}}{{concentration_suffix}}_shared_concentration; + std::vector < std::pair< std::size_t, int > > {{concentration_name}}{{concentration_suffix}}_con_area; + {% endfor -%} +{% endwith %} + // synapses +{% with %} + {%- for synapse_name, synapse_info in syns_info.items() %} + std::vector < double > {{synapse_name}}{{synapse_suffix}}_shared_current; + std::vector < std::pair< std::size_t, int > > {{synapse_name}}{{synapse_suffix}}_con_area; + {% endfor -%} +{% endwith %} +// continuous inputs +{% with %} + {%- for continuous_name, continuous_info in con_in_info.items() %} + std::vector < double > {{continuous_name}}{{continuous_suffix}}_shared_current; + std::vector < std::pair< std::size_t, int > > {{continuous_name}}{{continuous_suffix}}_con_area; + {% endfor -%} +{% endwith %} + + //compartment gi states + std::vector < std::pair < double, double > > comps_gi; + +public: + NeuronCurrents{{cm_unique_suffix}}(){}; + ~NeuronCurrents{{cm_unique_suffix}}(){}; + +{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} + void calibrate() { +{%- else %} + void pre_run_hook() { +{%- endif %} + // initialization of ion channels + {%- with %} +{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} + {%- for ion_channel_name, channel_info in chan_info.items() %} + {{ion_channel_name}}{{channel_suffix}}.calibrate(); + {% endfor -%} + {%- for concentration_name, concentration_info in conc_info.items() %} + {{concentration_name}}{{concentration_suffix}}.calibrate(); + {% endfor -%} + {%- for synapse_name, synapse_info in syns_info.items() %} + {{synapse_name}}{{synapse_suffix}}.calibrate(); + {% endfor -%} + {%- for continuous_name, continuous_info in con_in_info.items() %} + {{continuous_name}}{{continuous_suffix}}.calibrate(); + {% endfor -%} +{%- else %} + {%- for ion_channel_name, channel_info in chan_info.items() %} + {{ion_channel_name}}{{channel_suffix}}.pre_run_hook(); + {% endfor -%} + {%- for concentration_name, concentration_info in conc_info.items() %} + {{concentration_name}}{{concentration_suffix}}.pre_run_hook(); + {% endfor -%} + {%- for synapse_name, synapse_info in syns_info.items() %} + {{synapse_name}}{{synapse_suffix}}.pre_run_hook(); + {% endfor -%} + {%- for continuous_name, continuous_info in con_in_info.items() %} + {{continuous_name}}{{continuous_suffix}}.pre_run_hook(); + {% endfor -%} +{%- endif %} + int con_end_index; + {%- for ion_channel_name, channel_info in chan_info.items() %} + if({{ion_channel_name}}{{channel_suffix}}.neuron_{{ ion_channel_name }}_channel_count){ + con_end_index = int({{ion_channel_name}}{{channel_suffix}}.compartment_association[0]); + {{ion_channel_name}}{{channel_suffix}}_con_area.push_back(std::pair< std::size_t, int >(0, con_end_index)); + } + for(std::size_t chan_id = 1; chan_id < {{ion_channel_name}}{{channel_suffix}}.neuron_{{ ion_channel_name }}_channel_count; chan_id++){ + if(!({{ion_channel_name}}{{channel_suffix}}.compartment_association[chan_id] == size_t(int(chan_id) + con_end_index))){ + con_end_index = int({{ion_channel_name}}{{channel_suffix}}.compartment_association[chan_id]) - int(chan_id); + {{ion_channel_name}}{{channel_suffix}}_con_area.push_back(std::pair< std::size_t, int >(chan_id, con_end_index)); + } + } + {% endfor -%} + {%- for concentration_name, concentration_info in conc_info.items() %} + if({{concentration_name}}{{concentration_suffix}}.neuron_{{ concentration_name }}_concentration_count){ + con_end_index = int({{concentration_name}}{{concentration_suffix}}.compartment_association[0]); + {{concentration_name}}{{concentration_suffix}}_con_area.push_back(std::pair< std::size_t, int >(0, con_end_index)); + } + for(std::size_t conc_id = 0; conc_id < {{concentration_name}}{{concentration_suffix}}.neuron_{{ concentration_name }}_concentration_count; conc_id++){ + if(!({{concentration_name}}{{concentration_suffix}}.compartment_association[conc_id] == size_t(int(conc_id) + con_end_index))){ + con_end_index = int({{concentration_name}}{{concentration_suffix}}.compartment_association[conc_id]) - int(conc_id); + {{concentration_name}}{{concentration_suffix}}_con_area.push_back(std::pair< std::size_t, int >(conc_id, con_end_index)); + } + } + {% endfor -%} + {%- for synapse_name, synapse_info in syns_info.items() %} + if({{synapse_name}}{{synapse_suffix}}.neuron_{{ synapse_name }}_synapse_count){ + con_end_index = int({{synapse_name}}{{synapse_suffix}}.compartment_association[0]); + {{synapse_name}}{{synapse_suffix}}_con_area.push_back(std::pair< std::size_t, int >(0, con_end_index)); + } + for(std::size_t syn_id = 0; syn_id < {{synapse_name}}{{synapse_suffix}}.neuron_{{ synapse_name }}_synapse_count; syn_id++){ + if(!({{synapse_name}}{{synapse_suffix}}.compartment_association[syn_id] == size_t(int(syn_id) + con_end_index))){ + con_end_index = int({{synapse_name}}{{synapse_suffix}}.compartment_association[syn_id]) - int(syn_id); + {{synapse_name}}{{synapse_suffix}}_con_area.push_back(std::pair< std::size_t, int >(syn_id, con_end_index)); + } + } + {% endfor -%} + {%- for continuous_name, continuous_info in con_in_info.items() %} + if({{continuous_name}}{{continuous_suffix}}.neuron_{{ continuous_name }}_continuous_input_count){ + con_end_index = int({{continuous_name}}{{continuous_suffix}}.compartment_association[0]); + {{continuous_name}}{{continuous_suffix}}_con_area.push_back(std::pair< std::size_t, int >(0, con_end_index)); + } + for(std::size_t cont_id = 0; cont_id < {{continuous_name}}{{continuous_suffix}}.neuron_{{ continuous_name }}_continuous_input_count; cont_id++){ + if(!({{continuous_name}}{{continuous_suffix}}.compartment_association[cont_id] == size_t(int(cont_id) + con_end_index))){ + con_end_index = int({{continuous_name}}{{continuous_suffix}}.compartment_association[cont_id]) - (cont_id); + {{continuous_name}}{{continuous_suffix}}_con_area.push_back(std::pair< std::size_t, int >(cont_id, con_end_index)); + } + } + {% endfor -%} + {% endwith -%} + }; + + void add_mechanism( const std::string& type, const std::size_t compartment_id, const long multi_mech_index = 0) + { + {%- with %} + bool mech_found = false; + {%- for ion_channel_name, channel_info in chan_info.items() %} + if ( type == "{{ion_channel_name}}" ) + { + {{ion_channel_name}}{{channel_suffix}}.new_channel(compartment_id); + mech_found = true; + } + {% endfor -%} + + {%- for concentration_name, concentration_info in conc_info.items() %} + if ( type == "{{concentration_name}}" ) + { + {{concentration_name}}{{concentration_suffix}}.new_concentration(compartment_id); + mech_found = true; + } + {% endfor -%} + + {%- for synapse_name, synapse_info in syns_info.items() %} + if ( type == "{{synapse_name}}" ) + { + {{synapse_name}}{{synapse_suffix}}.new_synapse(compartment_id, multi_mech_index); + mech_found = true; + } + {% endfor -%} + + {%- for continuous_name, continuous_info in con_in_info.items() %} + if ( type == "{{continuous_name}}" ) + { + {{continuous_name}}{{continuous_suffix}}.new_continuous_input(compartment_id, multi_mech_index); + mech_found = true; + } + {% endfor -%} + + {% endwith -%} + if(!mech_found) + { + assert( false ); + } + }; + + void add_mechanism( const std::string& type, const std::size_t compartment_id, const DictionaryDatum& mechanism_params, const long multi_mech_index = 0) + { + {%- with %} + bool mech_found = false; + {%- for ion_channel_name, channel_info in chan_info.items() %} + if ( type == "{{ion_channel_name}}" ) + { + {{ion_channel_name}}{{channel_suffix}}.new_channel(compartment_id, mechanism_params); + mech_found = true; + } + {% endfor -%} + + {%- for concentration_name, concentration_info in conc_info.items() %} + if ( type == "{{concentration_name}}" ) + { + {{concentration_name}}{{concentration_suffix}}.new_concentration(compartment_id, mechanism_params); + mech_found = true; + } + {% endfor -%} + + {%- for synapse_name, synapse_info in syns_info.items() %} + if ( type == "{{synapse_name}}" ) + { + {{synapse_name}}{{synapse_suffix}}.new_synapse(compartment_id, multi_mech_index, mechanism_params); + mech_found = true; + } + {% endfor -%} + + {%- for continuous_name, continuous_info in con_in_info.items() %} + if ( type == "{{continuous_name}}" ) + { + {{continuous_name}}{{continuous_suffix}}.new_continuous_input(compartment_id, multi_mech_index, mechanism_params); + mech_found = true; + } + {% endfor -%} + {% endwith -%} + if(!mech_found) + { + assert( false ); + } + }; + + void add_compartment(){ + {%- for ion_channel_name, channel_info in chan_info.items() %} + this->add_mechanism("{{ ion_channel_name }}", compartment_number); + {% endfor -%} + + {%- for concentration_name, concentration_info in conc_info.items() %} + this->add_mechanism("{{ concentration_name }}", compartment_number); + {% endfor -%} + + compartment_number++; + + {%- for ion_channel_name, channel_info in chan_info.items() %} + this->{{ion_channel_name}}{{channel_suffix}}_shared_current.push_back(0.0); + {% endfor -%} + + {%- for concentration_name, concentration_info in conc_info.items() %} + this->{{concentration_name}}{{concentration_suffix}}_shared_concentration.push_back(0.0); + {% endfor -%} + + {%- for synapse_name, synapse_info in syns_info.items() %} + this->{{synapse_name}}{{synapse_suffix}}_shared_current.push_back(0.0); + {% endfor -%} + + {%- for continuous_name, continuous_info in con_in_info.items() %} + this->{{continuous_name}}{{continuous_suffix}}_shared_current.push_back(0.0); + {% endfor -%} + + }; + + void add_compartment(const DictionaryDatum& compartment_params){ + {%- for ion_channel_name, channel_info in chan_info.items() %} + this->add_mechanism("{{ ion_channel_name }}", compartment_number, compartment_params); + {% endfor -%} + + {%- for concentration_name, concentration_info in conc_info.items() %} + this->add_mechanism("{{ concentration_name }}", compartment_number, compartment_params); + {% endfor -%} + + compartment_number++; + + {%- for ion_channel_name, channel_info in chan_info.items() %} + this->{{ion_channel_name}}{{channel_suffix}}_shared_current.push_back(0.0); + {% endfor -%} + + {%- for concentration_name, concentration_info in conc_info.items() %} + this->{{concentration_name}}{{concentration_suffix}}_shared_concentration.push_back(0.0); + {% endfor -%} + + {%- for synapse_name, synapse_info in syns_info.items() %} + this->{{synapse_name}}{{synapse_suffix}}_shared_current.push_back(0.0); + {% endfor -%} + + {%- for continuous_name, continuous_info in con_in_info.items() %} + this->{{continuous_name}}{{continuous_suffix}}_shared_current.push_back(0.0); + {% endfor -%} + }; + + void add_receptor_info( ArrayDatum& ad, long compartment_index ) + { + {%- with %} + {%- for synapse_name, synapse_info in syns_info.items() %} + for( std::size_t syn_it = 0; syn_it != {{synapse_name}}{{synapse_suffix}}.neuron_{{synapse_name}}_synapse_count; syn_it++) + { + DictionaryDatum dd = DictionaryDatum( new Dictionary ); + def< long >( dd, names::receptor_idx, syn_it ); + def< long >( dd, names::comp_idx, compartment_index ); + def< std::string >( dd, names::receptor_type, "{{synapse_name}}" ); + ad.push_back( dd ); + } + {% endfor -%} + + {%- for continuous_name, continuous_info in con_in_info.items() %} + for( std::size_t con_it = 0; con_it != {{continuous_name}}{{continuous_suffix}}.neuron_{{continuous_name}}_continuous_input_count; con_it++) + { + DictionaryDatum dd = DictionaryDatum( new Dictionary ); + def< long >( dd, names::receptor_idx, con_it ); + def< long >( dd, names::comp_idx, compartment_index ); + def< std::string >( dd, names::receptor_type, "{{continuous_name}}" ); + ad.push_back( dd ); + } + {% endfor -%} + + {% endwith -%} + }; + + void set_buffers( std::vector< RingBuffer >& buffers) + { + // spike and continuous buffers for synapses and continuous inputs + + {%- with %} + {%- for synapse_name, synapse_info in syns_info.items() %} + {{synapse_name}}{{ synapse_suffix }}.set_buffer_ptr( buffers ); + {% endfor -%} + {%- for continuous_name, continuous_info in con_in_info.items() %} + {{continuous_name}}{{ continuous_suffix }}.set_buffer_ptr( buffers ); + {% endfor -%} + {% endwith %} + + }; + + std::map< Name, double* > get_recordables( const long compartment_idx ) + { + std::map< Name, double* > recordables; + + // append ion channel state variables to recordables + {%- with %} + {%- for ion_channel_name, channel_info in chan_info.items() %} + {{ion_channel_name}}{{channel_suffix}}.append_recordables( &recordables, compartment_idx ); + {% endfor %} + {% endwith %} + + // append concentration state variables to recordables + {%- with %} + {%- for concentration_name, concentration_info in conc_info.items() %} + {{concentration_name}}{{concentration_suffix}}.append_recordables( &recordables, compartment_idx ); + {% endfor %} + {% endwith %} + + // append synapse state variables to recordables + {%- with %} + {%- for synapse_name, synapse_info in syns_info.items() %} + {{synapse_name}}{{synapse_suffix}}.append_recordables( &recordables, compartment_idx ); + {% endfor %} + {% endwith %} + + // append continuous input state variables to recordables + {%- with %} + {%- for continuous_name, continuous_info in con_in_info.items() %} + {{continuous_name}}{{continuous_suffix}}.append_recordables( &recordables, compartment_idx ); + {% endfor %} + {% endwith %} + + return recordables; + }; + + std::vector< std::pair< double, double > > f_numstep( std::vector< double > v_comp_vec, const long lag ) + { + std::vector< std::pair< double, double > > comp_to_gi(compartment_number, std::make_pair(0., 0.)); +{%- for synapse_name, synapse_info in syns_info.items() %} + {{synapse_name}}{{synapse_suffix}}.get_currents_per_compartment({{synapse_name}}{{synapse_suffix}}_shared_current); +{% endfor %} +{%- for continuous_name, continuous_info in con_in_info.items() %} + {{continuous_name}}{{continuous_suffix}}.get_currents_per_compartment({{continuous_name}}{{continuous_suffix}}_shared_current); +{% endfor %} +{%- for concentration_name, concentration_info in conc_info.items() %} + {{ concentration_name }}{{concentration_suffix}}.get_concentrations_per_compartment({{concentration_name}}{{concentration_suffix}}_shared_concentration); +{% endfor -%} +{%- for ion_channel_name, channel_info in chan_info.items() %} + {{ion_channel_name}}{{channel_suffix}}.get_currents_per_compartment({{ion_channel_name}}{{channel_suffix}}_shared_current); +{% endfor -%} + + {%- with %} + {%- for concentration_name, concentration_info in conc_info.items() %} + // computation of {{ concentration_name }} concentration + {{ concentration_name }}{{concentration_suffix}}.f_numstep( {{ concentration_name }}{{concentration_suffix}}.distribute_shared_vector(v_comp_vec){% for ode in concentration_info["Dependencies"]["concentrations"] %}, {{ concentration_name }}{{concentration_suffix}}.distribute_shared_vector({{ode.lhs.name}}{{concentration_suffix}}_shared_concentration){% endfor %}{% if concentration_info["Dependencies"]["receptors"]|length %} + {% endif %}{% for inline in concentration_info["Dependencies"]["receptors"] %}, {{ concentration_name }}{{concentration_suffix}}.distribute_shared_vector({{inline.variable_name}}{{synapse_suffix}}_shared_current){% endfor %}{% if concentration_info["Dependencies"]["channels"]|length %} + {% endif %}{% for inline in concentration_info["Dependencies"]["channels"] %}, {{ concentration_name }}{{concentration_suffix}}.distribute_shared_vector({{inline.variable_name}}{{channel_suffix}}_shared_current){% endfor %}{% if concentration_info["Dependencies"]["continuous"]|length %} + {% endif %}{% for inline in concentration_info["Dependencies"]["continuous"] %}, {{ concentration_name }}{{concentration_suffix}}.distribute_shared_vector({{inline.variable_name}}{{continuous_suffix}}_shared_current){% endfor %}); + + {% endfor -%} + {% endwith -%} + + std::pair< std::vector< double >, std::vector< double > > gi_mech; + std::size_t con_area_count; + + {%- with %} + {%- for ion_channel_name, channel_info in chan_info.items() %} + // contribution of {{ion_channel_name}} channel + gi_mech = {{ion_channel_name}}{{channel_suffix}}.f_numstep( {{ion_channel_name}}{{channel_suffix}}.distribute_shared_vector(v_comp_vec){% for ode in channel_info["Dependencies"]["concentrations"] %}, {{ion_channel_name}}{{channel_suffix}}.distribute_shared_vector({{ode.lhs.name}}{{concentration_suffix}}_shared_concentration){% endfor %}{% if channel_info["Dependencies"]["receptors"]|length %} + {% endif %}{% for inline in channel_info["Dependencies"]["receptors"] %}, {{ion_channel_name}}{{channel_suffix}}.distribute_shared_vector({{inline.variable_name}}{{synapse_suffix}}_shared_current){% endfor %}{% if channel_info["Dependencies"]["channels"]|length %} + {% endif %}{% for inline in channel_info["Dependencies"]["channels"] %}, {{ion_channel_name}}{{channel_suffix}}.distribute_shared_vector({{inline.variable_name}}{{channel_suffix}}_shared_current){% endfor %}{% if channel_info["Dependencies"]["continuous"]|length %} + {% endif %}{% for inline in channel_info["Dependencies"]["continuous"] %}, {{ion_channel_name}}{{channel_suffix}}.distribute_shared_vector({{inline.variable_name}}{{continuous_suffix}}_shared_current){% endfor %}); + + con_area_count = {{ion_channel_name}}{{channel_suffix}}_con_area.size(); + if(con_area_count > 0){ + for(std::size_t con_area_index = 0; con_area_index < con_area_count-1; con_area_index++){ + std::size_t con_area = {{ion_channel_name}}{{channel_suffix}}_con_area[con_area_index].first; + std::size_t next_con_area = {{ion_channel_name}}{{channel_suffix}}_con_area[con_area_index+1].first; + int offset = {{ion_channel_name}}{{channel_suffix}}_con_area[con_area_index].second; + + #pragma omp simd + for(std::size_t chan_id = con_area; chan_id < next_con_area; chan_id++){ + comp_to_gi[chan_id+offset].first += gi_mech.first[chan_id]; + comp_to_gi[chan_id+offset].second += gi_mech.second[chan_id]; + } + } + + std::size_t con_area = {{ion_channel_name}}{{channel_suffix}}_con_area[con_area_count-1].first; + int offset = {{ion_channel_name}}{{channel_suffix}}_con_area[con_area_count-1].second; + + #pragma omp simd + for(std::size_t chan_id = con_area; chan_id < {{ion_channel_name}}{{channel_suffix}}.neuron_{{ ion_channel_name }}_channel_count; chan_id++){ + comp_to_gi[chan_id+offset].first += gi_mech.first[chan_id]; + comp_to_gi[chan_id+offset].second += gi_mech.second[chan_id]; + } + } + {% endfor -%} + {% endwith -%} + + {%- with %} + {%- for synapse_name, synapse_info in syns_info.items() %} + // contribution of {{synapse_name}} synapses + gi_mech = {{synapse_name}}{{synapse_suffix}}.f_numstep( {{synapse_name}}{{synapse_suffix}}.distribute_shared_vector(v_comp_vec), lag {% for ode in synapse_info["Dependencies"]["concentrations"] %}, {{synapse_name}}{{synapse_suffix}}.distribute_shared_vector({{ode.lhs.name}}{{concentration_suffix}}_shared_concentration){% endfor %}{% if synapse_info["Dependencies"]["receptors"]|length %} + {% endif %}{% for inline in synapse_info["Dependencies"]["receptors"] %}, {{synapse_name}}{{synapse_suffix}}.distribute_shared_vector({{inline.variable_name}}{{synapse_suffix}}_shared_current){% endfor %}{% if synapse_info["Dependencies"]["channels"]|length %} + {% endif %}{% for inline in synapse_info["Dependencies"]["channels"] %}, {{synapse_name}}{{synapse_suffix}}.distribute_shared_vector({{inline.variable_name}}{{channel_suffix}}_shared_current){% endfor %}{% if synapse_info["Dependencies"]["continuous"]|length %} + {% endif %}{% for inline in synapse_info["Dependencies"]["continuous"] %}, {{synapse_name}}{{synapse_suffix}}.distribute_shared_vector({{inline.variable_name}}{{continuous_suffix}}_shared_current){% endfor %}); + + con_area_count = {{synapse_name}}{{synapse_suffix}}_con_area.size(); + if(con_area_count > 0){ + for(std::size_t con_area_index = 0; con_area_index < con_area_count-1; con_area_index++){ + std::size_t con_area = {{synapse_name}}{{synapse_suffix}}_con_area[con_area_index].first; + std::size_t next_con_area = {{synapse_name}}{{synapse_suffix}}_con_area[con_area_index+1].first; + int offset = {{synapse_name}}{{synapse_suffix}}_con_area[con_area_index].second; + + #pragma omp simd + for(std::size_t syn_id = con_area; syn_id < next_con_area; syn_id++){ + comp_to_gi[syn_id+offset].first += gi_mech.first[syn_id]; + comp_to_gi[syn_id+offset].second += gi_mech.second[syn_id]; + } + } + + std::size_t con_area = {{synapse_name}}{{synapse_suffix}}_con_area[con_area_count-1].first; + int offset = {{synapse_name}}{{synapse_suffix}}_con_area[con_area_count-1].second; + + #pragma omp simd + for(std::size_t syn_id = con_area; syn_id < {{synapse_name}}{{synapse_suffix}}.neuron_{{ synapse_name }}_synapse_count; syn_id++){ + comp_to_gi[syn_id+offset].first += gi_mech.first[syn_id]; + comp_to_gi[syn_id+offset].second += gi_mech.second[syn_id]; + } + } + {% endfor -%} + {% endwith -%} + + {%- with %} + {%- for continuous_name, continuous_info in con_in_info.items() %} + // contribution of {{continuous_name}} continuous inputs + gi_mech = {{continuous_name}}{{continuous_suffix}}.f_numstep( {{continuous_name}}{{continuous_suffix}}.distribute_shared_vector(v_comp_vec), lag {% for ode in continuous_info["Dependencies"]["concentrations"] %}, {{continuous_name}}{{continuous_suffix}}.distribute_shared_vector({{ode.lhs.name}}{{concentration_suffix}}_shared_concentration){% endfor %}{% if continuous_info["Dependencies"]["receptors"]|length %} + {% endif %}{% for inline in continuous_info["Dependencies"]["receptors"] %}, {{continuous_name}}{{continuous_suffix}}.distribute_shared_vector({{inline.variable_name}}{{synapse_suffix}}_shared_current){% endfor %}{% if continuous_info["Dependencies"]["channels"]|length %} + {% endif %}{% for inline in continuous_info["Dependencies"]["channels"] %}, {{continuous_name}}{{continuous_suffix}}.distribute_shared_vector({{inline.variable_name}}{{channel_suffix}}_shared_current){% endfor %}{% if continuous_info["Dependencies"]["continuous"]|length %} + {% endif %}{% for inline in continuous_info["Dependencies"]["continuous"] %}, {{continuous_name}}{{continuous_suffix}}.distribute_shared_vector({{inline.variable_name}}{{continuous_suffix}}_shared_current){% endfor %}); + + con_area_count = {{continuous_name}}{{continuous_suffix}}_con_area.size(); + if(con_area_count > 0){ + for(std::size_t con_area_index = 0; con_area_index < con_area_count-1; con_area_index++){ + std::size_t con_area = {{continuous_name}}{{continuous_suffix}}_con_area[con_area_index].first; + std::size_t next_con_area = {{continuous_name}}{{continuous_suffix}}_con_area[con_area_index+1].first; + int offset = {{continuous_name}}{{continuous_suffix}}_con_area[con_area_index].second; + + #pragma omp simd + for(std::size_t cont_id = con_area; cont_id < next_con_area; cont_id++){ + comp_to_gi[cont_id+offset].first += gi_mech.first[cont_id]; + comp_to_gi[cont_id+offset].second += gi_mech.second[cont_id]; + } + } + + std::size_t con_area = {{continuous_name}}{{continuous_suffix}}_con_area[con_area_count-1].first; + int offset = {{continuous_name}}{{continuous_suffix}}_con_area[con_area_count-1].second; + + #pragma omp simd + for(std::size_t cont_id = con_area; cont_id < {{continuous_name}}{{continuous_suffix}}.neuron_{{ continuous_name }}_continuous_input_count; cont_id++){ + comp_to_gi[cont_id+offset].first += gi_mech.first[cont_id]; + comp_to_gi[cont_id+offset].second += gi_mech.second[cont_id]; + } + } + {% endfor -%} + {% endwith -%} + + return comp_to_gi; + }; +}; + +} // namespace + +#endif /* #ifndef SYNAPSES_NEAT_H_{{cm_unique_suffix | upper }} */ diff --git a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_tree_@NEURON_NAME@.cpp.jinja2 b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_tree_@NEURON_NAME@.cpp.jinja2 index 38bf6d446..a70f43213 100644 --- a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_tree_@NEURON_NAME@.cpp.jinja2 +++ b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_tree_@NEURON_NAME@.cpp.jinja2 @@ -21,52 +21,69 @@ */ #include "{{neuronSpecificFileNamesCmSyns["tree"]}}.h" - -nest::Compartment{{cm_unique_suffix}}::Compartment{{cm_unique_suffix}}( const long compartment_index, const long parent_index ) +nest::Compartment{{cm_unique_suffix}}::Compartment{{cm_unique_suffix}}( const long compartment_index, const long parent_index) : xx_( 0.0 ) , yy_( 0.0 ) , comp_index( compartment_index ) , p_index( parent_index ) , parent( nullptr ) - , v_comp( 0.0 ) + , v_comp( new double(0.0) ) , ca( 1.0 ) , gc( 0.01 ) , gl( 0.1 ) , el( -70. ) - , gg0( 0.0 ) - , ca__div__dt( 0.0 ) - , gl__div__2( 0.0 ) - , gc__div__2( 0.0 ) - , gl__times__el( 0.0 ) - , ff( 0.0 ) - , gg( 0.0 ) + , gg0( nullptr ) + , ca__div__dt( nullptr ) + , gl__times__el( nullptr ) + , ff( nullptr ) + , gg( nullptr ) , hh( 0.0 ) , n_passed( 0 ) { - v_comp = el; + *v_comp = el; +} - compartment_currents = CompartmentCurrents{{cm_unique_suffix}}(); +nest::Compartment{{cm_unique_suffix}}::Compartment{{cm_unique_suffix}}( const long compartment_index, const long parent_index, + double* v_comp_ref, double* ca__div__dt_ref, double* gl__times__el_ref, double* gg0_ref, double* gg_ref, double* ff_ref ) + : xx_( 0.0 ) + , yy_( 0.0 ) + , comp_index( compartment_index ) + , p_index( parent_index ) + , parent( nullptr ) + , v_comp( v_comp_ref ) + , ca( 1.0 ) + , gc( 0.01 ) + , gl( 0.1 ) + , el( -70. ) + , gg0( gg0_ref ) + , ca__div__dt( ca__div__dt_ref ) + , gl__times__el( gl__times__el_ref ) + , ff( ff_ref ) + , gg( gg_ref ) + , hh( 0.0 ) + , n_passed( 0 ) +{ + *v_comp = el; } nest::Compartment{{cm_unique_suffix}}::Compartment{{cm_unique_suffix}}( const long compartment_index, const long parent_index, - const DictionaryDatum& compartment_params ) + const DictionaryDatum& compartment_params, + double* v_comp_ref, double* ca__div__dt_ref, double* gl__times__el_ref, double* gg0_ref, double* gg_ref, double* ff_ref) : xx_( 0.0 ) , yy_( 0.0 ) , comp_index( compartment_index ) , p_index( parent_index ) , parent( nullptr ) - , v_comp( 0.0 ) + , v_comp( v_comp_ref ) , ca( 1.0 ) , gc( 0.01 ) , gl( 0.1 ) , el( -70. ) - , gg0( 0.0 ) - , ca__div__dt( 0.0 ) - , gl__div__2( 0.0 ) - , gc__div__2( 0.0 ) - , gl__times__el( 0.0 ) - , ff( 0.0 ) - , gg( 0.0 ) + , gg0( gg0_ref ) + , ca__div__dt( ca__div__dt_ref ) + , gl__times__el( gl__times__el_ref ) + , ff( ff_ref ) + , gg( gg_ref ) , hh( 0.0 ) , n_passed( 0 ) { @@ -75,10 +92,10 @@ nest::Compartment{{cm_unique_suffix}}::Compartment{{cm_unique_suffix}}( const lo updateValue< double >( compartment_params, names::g_C, gc ); updateValue< double >( compartment_params, names::g_L, gl ); updateValue< double >( compartment_params, names::e_L, el ); + double v_comp_update = el; + if( compartment_params->known( "v_comp" ) ) updateValue< double >( compartment_params, "v_comp", v_comp_update); - v_comp = el; - - compartment_currents = CompartmentCurrents{{cm_unique_suffix}}( compartment_params ); + *v_comp = v_comp_update; } void @@ -88,73 +105,39 @@ nest::Compartment{{cm_unique_suffix}}::calibrate() nest::Compartment{{cm_unique_suffix}}::pre_run_hook() {%- endif %} { -{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} - compartment_currents.calibrate(); -{%- else %} - compartment_currents.pre_run_hook(); -{%- endif %} const double dt = Time::get_resolution().get_ms(); - ca__div__dt = ca / dt; - gl__div__2 = gl / 2.; - gg0 = ca__div__dt + gl__div__2; - gc__div__2 = gc / 2.; - gl__times__el = gl * el; - - // initialize the buffer - currents.clear(); + // used in vectorized passive loop + *ca__div__dt = ca / dt; + *gl__times__el = gl * el; + *gg0 = *ca__div__dt + gl; } std::map< Name, double* > nest::Compartment{{cm_unique_suffix}}::get_recordables() { - std::map< Name, double* > recordables = compartment_currents.get_recordables( comp_index ); + std::map< Name, double* > recordables; - recordables.insert( recordables.begin(), recordables.end() ); - recordables[ Name( "v_comp" + std::to_string( comp_index ) ) ] = &v_comp; + recordables[ Name( "v_comp" + std::to_string( comp_index ) ) ] = v_comp; return recordables; } // for matrix construction void -nest::Compartment{{cm_unique_suffix}}::construct_matrix_element( const long lag ) +nest::Compartment{{cm_unique_suffix}}::construct_matrix_coupling_elements() { - // matrix diagonal element - gg = gg0; - if ( parent != nullptr ) { - gg += gc__div__2; + *gg += gc; // matrix off diagonal element - hh = -gc__div__2; - } - - for ( auto child_it = children.begin(); child_it != children.end(); ++child_it ) - { - gg += ( *child_it ).gc__div__2; - } - - // right hand side - ff = ( ca__div__dt - gl__div__2 ) * v_comp + gl__times__el; - - if ( parent != nullptr ) - { - ff -= gc__div__2 * ( v_comp - parent->v_comp ); + hh = -gc; } for ( auto child_it = children.begin(); child_it != children.end(); ++child_it ) { - ff -= ( *child_it ).gc__div__2 * ( v_comp - ( *child_it ).v_comp ); + *gg += ( *child_it ).gc; } - - // add all currents to compartment - std::pair< double, double > gi = compartment_currents.f_numstep( v_comp, lag ); - gg += gi.first; - ff += gi.second; - - // add input current - ff += currents.get_value( lag ); } @@ -174,14 +157,116 @@ nest::CompTree{{cm_unique_suffix}}::CompTree{{cm_unique_suffix}}() void nest::CompTree{{cm_unique_suffix}}::add_compartment( const long parent_index ) { - Compartment{{cm_unique_suffix}}* compartment = new Compartment{{cm_unique_suffix}}( size_, parent_index ); + v_comp_vec.push_back(0.0); + ca__div__dt_vec.push_back(0.0); + gl__times__el_vec.push_back(0.0); + gg0_vec.push_back(0.0); + gg_vec.push_back(0.0); + ff_vec.push_back(0.0); + size_t comp_index = v_comp_vec.size()-1; + + Compartment{{cm_unique_suffix}}* compartment = new Compartment{{cm_unique_suffix}}( + size_, parent_index, &(v_comp_vec[comp_index]), &(ca__div__dt_vec[comp_index]), + &(gl__times__el_vec[comp_index]), &(gg0_vec[comp_index]), &(gg_vec[comp_index]), &(ff_vec[comp_index]) + ); + + neuron_currents.add_compartment(); add_compartment( compartment, parent_index ); } void nest::CompTree{{cm_unique_suffix}}::add_compartment( const long parent_index, const DictionaryDatum& compartment_params ) { - Compartment{{cm_unique_suffix}}* compartment = new Compartment{{cm_unique_suffix}}( size_, parent_index, compartment_params ); + //Check whether all passed parameters exist within the used neuron model: + Dictionary* comp_param_copy = new Dictionary(*compartment_params); + + comp_param_copy->remove(names::C_m); + comp_param_copy->remove(names::g_C); + comp_param_copy->remove(names::g_L); + comp_param_copy->remove(names::e_L); + comp_param_copy->remove("v_comp"); +{%- for ion_channel_name, channel_info in chan_info.items() %} + {%- for variable_type, variable_info in channel_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + if( comp_param_copy->known( "{{variable.name}}" ) ) comp_param_copy->remove("{{variable.name}}"); + {%- endfor %} + {%- for variable_type, variable_info in channel_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + if( comp_param_copy->known( "{{variable.name}}" ) ) comp_param_copy->remove("{{variable.name}}"); + {%- endfor %} + {%- for variable_type, variable_info in channel_info["ODEs"].items() %} + if( comp_param_copy->known( "{{variable_type}}" ) ) comp_param_copy->remove("{{variable_type}}"); + {%- endfor %} +{%- endfor %} +{%- for concentration_name, concentration_info in conc_info.items() %} + {%- for variable_type, variable_info in concentration_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + if( comp_param_copy->known( "{{variable.name}}" ) ) comp_param_copy->remove("{{variable.name}}"); + {%- endfor %} +{%- endfor %} +{%- for concentration_name, concentration_info in conc_info.items() %} + {%- for variable_type, variable_info in concentration_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + if( comp_param_copy->known( "{{variable.name}}" ) ) comp_param_copy->remove("{{variable.name}}"); + {%- endfor %} + {%- for variable_type, variable_info in concentration_info["ODEs"].items() %} + if( comp_param_copy->known( "{{variable_type}}" ) ) comp_param_copy->remove("{{variable_type}}"); + {%- endfor %} +{%- endfor %} +{%- for synapse_name, synapse_info in syns_info.items() %} + {%- for variable_type, variable_info in synapse_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + if( comp_param_copy->known( "{{variable.name}}" ) ) comp_param_copy->remove("{{variable.name}}"); + {%- endfor %} +{%- endfor %} +{%- for synapse_name, synapse_info in syns_info.items() %} + {%- for variable_type, variable_info in synapse_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + if( comp_param_copy->known( "{{variable.name}}" ) ) comp_param_copy->remove("{{variable.name}}"); + {%- endfor %} + {%- for variable_type, variable_info in synapse_info["ODEs"].items() %} + if( comp_param_copy->known( "{{variable_type}}" ) ) comp_param_copy->remove("{{variable_type}}"); + {%- endfor %} +{%- endfor %} +{%- for continuous_name, continuous_info in con_in_info.items() %} + {%- for variable_type, variable_info in continuous_info["Parameters"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + if( comp_param_copy->known( "{{variable.name}}" ) ) comp_param_copy->remove("{{variable.name}}"); + {%- endfor %} + {%- for variable_type, variable_info in continuous_info["ODEs"].items() %} + if( comp_param_copy->known( "{{variable_type}}" ) ) comp_param_copy->remove("{{variable_type}}"); + {%- endfor %} +{%- endfor %} +{%- for continuous_name, continuous_info in con_in_info.items() %} + {%- for variable_type, variable_info in continuous_info["States"].items() %} + {%- set variable = variable_info["ASTVariable"] %} + if( comp_param_copy->known( "{{variable.name}}" ) ) comp_param_copy->remove("{{variable.name}}"); + {%- endfor %} +{%- endfor %} + + if(!comp_param_copy->empty()){ + std::string msg = "Following parameters are invalid: "; + for(auto& param : *comp_param_copy){ + msg += param.first.toString(); + msg += "\n"; + } + throw BadParameter(msg); + } + + v_comp_vec.push_back(0.0); + ca__div__dt_vec.push_back(0.0); + gl__times__el_vec.push_back(0.0); + gg0_vec.push_back(0.0); + gg_vec.push_back(0.0); + ff_vec.push_back(0.0); + size_t comp_index = v_comp_vec.size()-1; + + Compartment{{cm_unique_suffix}}* compartment = new Compartment{{cm_unique_suffix}}( + size_, parent_index, compartment_params, &(v_comp_vec[comp_index]), &(ca__div__dt_vec[comp_index]), + &(gl__times__el_vec[comp_index]), &(gg0_vec[comp_index]), &(gg_vec[comp_index]), &(ff_vec[comp_index]) + ); + + neuron_currents.add_compartment(compartment_params); add_compartment( compartment, parent_index ); } @@ -285,6 +370,7 @@ nest::CompTree{{cm_unique_suffix}}::init_pointers() { set_parents(); set_compartments(); + set_compartment_variables(); set_leafs(); } @@ -320,6 +406,23 @@ nest::CompTree{{cm_unique_suffix}}::set_compartments() } } +/** + * Set pointer variables within a compartment + */ +void +nest::CompTree{{cm_unique_suffix}}::set_compartment_variables() +{ + //reset compartment pointers due to unsafe pointers in vectors when resizing during compartment creation + for( size_t i = 0; i < v_comp_vec.size(); i++){ + compartments_[i]->v_comp = &(v_comp_vec[i]); + compartments_[i]->ca__div__dt = &(ca__div__dt_vec[i]); + compartments_[i]->gl__times__el = &(gl__times__el_vec[i]); + compartments_[i]->gg0 = &(gg0_vec[i]); + compartments_[i]->gg = &(gg_vec[i]); + compartments_[i]->ff = &(ff_vec[i]); + } +} + /** * Creates a vector of compartment pointers of compartments that are also leafs of the tree. */ @@ -342,10 +445,7 @@ nest::CompTree{{cm_unique_suffix}}::set_leafs() void nest::CompTree{{cm_unique_suffix}}::set_syn_buffers( std::vector< RingBuffer >& syn_buffers ) { - for ( auto compartment_it = compartments_.begin(); compartment_it != compartments_.end(); ++compartment_it ) - { - ( *compartment_it )->compartment_currents.set_syn_buffers( syn_buffers ); - } + neuron_currents.set_buffers( syn_buffers ); } /** @@ -362,7 +462,11 @@ nest::CompTree{{cm_unique_suffix}}::get_recordables() */ for ( auto compartment_it = compartments_.begin(); compartment_it != compartments_.end(); ++compartment_it ) { - std::map< Name, double* > recordables_comp = ( *compartment_it )->get_recordables(); + long comp_index = (*compartment_it)->comp_index; + std::map< Name, double* > recordables_comp = neuron_currents.get_recordables( comp_index ); + recordables.insert( recordables_comp.begin(), recordables_comp.end() ); + + recordables_comp = ( *compartment_it )->get_recordables(); recordables.insert( recordables_comp.begin(), recordables_comp.end() ); } return recordables; @@ -393,6 +497,11 @@ nest::CompTree{{cm_unique_suffix}}::pre_run_hook() ( *compartment_it )->pre_run_hook(); {%- endif %} } +{%- if nest_version.startswith("v2") or nest_version.startswith("v3.1") or nest_version.startswith("v3.2") or nest_version.startswith("v3.3") %} + neuron_currents.calibrate(); +{%- else %} + neuron_currents.pre_run_hook(); +{%- endif %} } /** @@ -401,12 +510,7 @@ nest::CompTree{{cm_unique_suffix}}::pre_run_hook() std::vector< double > nest::CompTree{{cm_unique_suffix}}::get_voltage() const { - std::vector< double > v_comps; - for ( auto compartment_it = compartments_.cbegin(); compartment_it != compartments_.cend(); ++compartment_it ) - { - v_comps.push_back( ( *compartment_it )->v_comp ); - } - return v_comps; + return v_comp_vec; } /** @@ -415,7 +519,7 @@ nest::CompTree{{cm_unique_suffix}}::get_voltage() const double nest::CompTree{{cm_unique_suffix}}::get_compartment_voltage( const long compartment_index ) { - return compartments_[ compartment_index ]->v_comp; + return *(compartments_[ compartment_index ]->v_comp); } /** @@ -424,9 +528,23 @@ nest::CompTree{{cm_unique_suffix}}::get_compartment_voltage( const long compartm void nest::CompTree{{cm_unique_suffix}}::construct_matrix( const long lag ) { + // compute all channel currents, receptor currents, and input currents + std::vector< std::pair< double, double > > comps_gi = neuron_currents.f_numstep( v_comp_vec, lag ); + + #pragma omp simd + for( size_t i = 0; i < v_comp_vec.size(); i++ ){ + // passive current left hand side + gg_vec[i] = gg0_vec[i]; + // passive currents right hand side + ff_vec[i] = ca__div__dt_vec[i] * v_comp_vec[i] + gl__times__el_vec[i]; + + // add all currents to compartment + gg_vec[i] += comps_gi[i].first; + ff_vec[i] += comps_gi[i].second; + } for ( auto compartment_it = compartments_.begin(); compartment_it != compartments_.end(); ++compartment_it ) { - ( *compartment_it )->construct_matrix_element( lag ); + ( *compartment_it )->construct_matrix_coupling_elements(); } } diff --git a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_tree_@NEURON_NAME@.h.jinja2 b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_tree_@NEURON_NAME@.h.jinja2 index fe1942c03..294970755 100644 --- a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_tree_@NEURON_NAME@.h.jinja2 +++ b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/cm_tree_@NEURON_NAME@.h.jinja2 @@ -56,7 +56,7 @@ along with NEST. If not, see . #include "ring_buffer.h" // compartmental model -#include "{{neuronSpecificFileNamesCmSyns["compartmentcurrents"]}}.h" +#include "{{neuronSpecificFileNamesCmSyns["neuroncurrents"]}}.h" // Includes from libnestutil: #include "dict_util.h" @@ -71,6 +71,8 @@ along with NEST. If not, see . #include "dict.h" #include "dictutils.h" +#include + namespace nest { @@ -90,34 +92,33 @@ public: // tree structure indices Compartment{{cm_unique_suffix}}* parent; std::vector< Compartment{{cm_unique_suffix}} > children; - // vector for synapses - CompartmentCurrents{{cm_unique_suffix}} compartment_currents; // buffer for currents RingBuffer currents; // voltage variable - double v_comp; + double* v_comp; // electrical parameters double ca; // compartment capacitance [uF] double gc; // coupling conductance with parent (meaningless if root) [uS] double gl; // leak conductance of compartment [uS] double el; // leak current reversal potential [mV] // auxiliary variables for efficienchy - double gg0; - double ca__div__dt; - double gl__div__2; - double gc__div__2; - double gl__times__el; + double* gg0; + double* ca__div__dt; + double* gl__times__el; // for numerical integration - double ff; - double gg; + double* ff; + double* gg; double hh; // passage counter for recursion int n_passed; // constructor, destructor - Compartment{{cm_unique_suffix}}( const long compartment_index, const long parent_index ); - Compartment{{cm_unique_suffix}}( const long compartment_index, const long parent_index, const DictionaryDatum& compartment_params ); + Compartment{{cm_unique_suffix}}( const long compartment_index, const long parent_index); + Compartment{{cm_unique_suffix}}( const long compartment_index, const long parent_index, + double* v_comp_ref, double* ca__div__dt_ref, double* gl__times__el_ref, double* gg0_ref, double* gg_ref, double* ff_ref); + Compartment{{cm_unique_suffix}}( const long compartment_index, const long parent_index, const DictionaryDatum& compartment_params, + double* v_comp_ref, double* ca__div__dt_ref, double* gl__times__el_ref, double* gg0_ref, double* gg_ref, double* ff_ref); ~Compartment{{cm_unique_suffix}}(){}; // initialization @@ -129,7 +130,7 @@ public: std::map< Name, double* > get_recordables(); // matrix construction - void construct_matrix_element( const long lag ); + void construct_matrix_coupling_elements(); // maxtrix inversion inline void gather_input( const std::pair< double, double >& in ); @@ -137,7 +138,6 @@ public: inline double calc_v( const double v_in ); }; // Compartment - /* Short helper functions for solving the matrix equation. Can hopefully be inlined */ @@ -151,12 +151,12 @@ inline std::pair< double, double > nest::Compartment{{cm_unique_suffix}}::io() { // include inputs from child compartments - gg -= xx_; - ff -= yy_; + *gg -= xx_; + *ff -= yy_; // output values - double g_val( hh * hh / gg ); - double f_val( ff * hh / gg ); + double g_val( hh * hh / *gg ); + double f_val( *ff * hh / *gg ); return std::make_pair( g_val, f_val ); } @@ -168,9 +168,9 @@ nest::Compartment{{cm_unique_suffix}}::calc_v( const double v_in ) yy_ = 0.0; // compute voltage - v_comp = ( ff - v_in * hh ) / gg; + *v_comp = ( *ff - v_in * hh ) / *gg; - return v_comp; + return *v_comp; } @@ -194,13 +194,24 @@ private: // functions for pointer initialization void set_parents(); void set_compartments(); + void set_compartment_variables(); void set_leafs(); + std::vector< double > v_comp_vec; + std::vector< double > ca__div__dt_vec; + std::vector< double > gl__div__2_vec; + std::vector< double > gl__times__el_vec; + std::vector< double > gg0_vec; + std::vector< double > gg_vec; + std::vector< double > ff_vec; + public: // constructor, destructor CompTree{{cm_unique_suffix}}(); ~CompTree{{cm_unique_suffix}}(){}; + NeuronCurrents{{cm_unique_suffix}} neuron_currents; + // initialization functions for tree structure void add_compartment( const long parent_index ); void add_compartment( const long parent_index, const DictionaryDatum& compartment_params ); diff --git a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/setup/CMakeLists.txt.jinja2 b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/setup/CMakeLists.txt.jinja2 index dc0c1f506..40529e635 100644 --- a/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/setup/CMakeLists.txt.jinja2 +++ b/pynestml/codegeneration/resources_nest_compartmental/cm_neuron/setup/CMakeLists.txt.jinja2 @@ -59,7 +59,7 @@ set( MODULE_NAME ${SHORT_NAME} ) set( MODULE_SOURCES {{moduleName}}.h {{moduleName}}.cpp {%- for neuron in neurons %} - {{perNeuronFileNamesCm[neuron.get_name()]["compartmentcurrents"]}}.cpp {{perNeuronFileNamesCm[neuron.get_name()]["compartmentcurrents"]}}.h + {{perNeuronFileNamesCm[neuron.get_name()]["neuroncurrents"]}}.cpp {{perNeuronFileNamesCm[neuron.get_name()]["neuroncurrents"]}}.h {{perNeuronFileNamesCm[neuron.get_name()]["main"]}}.cpp {{perNeuronFileNamesCm[neuron.get_name()]["main"]}}.h {{perNeuronFileNamesCm[neuron.get_name()]["tree"]}}.cpp {{perNeuronFileNamesCm[neuron.get_name()]["tree"]}}.h {% endfor -%} @@ -255,7 +255,7 @@ if ( BUILD_SHARED_LIBS ) add_library( ${MODULE_NAME}_module MODULE ${MODULE_SOURCES} ) set_target_properties( ${MODULE_NAME}_module PROPERTIES - COMPILE_FLAGS "${NEST_CXXFLAGS} -DLTX_MODULE" + COMPILE_FLAGS "${NEST_CXXFLAGS} -DLTX_MODULE -march=native -ffast-math" LINK_FLAGS "${NEST_LIBS}" PREFIX "" OUTPUT_NAME ${MODULE_NAME} ) diff --git a/pynestml/utils/ast_mechanism_information_collector.py b/pynestml/utils/ast_mechanism_information_collector.py index eed37601d..bf38df32f 100644 --- a/pynestml/utils/ast_mechanism_information_collector.py +++ b/pynestml/utils/ast_mechanism_information_collector.py @@ -101,6 +101,7 @@ def extend_variables_with_initialisations(cls, neuron, mechs_info): neuron.accept(var_init_visitor) mechs_info[mechanism_name]["States"] = var_init_visitor.states mechs_info[mechanism_name]["Parameters"] = var_init_visitor.parameters + mechs_info[mechanism_name]["Internals"] = var_init_visitor.internals return mechs_info @@ -116,6 +117,7 @@ def collect_mechanism_related_definitions(cls, neuron, mechs_info): neuron.accept(variable_collector) global_states = variable_collector.all_states global_parameters = variable_collector.all_parameters + global_internals = variable_collector.all_internals function_collector = ASTFunctionCollectorVisitor() neuron.accept(function_collector) @@ -133,16 +135,23 @@ def collect_mechanism_related_definitions(cls, neuron, mechs_info): neuron.accept(kernel_collector) global_kernels = kernel_collector.all_kernels + continuous_input_collector = ASTContinuousInputDeclarationVisitor() + neuron.accept(continuous_input_collector) + global_continuous_inputs = continuous_input_collector.ports + mechanism_states = list() mechanism_parameters = list() + mechanism_internals = list() mechanism_functions = list() mechanism_inlines = list() mechanism_odes = list() synapse_kernels = list() + mechanism_continuous_inputs = list() mechanism_dependencies = defaultdict() mechanism_dependencies["concentrations"] = list() mechanism_dependencies["channels"] = list() mechanism_dependencies["receptors"] = list() + mechanism_dependencies["continuous"] = list() mechanism_inlines.append(mechs_info[mechanism_name]["root_expression"]) @@ -200,6 +209,10 @@ def collect_mechanism_related_definitions(cls, neuron, mechs_info): if not inline.variable_name in [i.variable_name for i in mechanism_dependencies["receptors"]]: mechanism_dependencies["receptors"].append(inline) + if "continuous" in [e.name for e in inline.get_decorators()]: + if not inline.variable_name in [i.variable_name for i in + mechanism_dependencies["continuous"]]: + mechanism_dependencies["continuous"].append(inline) if not is_dependency: mechanism_inlines.append(inline) @@ -252,6 +265,10 @@ def collect_mechanism_related_definitions(cls, neuron, mechs_info): if variable.name == parameter.name: mechanism_parameters.append(parameter) + for internal in global_internals: + if variable.name == internal.name: + mechanism_internals.append(internal) + for kernel in global_kernels: if variable.name == kernel.get_variables()[0].name: synapse_kernels.append(kernel) @@ -268,15 +285,21 @@ def collect_mechanism_related_definitions(cls, neuron, mechs_info): local_function_call_collector.all_function_calls, search_functions + found_functions) + for input in global_continuous_inputs: + if variable.name == input.name: + mechanism_continuous_inputs.append(input) + search_variables.remove(variable) found_variables.append(variable) # IMPLEMENT CATCH NONDEFINED!!! mechs_info[mechanism_name]["States"] = mechanism_states mechs_info[mechanism_name]["Parameters"] = mechanism_parameters + mechs_info[mechanism_name]["Internals"] = mechanism_internals mechs_info[mechanism_name]["Functions"] = mechanism_functions mechs_info[mechanism_name]["SecondaryInlineExpressions"] = mechanism_inlines mechs_info[mechanism_name]["ODEs"] = mechanism_odes + mechs_info[mechanism_name]["Continuous"] = mechanism_continuous_inputs mechs_info[mechanism_name]["Dependencies"] = mechanism_dependencies return mechs_info @@ -312,9 +335,11 @@ def __init__(self, channel_info): self.inside_declaration = False self.inside_parameter_block = False self.inside_state_block = False + self.inside_internal_block = False self.current_declaration = None self.states = defaultdict() self.parameters = defaultdict() + self.internals = defaultdict() self.channel_info = channel_info def visit_declaration(self, node): @@ -330,10 +355,13 @@ def visit_block_with_variables(self, node): self.inside_state_block = True if node.is_parameters: self.inside_parameter_block = True + if node.is_internals: + self.inside_internal_block = True def endvisit_block_with_variables(self, node): self.inside_state_block = False self.inside_parameter_block = False + self.inside_internal_block = False def visit_variable(self, node): self.inside_variable = True @@ -349,6 +377,12 @@ def visit_variable(self, node): self.parameters[node.name]["ASTVariable"] = node.clone() self.parameters[node.name]["rhs_expression"] = self.current_declaration.get_expression() + if self.inside_internal_block and self.inside_declaration: + if any(node.name == variable.name for variable in self.channel_info["Internals"]): + self.internals[node.name] = defaultdict() + self.internals[node.name]["ASTVariable"] = node.clone() + self.internals[node.name]["rhs_expression"] = self.current_declaration.get_expression() + def endvisit_variable(self, node): self.inside_variable = False @@ -374,8 +408,10 @@ def __init__(self): self.inside_block_with_variables = False self.all_states = list() self.all_parameters = list() + self.all_internals = list() self.inside_states_block = False self.inside_parameters_block = False + self.inside_internals_block = False self.all_variables = list() def visit_block_with_variables(self, node): @@ -384,10 +420,13 @@ def visit_block_with_variables(self, node): self.inside_states_block = True if node.is_parameters: self.inside_parameters_block = True + if node.is_internals: + self.inside_internals_block = True def endvisit_block_with_variables(self, node): self.inside_states_block = False self.inside_parameters_block = False + self.inside_internals_block = False self.inside_block_with_variables = False def visit_variable(self, node): @@ -397,6 +436,8 @@ def visit_variable(self, node): self.all_states.append(node.clone()) if self.inside_parameters_block: self.all_parameters.append(node.clone()) + if self.inside_internals_block: + self.all_internals.append(node.clone()) def endvisit_variable(self, node): self.inside_variable = False @@ -456,3 +497,20 @@ def visit_kernel(self, node): def endvisit_kernel(self, node): self.inside_kernel = False + + +class ASTContinuousInputDeclarationVisitor(ASTVisitor): + def __init__(self): + super(ASTContinuousInputDeclarationVisitor, self).__init__() + self.inside_port = False + self.current_port = None + self.ports = list() + + def visit_input_port(self, node): + self.inside_port = True + self.current_port = node + if self.current_port.is_continuous(): + self.ports.append(node.clone()) + + def endvisit_input_port(self, node): + self.inside_port = False diff --git a/pynestml/utils/ast_vector_parameter_setter_and_printer.py b/pynestml/utils/ast_vector_parameter_setter_and_printer.py new file mode 100644 index 000000000..da0ee5076 --- /dev/null +++ b/pynestml/utils/ast_vector_parameter_setter_and_printer.py @@ -0,0 +1,50 @@ +# -*- coding: utf-8 -*- +# +# ast_vector_parameter_setter_and_printer.py +# +# This file is part of NEST. +# +# Copyright (C) 2004 The NEST Initiative +# +# NEST is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 2 of the License, or +# (at your option) any later version. +# +# NEST is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with NEST. If not, see . + +from pynestml.codegeneration.printers.ast_printer import ASTPrinter +from pynestml.codegeneration.printers.nest_variable_printer import NESTVariablePrinter + + +class ASTVectorParameterSetterAndPrinter(ASTPrinter): + def __init__(self): + super(ASTVectorParameterSetterAndPrinter, self).__init__() + self.inside_variable = False + self.vector_parameter = "" + self.printer = None + self.model = None + + def set_vector_parameter(self, node, vector_parameter=None): + self.vector_parameter = vector_parameter + node.accept(self) + + def print(self, node): + assert isinstance(self.printer._simple_expression_printer._variable_printer, NESTVariablePrinter) + + self.printer._simple_expression_printer._variable_printer.cpp_variable_suffix = "" + + if self.vector_parameter: + self.printer._simple_expression_printer._variable_printer.cpp_variable_suffix = "[" + self.vector_parameter + "]" + + text = self.printer.print(node) + + self.printer._simple_expression_printer._variable_printer.cpp_variable_suffix = "" + + return text diff --git a/pynestml/utils/ast_vector_parameter_setter_and_printer_factory.py b/pynestml/utils/ast_vector_parameter_setter_and_printer_factory.py new file mode 100644 index 000000000..13c3b08d5 --- /dev/null +++ b/pynestml/utils/ast_vector_parameter_setter_and_printer_factory.py @@ -0,0 +1,42 @@ +# -*- coding: utf-8 -*- +# +# ast_vector_parameter_setter_and_printer_factory.py +# +# This file is part of NEST. +# +# Copyright (C) 2004 The NEST Initiative +# +# NEST is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 2 of the License, or +# (at your option) any later version. +# +# NEST is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with NEST. If not, see . + +from pynestml.utils.ast_vector_parameter_setter_and_printer import ASTVectorParameterSetterAndPrinter +from pynestml.visitors.ast_visitor import ASTVisitor + +from pynestml.utils.model_parser import ModelParser +from pynestml.visitors.ast_symbol_table_visitor import ASTSymbolTableVisitor +from pynestml.symbol_table.scope import Scope, ScopeType, Symbol, SymbolKind +from pynestml.symbols.variable_symbol import VariableSymbol + + +class ASTVectorParameterSetterAndPrinterFactory: + + def __init__(self, model, printer): + self.printer = printer + self.model = model + + def create_ast_vector_parameter_setter_and_printer(self, vector_parameter=None): + my_printer = ASTVectorParameterSetterAndPrinter() + my_printer.printer = self.printer + my_printer.model = self.model + my_printer.vector_parameter = vector_parameter + return my_printer diff --git a/pynestml/utils/con_in_info_enricher.py b/pynestml/utils/con_in_info_enricher.py new file mode 100644 index 000000000..fa8228ac8 --- /dev/null +++ b/pynestml/utils/con_in_info_enricher.py @@ -0,0 +1,56 @@ +# -*- coding: utf-8 -*- +# +# con_in_info_enricher.py +# +# This file is part of NEST. +# +# Copyright (C) 2004 The NEST Initiative +# +# NEST is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 2 of the License, or +# (at your option) any later version. +# +# NEST is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with NEST. If not, see . + +from pynestml.visitors.ast_symbol_table_visitor import ASTSymbolTableVisitor +from pynestml.utils.mechs_info_enricher import MechsInfoEnricher +from pynestml.utils.model_parser import ModelParser + +import sympy + + +class ConInInfoEnricher(MechsInfoEnricher): + """Class extends MechsInfoEnricher by the computation of the inline derivative. This hasn't been done in the + channel processing because it would cause a circular dependency through the coco checks used by the ModelParser + which we need to use.""" + def __init__(self, params): + super(MechsInfoEnricher, self).__init__(params) + + @classmethod + def enrich_mechanism_specific(cls, neuron, mechs_info): + mechs_info = cls.compute_expression_derivative(mechs_info) + return mechs_info + + @classmethod + def compute_expression_derivative(cls, chan_info): + for ion_channel_name, ion_channel_info in chan_info.items(): + inline_expression = chan_info[ion_channel_name]["root_expression"] + expr_str = str(inline_expression.get_expression()) + sympy_expr = sympy.parsing.sympy_parser.parse_expr(expr_str) + sympy_expr = sympy.diff(sympy_expr, "v_comp") + + ast_expression_d = ModelParser.parse_expression(str(sympy_expr)) + # copy scope of the original inline_expression into the the derivative + ast_expression_d.update_scope(inline_expression.get_scope()) + ast_expression_d.accept(ASTSymbolTableVisitor()) + + chan_info[ion_channel_name]["inline_derivative"] = ast_expression_d + + return chan_info diff --git a/pynestml/utils/conc_info_enricher.py b/pynestml/utils/conc_info_enricher.py index e4ed0507d..f62d4989e 100644 --- a/pynestml/utils/conc_info_enricher.py +++ b/pynestml/utils/conc_info_enricher.py @@ -23,6 +23,7 @@ class ConcInfoEnricher(MechsInfoEnricher): - """Just created for consistency. No more than the base-class enriching needs to be done""" + """Just created for consistency with the rest of the mechanism generation process. No more than the base-class + enriching needs to be done""" def __init__(self, params): super(MechsInfoEnricher, self).__init__(params) diff --git a/pynestml/utils/continuous_input_processing.py b/pynestml/utils/continuous_input_processing.py new file mode 100644 index 000000000..980217fc1 --- /dev/null +++ b/pynestml/utils/continuous_input_processing.py @@ -0,0 +1,43 @@ +# -*- coding: utf-8 -*- +# +# continuous_input_processing.py +# +# This file is part of NEST. +# +# Copyright (C) 2004 The NEST Initiative +# +# NEST is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 2 of the License, or +# (at your option) any later version. +# +# NEST is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with NEST. If not, see . + +import copy + +from pynestml.utils.mechanism_processing import MechanismProcessing + +from collections import defaultdict + + +class ContinuousInputProcessing(MechanismProcessing): + mechType = "continuous_input" + + def __init__(self, params): + super(MechanismProcessing, self).__init__(params) + + @classmethod + def collect_information_for_specific_mech_types(cls, neuron, mechs_info): + for continuous_name, continuous_info in mechs_info.items(): + continuous = defaultdict() + for port in continuous_info["Continuous"]: + continuous[port.name] = copy.deepcopy(port) + mechs_info[continuous_name]["Continuous"] = continuous + + return mechs_info diff --git a/pynestml/utils/mechs_info_enricher.py b/pynestml/utils/mechs_info_enricher.py index c5514bff8..456ece178 100644 --- a/pynestml/utils/mechs_info_enricher.py +++ b/pynestml/utils/mechs_info_enricher.py @@ -22,13 +22,13 @@ from collections import defaultdict from pynestml.meta_model.ast_model import ASTModel -from pynestml.symbols.predefined_functions import PredefinedFunctions -from pynestml.symbols.symbol import SymbolKind -from pynestml.utils.ast_utils import ASTUtils -from pynestml.utils.model_parser import ModelParser from pynestml.visitors.ast_parent_visitor import ASTParentVisitor from pynestml.visitors.ast_symbol_table_visitor import ASTSymbolTableVisitor +from pynestml.utils.ast_utils import ASTUtils from pynestml.visitors.ast_visitor import ASTVisitor +from pynestml.utils.model_parser import ModelParser +from pynestml.symbols.predefined_functions import PredefinedFunctions +from pynestml.symbols.symbol import SymbolKind class MechsInfoEnricher: diff --git a/pynestml/utils/messages.py b/pynestml/utils/messages.py index f39d9ecdb..3d5673081 100644 --- a/pynestml/utils/messages.py +++ b/pynestml/utils/messages.py @@ -1304,10 +1304,11 @@ def get_integrate_odes_wrong_arg(cls, arg: str): return MessageCode.INTEGRATE_ODES_WRONG_ARG, message @classmethod - def get_mechs_dictionary_info(cls, chan_info, syns_info, conc_info): + def get_mechs_dictionary_info(cls, chan_info, syns_info, conc_info, con_in_info): message = "" message += "chan_info:\n" + chan_info + "\n" message += "syns_info:\n" + syns_info + "\n" message += "conc_info:\n" + conc_info + "\n" + message += "con_in_info:\n" + con_in_info + "\n" return MessageCode.MECHS_DICTIONARY_INFO, message diff --git a/tests/nest_compartmental_tests/concmech_model_test.py b/tests/nest_compartmental_tests/concmech_model_test.py deleted file mode 100644 index 76f3436b1..000000000 --- a/tests/nest_compartmental_tests/concmech_model_test.py +++ /dev/null @@ -1,113 +0,0 @@ -# -*- coding: utf-8 -*- -# -# concmech_model_test.py -# -# This file is part of NEST. -# -# Copyright (C) 2004 The NEST Initiative -# -# NEST is free software: you can redistribute it and/or modify -# it under the terms of the GNU General Public License as published by -# the Free Software Foundation, either version 2 of the License, or -# (at your option) any later version. -# -# NEST is distributed in the hope that it will be useful, -# but WITHOUT ANY WARRANTY; without even the implied warranty of -# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the -# GNU General Public License for more details. -# -# You should have received a copy of the GNU General Public License -# along with NEST. If not, see . - -import os -import pytest - -import nest - -from pynestml.codegeneration.nest_tools import NESTTools -from pynestml.frontend.pynestml_frontend import generate_nest_compartmental_target - -# set to `True` to plot simulation traces -TEST_PLOTS = True -try: - import matplotlib - import matplotlib.pyplot as plt -except BaseException as e: - # always set TEST_PLOTS to False if matplotlib can not be imported - TEST_PLOTS = False - - -class TestCompartmentalConcmech: - @pytest.fixture(scope="module", autouse=True) - def setup(self): - nest.ResetKernel() - nest.SetKernelStatus(dict(resolution=.1)) - - generate_nest_compartmental_target(input_path=os.path.join(os.path.realpath(os.path.dirname(__file__)), "resources", "concmech.nestml"), - suffix="_nestml", - logging_level="DEBUG", - module_name="concmech_mockup_module") - nest.Install("concmech_mockup_module") - - def test_concmech(self): - cm = nest.Create('multichannel_test_model_nestml') - - soma_params = {'C_m': 10.0, 'g_c': 0.0, 'g_L': 1.5, 'e_L': -70.0, 'gbar_Ca_HVA': 1.0, 'gbar_Ca_LVAst': 0.0} - dend_params = {'C_m': 0.1, 'g_c': 0.1, 'g_L': 0.1, 'e_L': -70.0} - - # nest.AddCompartment(cm, 0, -1, soma_params) - cm.compartments = [ - {"parent_idx": -1, "params": soma_params} - # {"parent_idx": 0, "params": dend_params}, - # {"parent_idx": 0, "params": dend_params} - ] - # nest.AddCompartment(cm, 1, 0, dend_params) - # nest.AddCompartment(cm, 2, 0, dend_params) - - # cm.V_th = -50. - - cm.receptors = [ - {"comp_idx": 0, "receptor_type": "AMPA"} - # {"comp_idx": 1, "receptor_type": "AMPA"}, - # {"comp_idx": 2, "receptor_type": "AMPA"} - ] - - # syn_idx_GABA = 0 - # syn_idx_AMPA = 1 - # syn_idx_NMDA = 2 - - # sg1 = nest.Create('spike_generator', 1, {'spike_times': [50., 100., 125., 137., 143., 146., 600.]}) - sg1 = nest.Create('spike_generator', 1, {'spike_times': [100., 1000., 1100., 1200., 1300., 1400., 1500., 1600., 1700., 1800., 1900., 2000., 5000.]}) - # sg1 = nest.Create('spike_generator', 1, {'spike_times': [(item*6000) for item in range(1, 20)]}) - # sg2 = nest.Create('spike_generator', 1, {'spike_times': [115., 155., 160., 162., 170., 254., 260., 272., 278.]}) - # sg3 = nest.Create('spike_generator', 1, {'spike_times': [250., 255., 260., 262., 270.]}) - - nest.Connect(sg1, cm, syn_spec={'synapse_model': 'static_synapse', 'weight': 4.0, 'delay': 0.5, 'receptor_type': 0}) - # nest.Connect(sg2, cm, syn_spec={'synapse_model': 'static_synapse', 'weight': .2, 'delay': 0.5, 'receptor_type': 1}) - # nest.Connect(sg3, cm, syn_spec={'synapse_model': 'static_synapse', 'weight': .3, 'delay': 0.5, 'receptor_type': 2}) - - mm = nest.Create('multimeter', 1, {'record_from': ['v_comp0', 'c_Ca0', 'i_tot_Ca_LVAst0', 'i_tot_Ca_HVA0'], 'interval': .1}) - - nest.Connect(mm, cm) - - nest.Simulate(6000.) - - res = nest.GetStatus(mm, 'events')[0] - - fig, axs = plt.subplots(5) - - axs[0].plot(res['times'], res['v_comp0'], c='b', label='V_m_0') - axs[1].plot(res['times'], res['i_tot_Ca_LVAst0'], c='r', label='i_Ca_LVAst_0') - axs[1].plot(res['times'], res['i_tot_Ca_HVA0'], c='g', label='i_Ca_HVA_0') - axs[2].plot(res['times'], res['c_Ca0'], c='r', label='c_Ca_0') - - axs[0].set_title('V_m_0') - axs[1].set_title('i_Ca_HVA/LVA_0') - axs[2].set_title('c_Ca_0') - # plt.plot(res['times'], res['v_comp2'], c='g', label='V_m_2') - - axs[0].legend() - axs[1].legend() - axs[2].legend() - - plt.savefig("concmech_test.png") diff --git a/tests/nest_compartmental_tests/resources/cm_default.nestml b/tests/nest_compartmental_tests/resources/cm_default.nestml index 83e4b5e44..dadf1dab3 100644 --- a/tests/nest_compartmental_tests/resources/cm_default.nestml +++ b/tests/nest_compartmental_tests/resources/cm_default.nestml @@ -4,7 +4,29 @@ Example compartmental model for NESTML Description +++++++++++ Corresponds to standard compartmental model implemented in NEST. + + +Copyright statement ++++++++++++++++++++ + +This file is part of NEST. + +Copyright (C) 2004 The NEST Initiative + +NEST is free software: you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation, either version 2 of the License, or +(at your option) any later version. + +NEST is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License +along with NEST. If not, see . """ + model cm_default: state: diff --git a/tests/nest_compartmental_tests/resources/concmech.nestml b/tests/nest_compartmental_tests/resources/concmech.nestml index 703c33fbf..486403661 100644 --- a/tests/nest_compartmental_tests/resources/concmech.nestml +++ b/tests/nest_compartmental_tests/resources/concmech.nestml @@ -1,3 +1,24 @@ +""" +Copyright statement ++++++++++++++++++++ + +This file is part of NEST. + +Copyright (C) 2004 The NEST Initiative + +NEST is free software: you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation, either version 2 of the License, or +(at your option) any later version. + +NEST is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License +along with NEST. If not, see . +""" model multichannel_test_model: parameters: @@ -31,7 +52,7 @@ model multichannel_test_model: inf_Ca real = 0.0001 state: - v_comp real = -7500.00000000 + v_comp real = -75.00000000 # state variables Ca_HVA h_Ca_HVA real = 0.69823671 @@ -69,17 +90,17 @@ model multichannel_test_model: h_Ca_LVAst' = ( h_inf_Ca_LVAst(v_comp) - h_Ca_LVAst ) / (tau_h_Ca_LVAst(v_comp)*1s) # equations NaTa_t - #inline NaTa_t real = gbar_NaTa_t * (h_NaTa_t*m_NaTa_t**3) * (e_NaTa_t - v_comp) @mechanism::channel - #m_NaTa_t' = ( m_inf_NaTa_t(v_comp) - m_NaTa_t ) / (tau_m_NaTa_t(v_comp)*1s) - #h_NaTa_t' = ( h_inf_NaTa_t(v_comp) - h_NaTa_t ) / (tau_h_NaTa_t(v_comp)*1s) + inline NaTa_t real = gbar_NaTa_t * (h_NaTa_t*m_NaTa_t**3) * (e_NaTa_t - v_comp) @mechanism::channel + m_NaTa_t' = ( m_inf_NaTa_t(v_comp) - m_NaTa_t ) / (tau_m_NaTa_t(v_comp)*1s) + h_NaTa_t' = ( h_inf_NaTa_t(v_comp) - h_NaTa_t ) / (tau_h_NaTa_t(v_comp)*1s) # equations SKv3_1 #inline SKv3_1 real = gbar_SKv3_1 * (z_SKv3_1) * (e_SKv3_1 - v_comp) @mechanism::channel #z_SKv3_1' = ( z_inf_SKv3_1(v_comp) - z_SKv3_1 ) / (tau_z_SKv3_1(v_comp)*1s) # equations SK_E2 - #inline SK_E2 real = gbar_SK_E2 * (z_SK_E2) * (e_SK_E2 - v_comp) @mechanism::channel - #z_SK_E2' = ( z_inf_SK_E2(c_Ca) - z_SK_E2) / 1.0s + inline SK_E2 real = gbar_SK_E2 * (z_SK_E2) * (e_SK_E2 - v_comp) @mechanism::channel + z_SK_E2' = ( z_inf_SK_E2(c_Ca) - z_SK_E2) / 1.0s # equations Ca concentration mechanism c_Ca' = (inf_Ca - c_Ca) / (tau_Ca*1s) + (gamma_Ca * (Ca_HVA + Ca_LVAst)) / 1s @mechanism::concentration diff --git a/tests/nest_compartmental_tests/resources/continuous_test.nestml b/tests/nest_compartmental_tests/resources/continuous_test.nestml new file mode 100644 index 000000000..2f1387d7e --- /dev/null +++ b/tests/nest_compartmental_tests/resources/continuous_test.nestml @@ -0,0 +1,44 @@ +""" +Copyright statement ++++++++++++++++++++ + +This file is part of NEST. + +Copyright (C) 2004 The NEST Initiative + +NEST is free software: you can redistribute it and/or modify +it under the terms of the GNU General Public License as published by +the Free Software Foundation, either version 2 of the License, or +(at your option) any later version. + +NEST is distributed in the hope that it will be useful, +but WITHOUT ANY WARRANTY; without even the implied warranty of +MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +GNU General Public License for more details. + +You should have received a copy of the GNU General Public License +along with NEST. If not, see . +""" + +model continuous_test_model: + state: + v_comp real = 0 + + parameters: + e_AMPA real = 0 mV + tau_r_AMPA real = 0.2 ms + tau_d_AMPA real = 3.0 ms + + equations: + inline con_in real = (I_stim*10) @mechanism::continuous_input + + kernel g_AMPA = g_norm_AMPA * ( - exp(-t / tau_r_AMPA) + exp(-t / tau_d_AMPA) ) + inline AMPA real = convolve(g_AMPA, spikes_AMPA) * (e_AMPA - v_comp) @mechanism::receptor + + internals: + tp_AMPA real = (tau_r_AMPA * tau_d_AMPA) / (tau_d_AMPA - tau_r_AMPA) * ln( tau_d_AMPA / tau_r_AMPA ) + g_norm_AMPA real = 1. / ( -exp( -tp_AMPA / tau_r_AMPA ) + exp( -tp_AMPA / tau_d_AMPA ) ) + + input: + I_stim real <- continuous + spikes_AMPA <- spike diff --git a/tests/nest_compartmental_tests/cocos_test.py b/tests/nest_compartmental_tests/test__cocos.py similarity index 66% rename from tests/nest_compartmental_tests/cocos_test.py rename to tests/nest_compartmental_tests/test__cocos.py index b3355d8e3..dc4daa28c 100644 --- a/tests/nest_compartmental_tests/cocos_test.py +++ b/tests/nest_compartmental_tests/test__cocos.py @@ -1,6 +1,6 @@ # -*- coding: utf-8 -*- # -# cocos_test.py +# test__cocos.py # # This file is part of NEST. # @@ -22,7 +22,8 @@ from __future__ import print_function import os -import unittest +import pytest + from pynestml.frontend.frontend_configuration import FrontendConfiguration from pynestml.utils.ast_source_location import ASTSourceLocation @@ -35,23 +36,25 @@ from pynestml.utils.model_parser import ModelParser -class CoCosTest(unittest.TestCase): +@pytest.fixture +def setUp(): + Logger.init_logger(LoggingLevel.INFO) + SymbolTable.initialize_symbol_table( + ASTSourceLocation( + start_line=0, + start_column=0, + end_line=0, + end_column=0)) + PredefinedUnits.register_units() + PredefinedTypes.register_types() + PredefinedVariables.register_variables() + PredefinedFunctions.register_functions() + FrontendConfiguration.target_platform = "NEST_COMPARTMENTAL" + - def setUp(self): - Logger.init_logger(LoggingLevel.INFO) - SymbolTable.initialize_symbol_table( - ASTSourceLocation( - start_line=0, - start_column=0, - end_line=0, - end_column=0)) - PredefinedUnits.register_units() - PredefinedTypes.register_types() - PredefinedVariables.register_variables() - PredefinedFunctions.register_functions() - FrontendConfiguration.target_platform = "NEST_COMPARTMENTAL" +class TestCoCos: - def test_invalid_cm_variables_declared(self): + def test_invalid_cm_variables_declared(self, setUp): model = ModelParser.parse_file( os.path.join( os.path.realpath( @@ -59,10 +62,10 @@ def test_invalid_cm_variables_declared(self): os.path.dirname(__file__), 'resources', 'invalid')), 'CoCoCmVariablesDeclared.nestml')) - self.assertEqual(len(Logger.get_all_messages_of_level_and_or_node( - model.get_model_list()[0], LoggingLevel.ERROR)), 5) + assert len(Logger.get_all_messages_of_level_and_or_node( + model.get_model_list()[0], LoggingLevel.ERROR)) == 5 - def test_valid_cm_variables_declared(self): + def test_valid_cm_variables_declared(self, setUp): Logger.set_logging_level(LoggingLevel.INFO) model = ModelParser.parse_file( os.path.join( @@ -71,12 +74,12 @@ def test_valid_cm_variables_declared(self): os.path.dirname(__file__), 'resources', 'valid')), 'CoCoCmVariablesDeclared.nestml')) - self.assertEqual(len(Logger.get_all_messages_of_level_and_or_node( - model.get_model_list()[0], LoggingLevel.ERROR)), 0) + assert len(Logger.get_all_messages_of_level_and_or_node( + model.get_model_list()[0], LoggingLevel.ERROR)) == 0 # it is currently not enforced for the non-cm parameter block, but cm # needs that - def test_invalid_cm_variable_has_rhs(self): + def test_invalid_cm_variable_has_rhs(self, setUp): model = ModelParser.parse_file( os.path.join( os.path.realpath( @@ -84,10 +87,10 @@ def test_invalid_cm_variable_has_rhs(self): os.path.dirname(__file__), 'resources', 'invalid')), 'CoCoCmVariableHasRhs.nestml')) - self.assertEqual(len(Logger.get_all_messages_of_level_and_or_node( - model.get_model_list()[0], LoggingLevel.ERROR)), 2) + assert len(Logger.get_all_messages_of_level_and_or_node( + model.get_model_list()[0], LoggingLevel.ERROR)) == 2 - def test_valid_cm_variable_has_rhs(self): + def test_valid_cm_variable_has_rhs(self, setUp): Logger.set_logging_level(LoggingLevel.INFO) model = ModelParser.parse_file( os.path.join( @@ -96,12 +99,12 @@ def test_valid_cm_variable_has_rhs(self): os.path.dirname(__file__), 'resources', 'valid')), 'CoCoCmVariableHasRhs.nestml')) - self.assertEqual(len(Logger.get_all_messages_of_level_and_or_node( - model.get_model_list()[0], LoggingLevel.ERROR)), 0) + assert len(Logger.get_all_messages_of_level_and_or_node( + model.get_model_list()[0], LoggingLevel.ERROR)) == 0 # it is currently not enforced for the non-cm parameter block, but cm # needs that - def test_invalid_cm_v_comp_exists(self): + def test_invalid_cm_v_comp_exists(self, setUp): model = ModelParser.parse_file( os.path.join( os.path.realpath( @@ -109,10 +112,10 @@ def test_invalid_cm_v_comp_exists(self): os.path.dirname(__file__), 'resources', 'invalid')), 'CoCoCmVcompExists.nestml')) - self.assertEqual(len(Logger.get_all_messages_of_level_and_or_node( - model.get_model_list()[0], LoggingLevel.ERROR)), 4) + assert len(Logger.get_all_messages_of_level_and_or_node( + model.get_model_list()[0], LoggingLevel.ERROR)) == 4 - def test_valid_cm_v_comp_exists(self): + def test_valid_cm_v_comp_exists(self, setUp): Logger.set_logging_level(LoggingLevel.INFO) model = ModelParser.parse_file( os.path.join( @@ -121,5 +124,5 @@ def test_valid_cm_v_comp_exists(self): os.path.dirname(__file__), 'resources', 'valid')), 'CoCoCmVcompExists.nestml')) - self.assertEqual(len(Logger.get_all_messages_of_level_and_or_node( - model.get_model_list()[0], LoggingLevel.ERROR)), 0) + assert len(Logger.get_all_messages_of_level_and_or_node( + model.get_model_list()[0], LoggingLevel.ERROR)) == 0 diff --git a/tests/nest_compartmental_tests/compartmental_model_test.py b/tests/nest_compartmental_tests/test__compartmental_model.py similarity index 96% rename from tests/nest_compartmental_tests/compartmental_model_test.py rename to tests/nest_compartmental_tests/test__compartmental_model.py index 9af75bb43..95e9ff9b7 100644 --- a/tests/nest_compartmental_tests/compartmental_model_test.py +++ b/tests/nest_compartmental_tests/test__compartmental_model.py @@ -1,6 +1,6 @@ # -*- coding: utf-8 -*- # -# compartmental_model_test.py +# test__compartmental_model.py # # This file is part of NEST. # @@ -23,7 +23,6 @@ import os import copy import pytest -import unittest import nest @@ -76,7 +75,7 @@ } -class CMTest(unittest.TestCase): +class TestCM(): def reset_nest(self): nest.ResetKernel() @@ -103,10 +102,10 @@ def install_nestml_model(self): generate_nest_compartmental_target( input_path=input_path, - target_path="/tmp/nestml-component/", + target_path=target_path, module_name="cm_defaultmodule", suffix="_nestml", - logging_level="DEBUG" + logging_level="ERROR" ) def get_model(self, reinstall_flag=True): @@ -262,6 +261,8 @@ def run_model(self): @pytest.mark.skipif(NESTTools.detect_nest_version().startswith("v2"), reason="This test does not support NEST 2") def test_compartmental_model(self): + """We numerically compare the output of the standard nest compartmental model to the equivalent nestml + compartmental model""" self.nestml_flag = False recordables_nest = self.get_rec_list() res_act_nest, res_pas_nest = self.run_model() @@ -528,38 +529,35 @@ def test_compartmental_model(self): for var_nest, var_nestml in zip( recordables_nest[:8], recordables_nestml[:8]): if var_nest == "v_comp0": - atol = 0.51 + atol = 1.0 elif var_nest == "v_comp1": - atol = 0.15 + atol = 0.3 else: - atol = 0.01 - self.assertTrue(np.allclose( + atol = 0.02 + assert (np.allclose( res_act_nest[var_nest], res_act_nestml[var_nestml], atol=atol )) for var_nest, var_nestml in zip( recordables_nest[:8], recordables_nestml[:8]): - if var_nest == "v_comp0": - atol = 0.51 - elif var_nest == "v_comp1": - atol = 0.15 - else: - atol = 0.01 - self.assertTrue(np.allclose( - res_pas_nest[var_nest], res_pas_nestml[var_nestml], atol=atol - )) + if not var_nest in ["h_Na_1", "m_Na_1", "n_K_1"]: + if var_nest == "v_comp0": + atol = 1.0 + elif var_nest == "v_comp1": + atol = 0.3 + else: + atol = 0.02 + assert (np.allclose( + res_pas_nest[var_nest], res_pas_nestml[var_nestml], atol=atol + )) # check if synaptic conductances are equal - self.assertTrue( + assert ( np.allclose( res_act_nest['g_r_AN_AMPA_1'] + res_act_nest['g_d_AN_AMPA_1'], res_act_nestml['g_AN_AMPA1'], 5e-3)) - self.assertTrue( + assert ( np.allclose( res_act_nest['g_r_AN_NMDA_1'] + res_act_nest['g_d_AN_NMDA_1'], res_act_nestml['g_AN_NMDA1'], 5e-3)) - - -if __name__ == "__main__": - unittest.main() diff --git a/tests/nest_compartmental_tests/test__concmech_model.py b/tests/nest_compartmental_tests/test__concmech_model.py new file mode 100644 index 000000000..7c4add105 --- /dev/null +++ b/tests/nest_compartmental_tests/test__concmech_model.py @@ -0,0 +1,131 @@ +# -*- coding: utf-8 -*- +# +# test__concmech_model.py +# +# This file is part of NEST. +# +# Copyright (C) 2004 The NEST Initiative +# +# NEST is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 2 of the License, or +# (at your option) any later version. +# +# NEST is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with NEST. If not, see . + +import os + +import pytest + +import nest + +from pynestml.codegeneration.nest_tools import NESTTools +from pynestml.frontend.pynestml_frontend import generate_nest_compartmental_target + +# set to `True` to plot simulation traces +TEST_PLOTS = True +try: + import matplotlib + import matplotlib.pyplot as plt +except BaseException as e: + # always set TEST_PLOTS to False if matplotlib can not be imported + TEST_PLOTS = False + + +class TestCompartmentalConcmech: + @pytest.fixture(scope="module", autouse=True) + def setup(self): + tests_path = os.path.realpath(os.path.dirname(__file__)) + input_path = os.path.join( + tests_path, + "resources", + "concmech.nestml" + ) + target_path = os.path.join( + tests_path, + "target/" + ) + + if not os.path.exists(target_path): + os.makedirs(target_path) + + print( + f"Compiled nestml model 'cm_main_cm_default_nestml' not found, installing in:" + f" {target_path}" + ) + + nest.ResetKernel() + nest.SetKernelStatus(dict(resolution=.1)) + + generate_nest_compartmental_target( + input_path=input_path, + target_path=target_path, + module_name="concmech_mockup_module", + suffix="_nestml", + logging_level="DEBUG" + ) + + nest.Install("concmech_mockup_module.so") + + def test_concmech(self): + """We test the concentration mechanism by comparing the concentration value at a certain critical point in + time to a previously achieved value at this point""" + cm = nest.Create('multichannel_test_model_nestml') + + params = {'C_m': 10.0, 'g_C': 0.0, 'g_L': 1.5, 'e_L': -70.0, 'gbar_Ca_HVA': 1.0, 'gbar_SK_E2': 1.0} + + cm.compartments = [ + {"parent_idx": -1, "params": params} + ] + + cm.receptors = [ + {"comp_idx": 0, "receptor_type": "AMPA"} + ] + + sg1 = nest.Create('spike_generator', 1, {'spike_times': [100.]}) + + nest.Connect(sg1, cm, syn_spec={'synapse_model': 'static_synapse', 'weight': 4.0, 'delay': 0.5, 'receptor_type': 0}) + + mm = nest.Create('multimeter', 1, {'record_from': ['v_comp0', 'c_Ca0', 'i_tot_Ca_LVAst0', 'i_tot_Ca_HVA0', 'i_tot_SK_E20', 'm_Ca_HVA0', 'h_Ca_HVA0'], 'interval': .1}) + + nest.Connect(mm, cm) + + nest.Simulate(1000.) + + res = nest.GetStatus(mm, 'events')[0] + + step_time_delta = res['times'][1] - res['times'][0] + data_array_index = int(200 / step_time_delta) + + expected_conc = 0.03559438228347359 + + fig, axs = plt.subplots(5) + + axs[0].plot(res['times'], res['v_comp0'], c='r', label='V_m_0') + axs[1].plot(res['times'], res['c_Ca0'], c='y', label='c_Ca_0') + axs[2].plot(res['times'], res['i_tot_Ca_HVA0'], c='b', label='i_tot_Ca_HVA0') + axs[3].plot(res['times'], res['i_tot_SK_E20'], c='b', label='i_tot_SK_E20') + axs[4].plot(res['times'], res['m_Ca_HVA0'], c='g', label='gating var m') + axs[4].plot(res['times'], res['h_Ca_HVA0'], c='r', label='gating var h') + + axs[0].set_title('V_m_0') + axs[1].set_title('c_Ca_0') + axs[2].set_title('i_Ca_HVA_0') + axs[3].set_title('i_tot_SK_E20') + axs[4].set_title('gating vars') + + axs[0].legend() + axs[1].legend() + axs[2].legend() + axs[3].legend() + axs[4].legend() + + plt.savefig("concmech test.png") + + assert res['c_Ca0'][data_array_index] == expected_conc, ("the concentration (left) is not as expected (right). (" + str(res['c_Ca0'][data_array_index]) + "!=" + str(expected_conc) + ")") diff --git a/tests/nest_compartmental_tests/test__continuous_input.py b/tests/nest_compartmental_tests/test__continuous_input.py new file mode 100644 index 000000000..6f8a60060 --- /dev/null +++ b/tests/nest_compartmental_tests/test__continuous_input.py @@ -0,0 +1,124 @@ +# -*- coding: utf-8 -*- +# +# test__continuous_input.py +# +# This file is part of NEST. +# +# Copyright (C) 2004 The NEST Initiative +# +# NEST is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 2 of the License, or +# (at your option) any later version. +# +# NEST is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with NEST. If not, see . + +import os + +import pytest + +import nest + +from pynestml.frontend.pynestml_frontend import generate_nest_compartmental_target + +# set to `True` to plot simulation traces +TEST_PLOTS = True +try: + import matplotlib + import matplotlib.pyplot as plt +except BaseException as e: + # always set TEST_PLOTS to False if matplotlib can not be imported + TEST_PLOTS = False + + +class TestContinuousInput: + @pytest.fixture(scope="module", autouse=True) + def setup(self): + tests_path = os.path.realpath(os.path.dirname(__file__)) + input_path = os.path.join( + tests_path, + "resources", + "continuous_test.nestml" + ) + target_path = os.path.join( + tests_path, + "target/" + ) + + if not os.path.exists(target_path): + os.makedirs(target_path) + + print( + f"Compiled nestml model 'cm_main_cm_default_nestml' not found, installing in:" + f" {target_path}" + ) + + nest.ResetKernel() + nest.SetKernelStatus(dict(resolution=.1)) + + generate_nest_compartmental_target( + input_path=input_path, + target_path=target_path, + module_name="continuous_test_module", + suffix="_nestml", + logging_level="DEBUG" + ) + + nest.Install("continuous_test_module.so") + + def test_continuous_input(self): + """We test the continuous input mechanism by just comparing the input current at a certain critical point in + time to a previously achieved value at this point""" + cm = nest.Create('continuous_test_model_nestml') + + soma_params = {'C_m': 10.0, 'g_C': 0.0, 'g_L': 1.5, 'e_L': -70.0} + + cm.compartments = [ + {"parent_idx": -1, "params": soma_params} + ] + + cm.receptors = [ + {"comp_idx": 0, "receptor_type": "con_in"}, + {"comp_idx": 0, "receptor_type": "AMPA"} + ] + + dcg = nest.Create("ac_generator", {"amplitude": 2.0, "start": 200, "stop": 800, "frequency": 20}) + + nest.Connect(dcg, cm, syn_spec={"synapse_model": "static_synapse", "weight": 1.0, "delay": 0.1, "receptor_type": 0}) + + sg1 = nest.Create('spike_generator', 1, {'spike_times': [205]}) + + nest.Connect(sg1, cm, syn_spec={'synapse_model': 'static_synapse', 'weight': 3.0, 'delay': 0.5, 'receptor_type': 1}) + + mm = nest.Create('multimeter', 1, {'record_from': ['v_comp0', 'i_tot_con_in0', 'i_tot_AMPA0'], 'interval': .1}) + + nest.Connect(mm, cm) + + nest.Simulate(1000.) + + res = nest.GetStatus(mm, 'events')[0] + + fig, axs = plt.subplots(2) + + axs[0].plot(res['times'], res['v_comp0'], c='b', label='V_m_0') + axs[1].plot(res['times'], res['i_tot_con_in0'], c='r', label='continuous') + axs[1].plot(res['times'], res['i_tot_AMPA0'], c='g', label='synapse') + + axs[0].set_title('V_m_0') + axs[1].set_title('inputs') + + axs[0].legend() + axs[1].legend() + + plt.savefig("continuous input test.png") + + step_time_delta = res['times'][1] - res['times'][0] + data_array_index = int(212 / step_time_delta) + + assert 19.9 < res['i_tot_con_in0'][data_array_index] < 20.1, ("the current (left) is not close enough to expected (right). (" + str(res['i_tot_con_in0'][data_array_index]) + " != " + "20.0 +- 0.1" + ")") diff --git a/tests/nest_compartmental_tests/test__interaction_with_disabled_mechanism.py b/tests/nest_compartmental_tests/test__interaction_with_disabled_mechanism.py new file mode 100644 index 000000000..8834c5c7f --- /dev/null +++ b/tests/nest_compartmental_tests/test__interaction_with_disabled_mechanism.py @@ -0,0 +1,128 @@ +# -*- coding: utf-8 -*- +# +# test__interaction_with_disabled_mechanism.py +# +# This file is part of NEST. +# +# Copyright (C) 2004 The NEST Initiative +# +# NEST is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 2 of the License, or +# (at your option) any later version. +# +# NEST is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with NEST. If not, see . + +import os + +import pytest + +import nest + +from pynestml.codegeneration.nest_tools import NESTTools +from pynestml.frontend.pynestml_frontend import generate_nest_compartmental_target + +# set to `True` to plot simulation traces +TEST_PLOTS = True +try: + import matplotlib + import matplotlib.pyplot as plt +except BaseException as e: + # always set TEST_PLOTS to False if matplotlib can not be imported + TEST_PLOTS = False + + +class TestCompartmentalMechDisabled(): + @pytest.fixture(scope="module", autouse=True) + def setup(self): + tests_path = os.path.realpath(os.path.dirname(__file__)) + input_path = os.path.join( + tests_path, + "resources", + "concmech.nestml" + ) + target_path = os.path.join( + tests_path, + "target/" + ) + + if not os.path.exists(target_path): + os.makedirs(target_path) + + print( + f"Compiled nestml model 'cm_main_cm_default_nestml' not found, installing in:" + f" {target_path}" + ) + + nest.ResetKernel() + nest.SetKernelStatus(dict(resolution=.1)) + + generate_nest_compartmental_target( + input_path=input_path, + target_path="/tmp/nestml-component/", + module_name="concmech_mockup_module", + suffix="_nestml", + logging_level="DEBUG" + ) + + nest.Install("concmech_mockup_module.so") + + def test_interaction_with_disabled(self): + """We test the interaction of active mechanisms (the concentration in this case) with disabled mechanisms + (zero key parameters) by just comparing the concentration value at a certain critical point in + time to a previously achieved value at this point""" + cm = nest.Create('multichannel_test_model_nestml') + + params = {'C_m': 10.0, 'g_C': 0.0, 'g_L': 1.5, 'e_L': -70.0, 'gbar_Ca_HVA': 0.0, 'gbar_SK_E2': 1.0} + + cm.compartments = [ + {"parent_idx": -1, "params": params} + ] + + cm.receptors = [ + {"comp_idx": 0, "receptor_type": "AMPA"} + ] + + sg1 = nest.Create('spike_generator', 1, {'spike_times': [100.]}) + + nest.Connect(sg1, cm, syn_spec={'synapse_model': 'static_synapse', 'weight': 4.0, 'delay': 0.5, 'receptor_type': 0}) + + mm = nest.Create('multimeter', 1, {'record_from': ['v_comp0', 'c_Ca0', 'i_tot_Ca_LVAst0', 'i_tot_Ca_HVA0', 'i_tot_SK_E20'], 'interval': .1}) + + nest.Connect(mm, cm) + + nest.Simulate(1000.) + + res = nest.GetStatus(mm, 'events')[0] + + step_time_delta = res['times'][1] - res['times'][0] + data_array_index = int(200 / step_time_delta) + + expected_conc = 2.8159902294145262e-05 + + fig, axs = plt.subplots(4) + + axs[0].plot(res['times'], res['v_comp0'], c='r', label='V_m_0') + axs[1].plot(res['times'], res['c_Ca0'], c='y', label='c_Ca_0') + axs[2].plot(res['times'], res['i_tot_Ca_HVA0'], c='b', label='i_tot_Ca_HVA0') + axs[3].plot(res['times'], res['i_tot_SK_E20'], c='b', label='i_tot_SK_E20') + + axs[0].set_title('V_m_0') + axs[1].set_title('c_Ca_0') + axs[2].set_title('i_Ca_HVA_0') + axs[3].set_title('i_tot_SK_E20') + + axs[0].legend() + axs[1].legend() + axs[2].legend() + axs[3].legend() + + plt.savefig("interaction with disabled mechanism test.png") + + assert res['c_Ca0'][data_array_index] == expected_conc, ("the concentration (left) is not as expected (right). (" + str(res['c_Ca0'][data_array_index]) + "!=" + str(expected_conc) + ")") diff --git a/tests/nest_compartmental_tests/test__model_variable_initialization.py b/tests/nest_compartmental_tests/test__model_variable_initialization.py new file mode 100644 index 000000000..3a8b313f7 --- /dev/null +++ b/tests/nest_compartmental_tests/test__model_variable_initialization.py @@ -0,0 +1,131 @@ +# -*- coding: utf-8 -*- +# +# test__model_variable_initialization.py +# +# This file is part of NEST. +# +# Copyright (C) 2004 The NEST Initiative +# +# NEST is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 2 of the License, or +# (at your option) any later version. +# +# NEST is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with NEST. If not, see . + +import os + +import pytest + +import nest + +from pynestml.codegeneration.nest_tools import NESTTools +from pynestml.frontend.pynestml_frontend import generate_nest_compartmental_target + +# set to `True` to plot simulation traces +TEST_PLOTS = True +try: + import matplotlib + import matplotlib.pyplot as plt +except BaseException as e: + # always set TEST_PLOTS to False if matplotlib can not be imported + TEST_PLOTS = False + + +class TestInitialization(): + @pytest.fixture(scope="module", autouse=True) + def setup(self): + nest.ResetKernel() + nest.SetKernelStatus(dict(resolution=0.1)) + + tests_path = os.path.realpath(os.path.dirname(__file__)) + input_path = os.path.join( + tests_path, + "resources", + "concmech.nestml" + ) + target_path = os.path.join( + tests_path, + "target/" + ) + + if not os.path.exists(target_path): + os.makedirs(target_path) + + print( + f"Compiled nestml model 'cm_main_cm_default_nestml' not found, installing in:" + f" {target_path}" + ) + + generate_nest_compartmental_target( + input_path=input_path, + target_path="/tmp/nestml-component/", + module_name="concmech_module", + suffix="_nestml", + logging_level="DEBUG" + ) + + nest.Install("concmech_module.so") + + def test_non_existing_param(self): + nest.ResetKernel() + nest.SetKernelStatus(dict(resolution=0.1)) + nest.Install("concmech_module.so") + + params = {'C_m': 10.0, 'g_C': 0.0, 'g_L': 1., 'e_L': -70.0, 'non_existing': 1.0} + + with pytest.raises(nest.NESTErrors.BadParameter): + cm = nest.Create('multichannel_test_model_nestml') + cm.compartments = [{"parent_idx": -1, "params": params}] + + def test_existing_states(self): + """Testing whether the python initialization of variables works by looking up the variables at the very start of + the simulation. Since the values change dramatically in the very first step, before which we can not sample them + we test whether they are still large enough and not whether they are the same""" + nest.ResetKernel() + nest.SetKernelStatus(dict(resolution=0.1)) + nest.Install("concmech_module.so") + + params = {'C_m': 10.0, 'g_C': 0.0, 'g_L': 1., 'e_L': -70.0, 'gbar_NaTa_t': 1.0, 'h_NaTa_t': 1000.0, 'c_Ca': 1000.0, 'v_comp': 1000.0} + + cm = nest.Create('multichannel_test_model_nestml') + cm.compartments = [{"parent_idx": -1, "params": params}] + + mm = nest.Create('multimeter', 1, { + 'record_from': ['v_comp0', 'c_Ca0', 'h_NaTa_t0'], 'interval': .1}) + + nest.Connect(mm, cm) + + nest.Simulate(1000.) + + res = nest.GetStatus(mm, 'events')[0] + + data_array_index = 0 + + fig, axs = plt.subplots(3) + + axs[0].plot(res['times'], res['v_comp0'], c='r', label='v_comp0') + axs[1].plot(res['times'], res['c_Ca0'], c='y', label='c_Ca0') + axs[2].plot(res['times'], res['h_NaTa_t0'], c='b', label='h_NaTa_t0') + + axs[0].set_title('v_comp0') + axs[1].set_title('c_Ca0') + axs[2].set_title('h_NaTa_t') + + axs[0].legend() + axs[1].legend() + axs[2].legend() + + plt.savefig("initialization test.png") + + assert res['v_comp0'][data_array_index] > 50.0, ("the voltage (left) is not as expected (right). (" + str(res['v_comp0'][data_array_index]) + "<" + str(50.0) + ")") + + assert res['c_Ca0'][data_array_index] > 900.0, ("the concentration (left) is not as expected (right). (" + str(res['c_Ca0'][data_array_index]) + "<" + str(900.0) + ")") + + assert res['h_NaTa_t0'][data_array_index] > 5.0, ("the gating variable state (left) is not as expected (right). (" + str(res['h_NaTa_t0'][data_array_index]) + "<" + str(5.0) + ")") diff --git a/tests/nest_tests/test_gap_junction.py b/tests/nest_tests/test_gap_junction.py index e94a34de6..25ced348d 100644 --- a/tests/nest_tests/test_gap_junction.py +++ b/tests/nest_tests/test_gap_junction.py @@ -23,6 +23,7 @@ import os import pytest import scipy +import scipy.signal import nest @@ -58,7 +59,7 @@ def generate_code(self, neuron_model: str): files = [os.path.join("models", "neurons", neuron_model + ".nestml")] input_path = [os.path.realpath(os.path.join(os.path.dirname(__file__), os.path.join(os.pardir, os.pardir, s))) for s in files] generate_nest_target(input_path=input_path, - logging_level="DEBUG", + logging_level="WARNING", module_name="nestml_gap_" + neuron_model + "_module", suffix="_nestml", codegen_opts=codegen_opts)