How to locate mech var index in the model exported by bbcorewrite?

[ArcheyChen]

Overview of the feature

I tried to ask this question on The NEURON Forum, but the Confirmation of registration blocked me, so I can only ask question here.

I'm trying to edit the coreneuron data exported using function bbcorewrite, but I found that it's impossible to locate the mech var index.

like I have 10 Iclamps, and I want to modify the 3rd Iclamp's amp value, but I can't just edit the 3rd var in the xx_2.dat, because it's already sorted in some magical way.

I tried using the object_id(mech_obj,1) command, but it only return the order that I created the mech objects, it doesn't help.

I tried setting coreneuron.cell_permute = 0 or h.cvode.cache_efficient(0), but the vars will be permuted anyway.

Adding a custom range var (like id) to mech will be a temporary solution, but if the model gets complicated, it's not an efficient way.

How to get the var export order in NEURON? I tried looking for the code but didn't find where the order was changed.

And it doesn't seem like it have a API to export the info.

BTW: I found some macro defines in coreneuron's code synapseID and selected_for_report , it' seems like coreneuron don't know how to find the vars either, it also requires user to

Foreseeable Impact

Add a single API in Python, which can tell the user mech var export order.

[nrnhines]

I don't understand the context. CoreNEURON generally gets the model and data from NEURON whether CoreNEURON is running in file or direct-memory mode. Thus the value of the 3rd IClamp.amp in CoreNEURON is the same as NEURON's value.when it created the files or copied its data into CoreNEURON. bbcorewrite in a mod file is essentially a helper for this activity for the case of mod files that use Random123 or Vector. In direct memory mode there are functions that map a NEURON reference to a CoreNEURON data location to allow copying of trajectory information gathered during a CoreNEURON run back into NEURON (for Vector.record and Graph trajectories). Also when a direct memory mode CoreNEURON run is complete, all model data (including the event queue) is copied back into NEURON. I'm likely missing your point and it may help to further describe operationally when you need to modify the value of mech_obj.amp in CoreNEURON. For example, if I may be permitted to invent a context, I could say I want to force iclamp.amp t to follow a set of step changes described by a pair of Vectors that describe a step function of amp3 discontinuities at time points. In that case I would get it working first in NEURON using Vector.play
and then run with CoreNEURON.

[ArcheyChen]

Oh, maybe I should tell you my use case.

Building neuron model takes time, especially when the model is large, or, building too many times.

I need to run simulation on a same model with different imputs, so building the model again and again is a waste of time.

So, if I use coreNeuron's file mode to load the model, it will save the time, and I can send the model file to multiple servers to run.

The only problem here is that I can't find a way to change the input in the bbcorewrite exported model,and the coreneuron doesn't seem to provide a var change function in filemode, it can only simulate as it is.

It makes me unable to use filemode, I can only use memory mode, which a lot of time wasted in building models.

So, I really need some way to hack the exported model file, so I can skip the building process.

[ArcheyChen]

The workflow I'm looking forward to is:

1. Build a model with NEURON, pre-set several input mechs, and then export the model
2. Make multiple copies of the model, and then modify the parameters of each model differently
3. Multiple simulations can be run at the same time, no need to build model multiple times.

[ArcheyChen]

Thanks for your reply!

First, I'd like to share my idea for a minimally invasive way to allow variable modification in CoreNEURON's file mode.

There are generally two types of variables we want to modify:
(1) Global variables, such as membrane potential v on a node;
(2) Mechanism variables, such as amp in an IClamp instance.

For global variables (e.g., voltage), NEURON already provides the node_index() function to retrieve a unique index. In my experiments, I noticed that in the exported CoreNEURON model, the ordering of voltage values matches exactly with node_index()—no reordering occurs during export. Assuming this remains true in future versions, we can use a simple syntax like:

global node_idx var_name new_val

For example, to set soma(0.5).v = -66.0, and assuming soma(0.5).node_index() == 2, we could write:

global 2 v -66.0

For mechanism variables, assuming the user does not add or remove mechanisms after export, we can introduce a similar API like mech_obj.mech_index(), which would return a stable index corresponding to the mechanism's position in the exported data. Then, to modify the amp of the 4th IClamp (e.g., iclamp.mech_index() == 3) to 0.01, we could use:

