Link errors when generating special-core with kinetic scheme NMOD files

[pkuzyc]

Hello. I can successfully generate "special-core" using nrnivmodl-core with some general NMOD files. However, when there are some NMOD files implemented with kinetic scheme, it failed.

I used nrnivmodl-core with command:
./nrnivmodl-core /path/to/nmodfiles
and the errors are:

coreneuron/coreneuron.cpp:
@E@nvlink error   : Undefined reference to '_ZN10coreneuron13sparse_threadEPNS_9SparseObjEiPiS2_PddiiiiS3_S2_PNS_11ThreadDatumEPNS_9NrnThreadEd' in 'x86_64/libcorenrnmech.a:ca_ch.o'
@E@nvlink error   : Undefined reference to '_ZN10coreneuron17_ss_sparse_threadEPNS_9SparseObjEiPiS2_PddiiiiS3_S2_PNS_11ThreadDatumEPNS_9NrnThreadEd' in 'x86_64/libcorenrnmech.a:hcn12_ch.o'

It seems that it's about the function sparse_thread? Could you please take a look at it? Thank you.

[cattabiani]

Hello there!

could be useful to check also the modfile. Do you have a toy example where this appears?

[pkuzyc]

Here are the code of two example mod files:

ca_ch.mod:

:Migliore file Modify by Maciej Lazarewicz (mailto:mlazarew@gmu.edu) May/16/2001

TITLE Calcium ion accumulation and diffusion
: The internal coordinate system is set up in PROCEDURE coord_cadifus()
: and must be executed before computing the concentrations.
: The scale factors set up in this procedure do not have to be recomputed
: when diam1 or DFree are changed.
: The amount of calcium in an annulus is ca[i]*diam1^2*vol[i] with
: ca[0] being the second order correct concentration at the exact edge
: and ca[NANN-1] being the concentration at the exact center

NEURON {
	SUFFIX ca_ch
	USEION ca READ cao, cai, ica WRITE cai, ica
	RANGE ipump,last_ipump,test
	RANGE DFree, k1buf, k2buf, k1, k2, k3, k4, totpump, totbuf
	RANGE vol, Buffer0
}

DEFINE NANN  4

UNITS {
        (mol)   = (1)
	(molar) = (1/liter)
	(mM)	= (millimolar)
	(um)	= (micron)
	(mA)	= (milliamp)
	FARADAY = (faraday)	(10000 coulomb)
	PI	= (pi) 		(1)
}

PARAMETER {
	DFree	= .6	(um2/ms)
	diam 	= 1	(um)
	cao		(mM)
	ica		(mA/cm2)
	k1buf 	= 500	(/mM-ms)
	k2buf 	= 0.5	(/ms)
        k1	= 1.e10 (um3/s)	
        k2	= 50.e7 (/s)	: k1*50.e-3
        k3	= 1.e10 (/s)	: k1
        k4	= 5.e6	(um3/s)	: k1*5.e-4
	totpump	= 0.2	(mol/cm2)
	totbuf	= 1.2	(mM)
} 

CONSTANT { volo=1  (liter)}

ASSIGNED {
	area		(um2)
	test
	cai		(mM)
	vol[NANN]	(1)	: gets extra cm2 when multiplied by diam^2
	ipump           (mA/cm2)
	last_ipump	(mA/cm2)
}

STATE {
	ca[NANN]	(mM) <1.e-5> : ca[0] is equivalent to cai
	CaBuffer[NANN]  (mM)
	Buffer[NANN]    (mM)
        pump            (mol/cm2) <1.e-3>
        pumpca          (mol/cm2) <1.e-15>

}

BREAKPOINT {
	SOLVE state METHOD sparse
	:SOLVE state METHOD cnexp
	last_ipump=ipump
	ica = ipump
	test = 0
}

LOCAL coord_done

