Skip to content
This repository was archived by the owner on Mar 20, 2023. It is now read-only.

Commit f72026d

Browse files
nrnhinesnrnhines
authored
Extend POINTER transfer to any RANGE variable in a NRN_THREAD (#772)
* Extend POINTER from voltage to any RANGE variable. * trajectory recording is after AFTER_SOLVE and for consistency with NEURON, after BEFORE_STEP as well. * update coreneuron ringtest integration data to version 1.5 * Handle the checkpoint for the POINTER * Initialize reporting interface only if there are reports
1 parent f3d4615 commit f72026d

37 files changed

+311
-77
lines changed

‎coreneuron/apps/main1.cpp‎

Lines changed: 5 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -608,8 +608,11 @@ extern "C" int run_solve_core(int argc, char** argv) {
608608
if (corenrn_param.report_buff_size != corenrn_param.report_buff_size_default) {
609609
set_report_buffer_size(corenrn_param.report_buff_size);
610610
}
611-
setup_report_engine(min_report_dt, delay);
612-
configs.clear();
611+
612+
if (!configs.empty()) {
613+
setup_report_engine(min_report_dt, delay);
614+
configs.clear();
615+
}
613616

614617
// call prcellstate for prcellgid
615618
call_prcellstate_for_prcellgid(corenrn_param.prcellgid, compute_gpu, 0);

‎coreneuron/io/nrn2core_direct.h‎

Lines changed: 3 additions & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -9,6 +9,7 @@
99
#pragma once
1010

1111
#include <iostream>
12+
#include <vector>
1213

1314
extern "C" {
1415
// The callbacks into nrn/src/nrniv/nrnbbcore_write.cpp to get
@@ -56,7 +57,8 @@ extern int (*nrn2core_get_dat2_mech_)(int tid,
5657
int dsz_inst,
5758
int*& nodeindices,
5859
double*& data,
59-
int*& pdata);
60+
int*& pdata,
61+
std::vector<int>& pointer2type);
6062

6163
extern int (*nrn2core_get_dat2_3_)(int tid,
6264
int nweight,

‎coreneuron/io/nrn_checkpoint.cpp‎

Lines changed: 74 additions & 29 deletions
Original file line numberDiff line numberDiff line change
@@ -88,6 +88,33 @@ void CheckPoints::write_checkpoint(NrnThread* nt, int nb_threads) const {
8888
#endif
8989
}
9090

91+
// Factor out the body of ion handling below as the same code
92+
// handles POINTER
93+
static int nrn_original_aos_index(int etype, int ix, NrnThread& nt, int** ml_pinv) {
94+
// Determine ei_instance and ei from etype and ix.
95+
// Deal with existing permutation and SoA.
96+
Memb_list* eml = nt._ml_list[etype];
97+
int ecnt = eml->nodecount;
98+
int esz = corenrn.get_prop_param_size()[etype];
99+
int elayout = corenrn.get_mech_data_layout()[etype];
100+
// current index into eml->data is a function
101+
// of elayout, eml._permute, ei_instance, ei, and
102+
// eml padding.
103+
int p = ix - (eml->data - nt._data);
104+
assert(p >= 0 && p < eml->_nodecount_padded * esz);
105+
int ei_instance, ei;
106+
nrn_inverse_i_layout(p, ei_instance, ecnt, ei, esz, elayout);
107+
if (elayout == Layout::SoA) {
108+
if (eml->_permute) {
109+
if (!ml_pinv[etype]) {
110+
ml_pinv[etype] = inverse_permute(eml->_permute, eml->nodecount);
111+
}
112+
ei_instance = ml_pinv[etype][ei_instance];
113+
}
114+
}
115+
return ei_instance * esz + ei;
116+
}
117+
91118
void CheckPoints::write_phase2(NrnThread& nt) const {
92119
FileHandler fh;
93120

@@ -186,7 +213,7 @@ void CheckPoints::write_phase2(NrnThread& nt) const {
186213
}
187214

188215
auto& memb_func = corenrn.get_memb_funcs();
189-
// will need the ml_pinv inverse permutation of ml._permute for ions
216+
// will need the ml_pinv inverse permutation of ml._permute for ions and POINTER
190217
int** ml_pinv = (int**) ecalloc(memb_func.size(), sizeof(int*));
191218

192219
for (NrnThreadMembList* current_tml = nt.tml; current_tml; current_tml = current_tml->next) {
@@ -224,9 +251,10 @@ void CheckPoints::write_phase2(NrnThread& nt) const {
224251

225252
sz = nrn_prop_dparam_size_[type];
226253
if (sz) {
227-
int* d = soa2aos(ml->pdata, cnt, sz, layout, ml->_permute);
228254
// need to update some values according to Datum semantics.
229-
if (!nrn_is_artificial_[type])
255+
int* d = soa2aos(ml->pdata, cnt, sz, layout, ml->_permute);
256+
std::vector<int> pointer2type; // voltage or mechanism type (starts empty)
257+
if (!nrn_is_artificial_[type]) {
230258
for (int i_instance = 0; i_instance < cnt; ++i_instance) {
231259
for (int i = 0; i < sz; ++i) {
232260
int ix = i_instance * sz + i;
@@ -238,32 +266,33 @@ void CheckPoints::write_phase2(NrnThread& nt) const {
238266
int p = pinv_nt[d[ix] - (nt._actual_diam - nt._data)];
239267

240268
d[ix] = p; // relative to _actual_diam
241-
} else if (s == -5) { // Assume pointer to membrane voltage
242-
int p = pinv_nt[d[ix] - (nt._actual_v - nt._data)];
243-
d[ix] = p; // relative to _actual_v
244-
} else if (s >= 0 && s < 1000) { // ion
245-
// determine ei_instance and ei
246-
int etype = s;
247-
Memb_list* eml = nt._ml_list[etype];
248-
int ecnt = eml->nodecount;
249-
int esz = nrn_prop_param_size_[etype];
250-
int elayout = corenrn.get_mech_data_layout()[etype];
251-
// current index into eml->data is a function
252-
// of elayout, eml._permute, ei_instance, ei, and
253-
// eml padding.
254-
int p = d[ix] - (eml->data - nt._data);
255-
int ei_instance, ei;
256-
nrn_inverse_i_layout(p, ei_instance, ecnt, ei, esz, elayout);
257-
if (elayout == Layout::SoA) {
258-
if (eml->_permute) {
259-
if (!ml_pinv[etype]) {
260-
ml_pinv[etype] = inverse_permute(eml->_permute,
261-
eml->nodecount);
262-
}
263-
ei_instance = ml_pinv[etype][ei_instance];
264-
}
269+
} else if (s == -5) { // POINTER
270+
// loop over instances, then sz, means that we
271+
// visit consistent with natural order of
272+
// pointer2type
273+
274+
// Relevant code that this has to invert
275+
// is permute/node_permute.cpp :: update_pdata_values with
276+
// respect to permutation, and
277+
// io/phase2.cpp :: Phase2::pdata_relocation
278+
// with respect to that AoS -> SoA
279+
280+
// Step 1: what mechanism is d[ix] pointing to
281+
int ptype = type_of_ntdata(nt, d[ix], i_instance == 0);
282+
pointer2type.push_back(ptype);
283+
284+
// Step 2: replace d[ix] with AoS index relative to type
285+
if (ptype == voltage) {
286+
int p = pinv_nt[d[ix] - (nt._actual_v - nt._data)];
287+
d[ix] = p; // relative to _actual_v
288+
} else {
289+
// Since we know ptype, the situation is
290+
// identical to ion below. (which was factored
291+
// out into the following function.
292+
d[ix] = nrn_original_aos_index(ptype, d[ix], nt, ml_pinv);
265293
}
266-
d[ix] = ei_instance * esz + ei;
294+
} else if (s >= 0 && s < 1000) { // ion
295+
d[ix] = nrn_original_aos_index(s, d[ix], nt, ml_pinv);
267296
}
268297
#if CHKPNTDEBUG
269298
if (s != -8) { // WATCH values change
@@ -273,8 +302,14 @@ void CheckPoints::write_phase2(NrnThread& nt) const {
273302
#endif
274303
}
275304
}
305+
}
276306
fh.write_array<int>(d, cnt * sz);
277307
delete[] d;
308+
size_t s = pointer2type.size();
309+
fh << s << " npointer\n";
310+
if (s) {
311+
fh.write_array<int>(pointer2type.data(), s);
312+
}
278313
}
279314
}
280315

@@ -298,7 +333,17 @@ void CheckPoints::write_phase2(NrnThread& nt) const {
298333
} else {
299334
Point_process* pnt = ps->pntsrc_;
300335
assert(pnt);
301-
output_vindex[i] = -(pnt->_i_instance * 1000 + pnt->_type);
336+
int type = pnt->_type;
337+
int ix = pnt->_i_instance;
338+
if (nt._ml_list[type]->_permute) {
339+
// pnt->_i_instance is the permuted index into pnt->_type
340+
if (!ml_pinv[type]) {
341+
Memb_list* ml = nt._ml_list[type];
342+
ml_pinv[type] = inverse_permute(ml->_permute, ml->nodecount);
343+
}
344+
ix = ml_pinv[type][ix];
345+
}
346+
output_vindex[i] = -(ix * 1000 + type);
302347
}
303348
}
304349
fh.write_array<int>(output_vindex, nt.n_presyn);

‎coreneuron/io/phase2.cpp‎

Lines changed: 78 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -48,7 +48,8 @@ int (*nrn2core_get_dat2_mech_)(int tid,
4848
int dsz_inst,
4949
int*& nodeindices,
5050
double*& data,
51-
int*& pdata);
51+
int*& pdata,
52+
std::vector<int>& pointer2type);
5253

5354
int (*nrn2core_get_dat2_3_)(int tid,
5455
int nweight,
@@ -172,7 +173,14 @@ void Phase2::read_file(FileHandler& F, const NrnThread& nt) {
172173
if (dsz > 0) {
173174
pdata = F.read_vector<int>(dsz * n);
174175
}
175-
tmls.emplace_back(TML{nodeindices, pdata, 0, {}, {}});
176+
tmls.emplace_back(TML{nodeindices, pdata, mech_types[i], {}, {}});
177+
if (dsz > 0) {
178+
int sz = F.read_int();
179+
if (sz) {
180+
auto& p2t = tmls.back().pointer2type;
181+
p2t = F.read_vector<int>(sz);
182+
}
183+
}
176184
}
177185
output_vindex = F.read_vector<int>(nt.n_presyn);
178186
output_threshold = F.read_vector<double>(n_real_output);
@@ -314,8 +322,13 @@ void Phase2::read_direct(int thread_id, const NrnThread& nt) {
314322
int* nodeindices_ = nullptr;
315323
double* data_ = _data + offset;
316324
int* pdata_ = const_cast<int*>(tml.pdata.data());
317-
(*nrn2core_get_dat2_mech_)(
318-
thread_id, i, dparam_sizes[type] > 0 ? dsz_inst : 0, nodeindices_, data_, pdata_);
325+
(*nrn2core_get_dat2_mech_)(thread_id,
326+
i,
327+
dparam_sizes[type] > 0 ? dsz_inst : 0,
328+
nodeindices_,
329+
data_,
330+
pdata_,
331+
tml.pointer2type);
319332
if (dparam_sizes[type] > 0)
320333
dsz_inst++;
321334
offset += nrn_soa_padded_size(nodecounts[i], layout) * param_sizes[type];
@@ -586,6 +599,16 @@ void Phase2::pdata_relocation(const NrnThread& nt, const std::vector<Memb_func>&
586599
// -9 (diam), // or 0-999 (ion variables).
587600
// Note that pdata has a layout and the // type block in nt.data into which
588601
// it indexes, has a layout.
602+
603+
// For faster search of tmls[i].type == type, use a map.
604+
// (perhaps would be better to replace tmls so that we can use tmls[type].
605+
std::map<int, size_t> type2itml;
606+
for (size_t i = 0; i < tmls.size(); ++i) {
607+
if (tmls[i].pointer2type.size()) {
608+
type2itml[tmls[i].type] = i;
609+
}
610+
}
611+
589612
for (auto tml = nt.tml; tml; tml = tml->next) {
590613
int type = tml->index;
591614
int layout = corenrn.get_mech_data_layout()[type];
@@ -615,8 +638,23 @@ void Phase2::pdata_relocation(const NrnThread& nt, const std::vector<Memb_func>&
615638
nt._actual_diam - nt._data, cnt, pdata, i, szdp, layout, nt.end);
616639
break;
617640
case -5: // pointer assumes a pointer to membrane voltage
641+
// or mechanism data in this thread. The value of the
642+
// pointer on the NEURON side was analyzed by
643+
// nrn_dblpntr2nrncore which returned the
644+
// mechanism index and type. At this moment the index
645+
// is in pdata and the type is in tmls[type].pointer2type.
646+
// However the latter order is according to the nested
647+
// iteration for nodecount { for szdp {}}
648+
// Also the nodecount POINTER instances of mechanism
649+
// might possibly point to differnt range variables.
650+
// Therefore it is not possible to use transform_int_data
651+
// and the transform must be done one at a time.
652+
// So we do nothing here and separately iterate
653+
// after this loop instead of the former voltage only
654+
/**
618655
transform_int_data(
619656
nt._actual_v - nt._data, cnt, pdata, i, szdp, layout, nt.end);
657+
**/
620658
break;
621659
default:
622660
if (s >= 0 && s < 1000) { // ion
@@ -637,6 +675,36 @@ void Phase2::pdata_relocation(const NrnThread& nt, const std::vector<Memb_func>&
637675
}
638676
}
639677
}
678+
// Handle case -5 POINTER transformation (see comment above)
679+
auto search = type2itml.find(type);
680+
if (search != type2itml.end()) {
681+
auto& ptypes = tmls[type2itml[type]].pointer2type;
682+
assert(ptypes.size());
683+
size_t iptype = 0;
684+
for (int iml = 0; iml < cnt; ++iml) {
685+
for (int i = 0; i < szdp; ++i) {
686+
int s = semantics[i];
687+
if (semantics[i] == -5) { // POINTER
688+
int* pd = pdata + nrn_i_layout(iml, cnt, i, szdp, layout);
689+
int ix = *pd; // relative to elem0
690+
int ptype = ptypes[iptype++];
691+
if (ptype == voltage) {
692+
nrn_assert((ix >= 0) && (ix < nt.end));
693+
int elem0 = nt._actual_v - nt._data;
694+
*pd = elem0 + ix;
695+
} else {
696+
Memb_list* pml = nt._ml_list[ptype];
697+
int pcnt = pml->nodecount;
698+
int psz = corenrn.get_prop_param_size()[ptype];
699+
nrn_assert((ix >= 0) && (ix < pcnt * psz));
700+
int elem0 = pml->data - nt._data;
701+
*pd = elem0 + nrn_param_layout(ix, ptype, pml);
702+
}
703+
}
704+
}
705+
}
706+
ptypes.clear();
707+
}
640708
}
641709
}
642710
}
@@ -1079,9 +1147,15 @@ void Phase2::populate(NrnThread& nt, const UserParams& userParams) {
10791147
#endif
10801148

10811149
// specify the ml->_permute and sort the nodeindices
1150+
// Have to calculate all the permute before updating pdata in case
1151+
// POINTER to data of other mechanisms exist.
10821152
for (auto tml = nt.tml; tml; tml = tml->next) {
10831153
if (tml->ml->nodeindices) { // not artificial
10841154
permute_nodeindices(tml->ml, p);
1155+
}
1156+
}
1157+
for (auto tml = nt.tml; tml; tml = tml->next) {
1158+
if (tml->ml->nodeindices) { // not artificial
10851159
permute_ml(tml->ml, tml->index, nt);
10861160
}
10871161
}

‎coreneuron/io/phase2.hpp‎

Lines changed: 1 addition & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -113,6 +113,7 @@ class Phase2 {
113113
int type;
114114
std::vector<int> iArray;
115115
std::vector<double> dArray;
116+
std::vector<int> pointer2type;
116117
};
117118
std::vector<TML> tmls;
118119
std::vector<int> output_vindex;

0 commit comments

Comments
 (0)