IClamp 3 amp 0.01

Here, "IClamp" comes from the mechanism's name string.

In summary, we can define a unified API format like:

global/mech_name  node_index/mech_index  var_name  new_val

Although retrieving the index might not be intuitive for all users, I believe this approach is sufficient for now and keeps implementation complexity low. We only need to maintain a simple mapping from mech_index to its corresponding variable index in the exported binary.

If we want to go further and support more intuitive syntax like soma(0.5).v, then we would need bbcorewrite to export section information into the model files. Currently, section and segment names are lost during export. This would require more substantial changes, and we would also need to address naming mismatches between NEURON's internal section/mechanism names and Python variable names.

[ArcheyChen]

Additionally, I’d like to ask: if I want to modify the NEURON source code myself to figure out the exact export order of mech variables, which part of the codebase should I look into?

I've already written a small program to parse and edit the model files exported by bbcorewrite. The only missing piece is a reliable way to determine the index of a specific variable in the exported binary data.

I’m currently in a bit of a rush to get large-scale simulations running in parallel, so I’d like to explore implementing this myself rather than waiting for an official NEURON update. Any pointers on where to start in the code would be greatly appreciated.

[nrnhines]

The NEURON side communication to and from CoreNEURON is in nrn/src/nrniv/nrncore_write of the master branch. nrn/src/nrniv/nrncore_write.cpp has a lengthy comment and is the seed that calls functions in files of the nrncore_write subfolder. The CoreNEURON side communication from and to NEURON is in nrn/src/coreneuron/io subfolder and the principle file is nrn_setup.cpp. That file has a long comment that describes the data transfer file formats. I would be happy to give you a tour of the code activities that are involved in file transfer. Let's organize a zoom call.

As much as possible I'd like to base this on data handles. Take a look at print(soma(0.5)._ref_v) and print(iclamp[3]._ref_amp) but we also need to bring in the concept of a CellGroup which in some ways is synonymous to a NrnThread, a file, and a set of gids.

Two functions are absolutely critical to comprehending the relationship between NEURON and CoreNEURON data (at least for direct memory transfer mode). On
the NEURON side it is nrn_dblpntr2nrncore and on the CoreNEURON side it is legacy_index2pointer

[ArcheyChen]

Thank you very much for the detailed guidance and your kind offer for a Zoom walkthrough!

Following your suggestion, I’ve been exploring the nrn/src/nrniv/nrncore_write folder and have located the function responsible for exporting mech data: nrnthread_dat2_mech. However, I noticed that this function uses a straightforward loop:

for (auto instance = 0, k = 0; instance < n; ++instance)

This suggests that it’s writing mech instances in order. Yet in the exported data files, the ordering doesn’t match the creation order of the mechanism objects in my NEURON script.

From this, I suspect that the mech object instances must have already been reordered at an earlier stage.
I also checked Memb_list, which, as the comment suggests, is just a view of the mech instances and assumes they are already sorted.

So now I’m trying to trace back further — specifically, I want to understand where the mech instances are created and how they are sorted. It’s clearly not a simple push_back, otherwise the ordering should match the creation sequence.

Would you be able to point me to the part of the NEURON codebase where new mech instances are allocated and sorted into the Memb_list? That would help me a lot in understanding and tracking the export order.

I also briefly looked into legacy_index2pointer and related functions. I understand now that they handle conversion between AOS and SOA memory layouts, which is very helpful for understanding memory access during simulation. However, for my current goal, I believe those are less directly relevant, as I’m only interested in the ordering of the mech instances themselves, not the layout of variables within each instance.

Regarding the Zoom call — thank you again for the offer! Since I’m located in China, the time zone difference may make scheduling a bit tricky, and I also don’t want to take up too much of your time. I believe that if I can get the relevant code pointers, I’ll be able to read through and understand the details on my own.

Thanks again for your support and for making NEURON so accessible!

[ArcheyChen]

Today I found what seems to be the relevant reordering logic in nrnoc/treeset.cpp, inside the nrn_sort_mech_data function. The array mech_data_permutation appears to represent the actual permutation used for reordering mechanism instances.

