Facilitating user-provided propensity calculations

[heejochoi]

A helpful feature for systems with alternative propensity calculations (such as when stochastic constants are unavailable) might be support for externally provided propensities. For example, each reaction:

X + Y -> Z

may have an associated propensity function described by the user that depends on the counts of molecules at each iteration:

def get_propensity(X, Y, c1, c2, ...):
    # computation here could involve several constants like c1, c2, ...
    ...

which returns a value that represents the propensity for this reaction, for this iteration. Since the inputs (X and Y) are dynamic during the simulation, the return value from get_propensity() also varies during the simulation.

Two candidate approaches to supporting this feature might be providing a way for the user:

1. to provide the propensity functions to arrow once at the initialization of the system - (get_propensity() is called inside arrow)
2. to input the propensity values at each iteration - (get_propensity() is called outside arrow by the user)

[jmason42]

This does seem like a desirable feature, even if these alternative propensity functions aren't as well founded as the typical mass-action-kinetic forms.

I think (1) would provide a richer compliment of features (i.e. I think (2) is easy to emulate with (1), but not vice-versa). Ideally it would be possible to call compiled functions from native C code rather than needing to pass through a Python layer, though short of writing a modeling sublanguage (like Theano and other tools) I'm not sure how this could be accomplished elegantly. It might be most practical to write the needed propensity function in C.

[jmason42]

Further thoughts: maybe the healthy middle ground would be to implement compiled 'metafunctions'. The current approach is like a metafunction (by which I mean a function that returns functions) in that we provide some information (which molecules participate in a reaction) and in turn it caches all the information needed to (re)evaluate that particular function efficiently (indices, dependencies, etc.).

[heejochoi]

Thanks, @jmason42! I agree:

It might be most practical to write the needed propensity function in C.

@prismofeverything and I got a chance to sketch the implementation of this route -- I'll start a branch for this feature.

[1fish2]

FYI: A “higher order function” is the term for a function that returns functions and/or takes functions as arguments.

[jmason42]

Yeah, I suppose a metafunction wouldn't be a function emitter so much as something that modifies how functions work. (Incidentally there is a concept called a 'metafunction' that has nothing to do with programming and has one of the most baffling Wikipedia entries :rofl: .)