BlueBrain
Bugs/Miscompilations
Hi,
I am stresstesting Arbor's NMODL compiler and ran into some cases that I think are mishandled by nmodl as well:
KINETIC under conditionals
NEURON { SUFFIX test_kinetic_alternating_reaction }
STATE { A B }
BREAKPOINT {
SOLVE foobar METHOD sparse
}
KINETIC foobar {
LOCAL x, y
x = 23*v
y = 42*v
if (v<0) {
~ A <-> B (x, y)
}
else {
~ A <-> B (y, x)
}
}
Arbor
Throws an error; which I think should be the expected behaviour
nmodl
After 5/local rename
NEURON {
SUFFIX test_kinetic_conditional
}
STATE {
A
B
}
BREAKPOINT {
SOLVE foobar METHOD sparse
}
DERIVATIVE foobar {
LOCAL x, y
x = 23*v
y = 42*v
IF (v<0) {
} ELSE {
}
A' = (-1*(x*A-y*B))+(-1*(y*A-x*B))
B' = (1*(x*A-y*B))+(1*(y*A-x*B))
}
which is likely not what a user expects.
Alternating between reaction and non-reaction statements
NEURON { SUFFIX test_kinetic_alternating_reaction }
STATE { A B C D }
BREAKPOINT { : NOTE we need a newline here :/
SOLVE foobar METHOD sparse
}
KINETIC foobar {
LOCAL x, y
x = 23*v
y = 42*v
~ A <-> B (x, y)
x = sin(y)
y = cos(x)
~ C <-> D (x, y)
}
Arbor
Produces the same code as nmodl does, see below. However, I think this should constitute an error and Arbor will adopt this in the next version.
nmodl
After 5/local rename
NEURON {
SUFFIX test_kinetic_alternating_reaction
}
STATE {
A
B
C
D
}
BREAKPOINT {
SOLVE foobar METHOD sparse
}
DERIVATIVE foobar {
LOCAL x, y
x = 23*v
y = 42*v
x = sin(y)
y = cos(x)
A' = (-1*(x*A-y*B))
B' = (1*(x*A-y*B))
C' = (-1*(x*C-y*D))
D' = (1*(x*C-y*D))
}
Again, I think this is unexpected and likely constitutes a user error.
thanks, @thorstenhater!
The #89 mentions conditional statements in the title but there is no specific example there. These are bugs and I vaguely remember discussion around these:
- For 1, if NEURON does handle such MOD file correctly, we were wondering if that should be supported in NMODL as well
- For 2, we were also saying it would be easy for user to rewrite mod file. Or, whether NMODL should do something smart.
Just tagging here @nrnhines to confirm if nocmodl from NEURON is handling above examples.
For 1. the definition of correct should be this, I assume
KINETIC foobar {
LOCAL x1, y1, x2, y2
x1 = 23*v
y1 = 42*v
x2 = sin(y1)
y2 = cos(x1)
~ A <-> B (x1, y1)
~ C <-> D (x2, y2)
}
which is readily doable.
For 2.) I can imagine a similar transformation, but is much more involved, as one would need to hoist sideeffects across the conditional. However, things might become infeasible ad/or confusing to the user. Essentially -- I assume, since I see no better way -- one would compute all rates for all combinations of the conditionals and set them before producing the ODE system. That however is subject to combinatorical explosion, $2^N$ with $N$ conditionals.
Essentially this
KINETIC foobar {
LOCAL x, y
if (a()) {
x = f()
y = g()
~ A <-> B (x, y)
} else {
x = h()
y = i()
~ 2B <-> C (y, x)
~ C <-> D (x, y)
}
}
should become
KINETIC foobar {
LOCAL rAB, iAB, r2BC, i2BC, rCD, iCD
if (a()) {
rAB = f()
iAB = g()
r2BC = 0
i2BC = 0
rCD = 0
iCD = 0
} else {
rAB = 0
iAB = 0
r2BC = i()
i2BC = h()
rCD = h()
iCD = i()
}
~ A <-> B (rAB, iAB)
~ 2B <-> C (r2BC, i2BC)
~ C <-> D (rCD, iCD)
}
which already gets taxing.
Can we document the definitions of the correct outputs somewhere?
The kinetic schemes shown with conditional and alternating reaction/non-reaction work properly with NEURON's nocmodl (and are stable as long as the rates are positive). I do believe, though, that an argument can be made that conditionals are a bad idea as there are implicit requirements in regard to number of states being the same. To get sensible results with simple tests. I modified the examples to be.
$ cat test1.mod
NEURON { SUFFIX test_kinetic_conditional }
STATE { A B }
INITIAL {
A = 1
B = 0
}
BREAKPOINT {
SOLVE foobar METHOD sparse
}
KINETIC foobar {
LOCAL x, y
x = 23*fabs(v)
y = 42*fabs(v)
if (v<0) {
~ A <-> B (x, y)
}
else {
~ A <-> B (y, x)
}
}
and
$ cat test2.mod
NEURON { SUFFIX test_kinetic_alternating_reaction }
STATE { A B C D }
INITIAL {
A = 1
B = 0
C = 1
D = 0
}
BREAKPOINT { : NOTE we need a newline here :/
SOLVE foobar METHOD sparse
}
KINETIC foobar {
LOCAL x, y, z
x = 1*fabs(v)
y = 1
~ A <-> B (x, y)
z = fabs(sin(y))
y = fabs(cos(x))
x = z
~ C <-> D (x, y)
}
Test1 with a square wave voltage gives
and Test2 with a step function voltage gives
Alternating Reaction Equations & Assignment
NMODL has changed how it handles KINETIC blocks. All rates are converted to LOCAL variables. Therefore, the above example turns into the following:
KINETIC foobar {
LOCAL x, y, kf0_, kb0_, kf1_, kb1_
x = 23*v
y = 42*v
kf0_ = x
kb0_ = y
~ A <-> B (kf0_, kb0_)
x = sin(y)
y = cos(x)
kf1_ = x
kb1_ = y
~ A <-> B (kf1_, kb1_)
}
which you'll note that it can be safely reordered:
KINETIC foobar {
LOCAL x, y, kf0_, kb0_, kf1_, kb1_
x = 23*v
y = 42*v
kf0_ = x
kb0_ = y
x = sin(y)
y = cos(x)
kf1_ = x
kb1_ = y
~ A <-> B (kf0_, kb0_)
~ A <-> B (kf1_, kb1_)
}
and from there converted to an ODE.
This was implemented to make MOD files like work: https://modeldb.science/98005?tab=2&file=D2modulation/ca_ch.mod which contains:
FROM i=0 TO NANN-1 {
dsqvol = dsq*vol[i]
~ ca[i] + Buffer[i] <-> CaBuffer[i] (k1buf*dsqvol,k2buf*dsqvol)
}
which after loop unrolling turns into the exact issue reported.
It seems reasonably clear that the intended behaviour is:
~ ca[i] + Buffer[i] <-> CaBuffer[i] (k1buf*dsq*vol[i],k2buf*dsq*vol[i])
and expressing it via dsqvol is an "optimization". Similarly, with the reported example, there's a sequence in which statements appear. The assignment of x and y happens at particular points in that sequence, before the change the old value must be used and after the assignment the new value must be used, until it's changed again.