However, I ran into a few questions:

I’m not sure how many times nrn_sort_mech_data is called during model building. In my test case, it was called only once, so the permutation result seemed clean and correct. But if it can be called multiple times in other scenarios, the final ordering might become difficult to predict or trace.

The mech_data_permutation array is currently temporary—it gets discarded after sorting. I think it would be helpful to persist this data somewhere in the model structure, but I’m not sure which structure would be the most appropriate. Would it make more sense to store it in struct storage, or in struct Soa?

On the Python side, users interact with mechanism instances as hoc_obj. I'm still trying to figure out how to go from a hoc_obj back to the underlying C++ structure (like storage) in order to retrieve the permutation index. That seems to be the key missing link if I want to create a reliable mech_index() API.

Any insights on these issues would be greatly appreciated!

[ArcheyChen]

Oh, I just found that Memb_list::legacy_index(double const* ptr) might be the 'mech index' I'm looking for, maybe I just need to expose this function to Python. But how to do it?

[nrnhines]

Memb_list::legacy_index(double const* ptr) might be the 'mech index' I'm looking for

Possibly. I'm not sure. The function was written by the person who did all of the programming transformation from AOS to SOA and DataHandle in NEURON. There may be more than one way to determine it. I continue to believe, but am not certain, that the mech index you are looking for is in fact given by a data handle reference, e.g. print(soma(0.5)._ref_v). Til this question, I have not had the need to delve into this part of the new code. I will need to do some experiments to see how a datahandle or Memb_list.legacy_index(datahandle) match up with the corresponding value in CoreNEURON.

[nrnhines]