INITIAL {
	if (coord_done == 0) {
		coord_done = 1
		coord()
	}
	: note Buffer gets set to Buffer0 automatically
	: and CaBuffer gets set to 0 (Default value of CaBuffer0) as well
	FROM i=0 TO NANN-1 {
		ca[i] = cai
	}
 
       	ipump 	= 0
        pump 	= totpump
        pumpca 	= (1e-18)*pump*cao*k4/k3

        FROM i=0 TO NANN-1 {
               	ca[i] = cai
		CaBuffer[i] =(totbuf*ca[i])/(k2buf/k1buf+ca[i])
		Buffer[i] = totbuf - CaBuffer[i]	
	}
}

LOCAL frat[NANN] 	: gets extra cm when multiplied by diam

PROCEDURE coord() {
	LOCAL r, dr2
	: cylindrical coordinate system  with constant annuli thickness to
	: center of cell. Note however that the first annulus is half thickness
	: so that the concentration is second order correct spatially at
	: the membrane or exact edge of the cell.
	: note ca[0] is at edge of cell
	:      ca[NANN-1] is at center of cell
	r = 1/2					:starts at edge (half diam)
	dr2 = r/(NANN-1)/2			:half thickness of annulus
	vol[0] = 0
	frat[0] = 2*r
	FROM i=0 TO NANN-2 {
		vol[i] = vol[i] + PI*(r-dr2/2)*2*dr2	:interior half
		r = r - dr2
		frat[i+1] = 2*PI*r/(2*dr2)	:exterior edge of annulus
					: divided by distance between centers
		r = r - dr2
		vol[i+1] = PI*(r+dr2/2)*2*dr2	:outer half of annulus
	}
}

LOCAL dsq, dsqvol : can't define local variable in KINETIC block or use
		:  in COMPARTMENT
KINETIC state {
	COMPARTMENT i, diam*diam*vol[i]*1(um) {ca CaBuffer Buffer}
        COMPARTMENT (1.e10)*area {pump pumpca}
        COMPARTMENT (1.e15)*volo {cao}

	~ ca[0] << (-(ica-last_ipump)*PI*diam*frat[0]*1(um)/(2*FARADAY))

	FROM i=0 TO NANN-2 {
		~ ca[i] <-> ca[i+1] 	(DFree*frat[i+1]*1(um), DFree*frat[i+1]*1(um))
	}

	dsq = diam*diam*1(um)
	FROM i=0 TO NANN-1 {
		dsqvol = dsq*vol[i]
		~ ca[i] + Buffer[i] <-> CaBuffer[i] (k1buf*dsqvol,k2buf*dsqvol)
	}

        ~ca[0] + pump <-> pumpca 	((1.e-11)*k1*area, (1.e7)*k2*area)
        ~pumpca       <-> pump + cao 	((1.e7)*k3*area, (1.e-11)*k4*area)

        ipump = 2*FARADAY*(f_flux-b_flux)/area

	cai = ca[0]
}
	
COMMENT
At this time, conductances (and channel states and currents are
calculated at the midpoint of a dt interval.  Membrane potential and
concentrations are calculated at the edges of a dt interval.  With
secondorder=2 everything turns out to be second order correct.
ENDCOMMENT

hcn12_ch.mod:

COMMENT

Josh Held's adaptation to suit HCN1+2.  12/22/2003

****
Kinetic model of HCN2 channel gating from Wang et al 2002.

In this model channel opening is coupled to a change in the affinity of the cyclic nucleotide binding domain for cAMP which is manifest as a shift in the activation curve toward more positive potentials.  This model explains the slow activation kinetics of Ih associated with low concentrations of cAMP.

For further details email Matt Nolan at mfnolan@fido.cpmc.columbia.edu.

Reference

Wang J., Chen S., Nolan M.F. and Siegelbaum S.A. (2002). Activity-dependent regulation of HCN pacemaker channels by cyclicAMP: signalling through dynamic allosteric coupling. Neuron 36, 1-20.
****
ENDCOMMENT

NEURON {
	SUFFIX hcn12_ch
	NONSPECIFIC_CURRENT i
	RANGE i, ehcn, g, gbar
	RANGE a0, b0, ah, bh, ac, bc, aa0, ba0
	RANGE aa0, ba0, aah, bah, aac, bac
	RANGE kon, koff, b, bf, ai, gca, shift
}

UNITS {
	(mV) = (millivolt)
	(molar) = (1/liter)
	(mM) = (millimolar)
	(mA) = (milliamp)
	(S) = (siemens)
}

PARAMETER {
	gbar    = 1		(S/cm2)
	ehcn    = -20		(mV)
	a0      = .006		(/ms)	: parameters for alpha and beta
	b0      = .0008		(/ms)
	ah      = -96		(mV)
	bh      = -51.7		(mV)
	ac      = -.155		(/mV)
	bc      = .144		(/mV)
	aa0     = .0006		(/ms)	: parameters for alphaa and betaa
	ba0     = .004		(/ms)
	aah     = -94.2		(mV)
	bah     = -35.5		(mV)
	aac     = -.075		(/mV)
	bac     = .144		(/mV)
	kon     = 30		(/mM-ms) : cyclic AMP binding parameters
	koff    = 4.5e-05	(/ms)
	b       = 80
	bf      = 8.94
	ai	= 1e-05		(mM)	:concentration cyclic AMP
	gca     = 1			: relative conductance of the bound state
	shift   = -17		(mV)	: shift in voltage dependence
	q10v    = 4                     : q10 value from Magee 1998
	q10a    = 1.5			: estimated q10 for the cAMP binding reaction
	celsius			(degC)
}

ASSIGNED {
	v	(mV)
	g	(S/cm2)
	i	(mA/cm2)
	alpha	(/ms)
	beta    (/ms)
	alphaa	(/ms)
	betaa	(/ms)
}

STATE {
	c
	cac
	o
	cao
}

INITIAL {
	SOLVE kin STEADYSTATE sparse
	:SOLVE kin STEADYSTATE cnexp
}

BREAKPOINT {
	SOLVE kin METHOD sparse
	:SOLVE kin METHOD cnexp
	g = gbar*(o + cao*gca)
	i = g*(v-ehcn)
}

KINETIC kin {
	LOCAL qa
	qa = q10a^((celsius-22 (degC))/10 (degC))
	rates(v)
	~ c <-> o       (alpha, beta)
	~ c <-> cac     (kon*qa*ai/bf,koff*qa*b/bf)
	~ o <-> cao     (kon*qa*ai, koff*qa)
	~ cac <-> cao   (alphaa, betaa)
	CONSERVE c + cac + o + cao = 1
}

PROCEDURE rates(v(mV)) {
	LOCAL qv
	qv = q10v^((celsius-22 (degC))/10 (degC))
	if (v > -200) {
		alpha = a0*qv / (1 + exp(-(v-ah-shift)*ac))
		beta = b0*qv / (1 + exp(-(v-bh-shift)*bc))
		alphaa = aa0*qv / (1 + exp(-(v-aah-shift)*aac))
		betaa = ba0*qv / (1 + exp(-(v-bah-shift)*bac))
	} else {
		alpha = a0*qv / (1 + exp(-((-200)-ah-shift)*ac))
		beta = b0*qv / (1 + exp(-((-200)-bh-shift)*bc))
		alphaa = aa0*qv / (1 + exp(-((-200)-aah-shift)*aac))
		betaa = ba0*qv / (1 + exp(-((-200)-bah-shift)*bac))
	}
}

[cattabiani]

sorry but I forgot to mention: do you use mod2c or nmodl in coreneuron?

[pkuzyc]

I used nrnivmodl-core to generate the file "special-core", and the command I used is: ./nrnivmodl-core /path/to/mod

[cattabiani]

yes but, coreneuron (and thus nrnivmodl-core) was built with nmodl or mod2c activated?

[pkuzyc]

With mod2c activated, here are some outputs in compiling:

[cattabiani]

Thank you very much! We will investigate on this

[pramodk]

Sorry for delay @pkuzyc ! #307 should fix the link error you have reported.

One of the thing I would like to mention here is the the default MOD file to C++ converter we are using (https://github.com/BlueBrain/mod2c) doesn't generate optimal code for GPU/CPU especially for the KINETIC blocks.

New translator we are working on (https://github.com/BlueBrain/nmodl) has much better support but currently has one issue for GPU. We are still working on a fix.

I would be curious to know the performance you get. If you have any issues, please let us know.