I created a hines/legindex branch to experiment with how to derive enough information from a data handle to allow a change to a variable value when coreneuron is launched in file mode. To get the Memb_list::legacy_index( I have to figure out how to determine the Memb_list associated with the data handle. I added a small test file at the top level source folder so from there, after you build, you can try it with python legacy.py . It is an example of how one can expose an internal c++ function to python. In this case, since ParallelContext presently contains the methods that activate coreneuron, I added it there.

[nrnhines]

NEURON version 8.2.x was fully AOS and there had to be a conversion to CoreNEURON SOA. NEURON version 9.0 is now also SOA and each variable name has it's own std::vector container so much less conversion is needed. A mechanism has many variables and thus many std::vector and also an associated permutation vector. All of this container and data handle stuff is implemented in nrn/src/neuron. For simulation we often try to permute so that access to values in each vector is sequential for each time step. Turns out that for multiple threads, we now just permute the vectors so each thread is contiguous. Also, our pointers (data handles) now have the nice property that they always point to the correct value regardless of what permutation is used and whether backing store has been freed and reallocated. I other words, cell.section(x)._ref_v[0] is forever equal to cell.section(x).v, cell.section(x).hh._ref_m[0] is forever equal to cell.section(x).hh.m, and ICLamp[3]._ref_amp[0] is forever equal to IClamp[3].amp.

I see another message from you has arrived. I'll stop now and continue after pondering it.

[ArcheyChen]

Thank you for your reply. Based on your suggestion, I examined data_handler and found that it indeed has a very useful property: it can be directly printed in Python, and the resulting information includes a row field — which is exactly the index I was looking for.

Leveraging this feature, I wrote a wrapper that is compatible with the current mainline NEURON simulator code. By enumerating the __dict__ of a mech_obj, I look for attributes that contain _ref_{name}. If such an attribute exists, it indicates that this is a variable within the mechanism. I then store this variable into a dictionary and simultaneously record its associated row. With this information, it becomes possible to locate and modify the variable within the files exported by bbcorewrite.

Next, using Python's __setattr__ and __getattr__, I implemented support for accessing variables via wrapper.var = xx. This makes it very user-friendly to write model-modifying code sequences.

The code I pasted afterward was only intended as a technical demonstration. At this stage, the wrapper only serves to create a backup of the variable. When wrapper.var is accessed, it simply accesses this backup variable.

With slight modifications, this technique can be extended so that accessing wrapper.var actually modifies the corresponding value in the exported model file.

import re

class MechWrapper:
    def __init__(self, mech):
        self.__mech = mech
        self.__row_map = {}
        self.__backup = {}
        self.__allowed_vars = set()

        for name in mech.__dict__:
            ref_name = f"_ref_{name}"
            try:
                ref = getattr(mech, ref_name)
                match = re.search(r'row=(\d+)/\d+', str(ref))
                if match:
                    row = int(match.group(1))
                    self.__row_map[name] = row
                    self.__backup[name] = getattr(mech, name)
                    self.__allowed_vars.add(name)
                    self._create_property(name)
            except Exception:
                continue

    def _create_property(self, name):
        def getter(self):
            print(f"Accessing '{name}' at row {self.__row_map[name]}")
            return self.__backup[name]

        def setter(self, value):
            print(f"Updating backup value of '{name}' to {value}")
            self.__backup[name] = value

        setattr(self.__class__, name, property(getter, setter))

    def __setattr__(self, name, value):
        if name.startswith('_'):
            super().__setattr__(name, value)
        elif name in getattr(self, '__allowed_vars', set()):
            object.__setattr__(self, name, value)
        else:
            raise AttributeError(f"Cannot add new attribute '{name}'")

    def __getattr__(self, name):
        if name in self.__allowed_vars:
            return getattr(self, name)
        raise AttributeError(f"Attribute '{name}' not found")

from neuron import h
soma = h.Section()
iclamp = h.IClamp(soma(0.5))
wrapper = MechWrapper(iclamp)

# Accessing a property
print(wrapper.amp)     # Prints row info + returns backup value

# Modifying a property
wrapper.amp = 5.0      # Only modifies backup, does not affect NEURON internal value
print(wrapper.amp)     # Prints row info + returns updated backup value

print(wrapper.__dict__)
wrapper.not_a_real_var = 42  # Attempt to add a nonexistent variable — should raise error
print(wrapper.not_a_real_var)  # Attempt to access nonexistent variable — should raise error

[ArcheyChen]

However, this current technique does come with some limitations.

Firstly, I discovered that in NEURON, changes to morphology only take full effect after calling finitialize(). Moreover, after multiple calls to finitialize() and simulation runs, the mech_obj may be reordered multiple times.

Therefore, it's crucial to strictly follow the workflow:
finitialize -> bbcorewrite -> wrap mech_obj.
Of course, we could further encapsulate this logic by creating a wrapper_mgr that first records object references to each NEURON mech_obj, and then delays initialization of each wrapper until after the model export. For example:

# Build NEURON model
iclamp_wrap = wrap_mgr.create_wrapper(neuron_iclamp_obj)
wrap_mgr.dump_model(path)  # internally calls bbcorewrite + iclamp_wrap.init() to get correct row value 
iclamp_wrap.amp = xx       # generate model modification sequence or other operations

These wrapper classes can then be serialized and deserialized using Python's pickle module, which is convenient for repeated simulation runs.

However, one issue is that this approach only supports modifying mechanisms (mech) that were already created during NEURON's modeling phase. It does not support adding new mechanisms after loading the exported model. Additionally, for operations like accessing or modifying voltage values, it is necessary to specify the target location (e.g., soma(0.5)) during modeling in NEURON, in order to retrieve the corresponding node_index(). After loading the exported model, information such as sections is lost, so it's no longer possible to determine locations like soma(0.8).

In other words, NEURON’s current model export process results in some loss of information. While the use of wrapper classes can partially mitigate this, enabling model editing after loading would require NEURON itself to enhance its model serialization and deserialization mechanisms.

Regarding the Zoom meeting invitation you previously mentioned — I was wondering if it's still valid? If you're available, I’d love to learn more about NEURON’s internal data structure design (such as section, data_handler, and object), and also discuss ways to improve NEURON's model export and import interfaces.

My email is: 464270342@qq.com.

[nrnhines]

I'm glad you've been able to determine sufficient information within NEURON to modify specific range variable instance values in the files read by CoreNEURON (under the limitation that the model structure (as opposed to values within the structure) does not change ). I will contact you by email to set up a zoom meeting.