4 ms·
> With DifferentialEquations.jl I also could just yep... I took a look at the DE packages in Julia today, and quite frankly they're much better than the situat
by marmaduke 8y ago
> With DifferentialEquations.jl I also could just
yep... I took a look at the DE packages in Julia today, and quite frankly they're much better than the situation in Python, perhaps because of one or more prolific applied mathematicians are making a concerted effort, which is lacking Python? I dunno, but I did recommend my colleagues look at Julia for DEs, for this reason.
That said,
> Pass a list of functions, and some parameters and I construct a right hand side function from it. Unfortunately this does not work with numba.
I'm pretty sure I've done this before with numba, so maybe getting concrete would help, e.g. an Euler step
def euler(f, dt, nopython=False):
@numba.jit(nopython=nopython)
def step(x):
return x + dt * f(x)
where user can provide regular Python function or a @numba.jit'd function. If a @numba.jit'd function is provided, and nopython=True, this should result in fused machine code. This sort of code gen through nest functions can be done repeatedly for e.g. the time stepping loop.
I've done this for CPU & GPU code for a complex model space (10 neural mass models, 8 coupling functions, 2 integration schemes, N output functions, ...) which, by the above pattern, results in flexible but fast code.
Is this a pattern that captures your use case or not yet?
> implement sparse matrix algorithms by hand because I can't call scipy.sparse.
agreed, this is a surprising omission, which I attribute to not much of the numerical Python community making use of Numba, but could be fixed rapidly.
> constrained my users to use numba compatible code in their right hand side functions
what did you run into that was problematic?
> I had to implement a non-linear solver from scratch in numba rather than being able to use scipys excellent set of solvers
I didn't follow; passing @numba.jit'd functions to scipy is in the Numba guide, so what exactly didn't work?
- Certhas 8y agoThis pattern is how I wrote the SDE solver in Python. That works great and is really useful and the reason why I teach closures. The library we're building now though does something different. Something like this: def network_rhs(fs, Network) def rhs(y,t) y2 = np.dot(Network, y) r = empty_like(y) for i, f in enumerate(fs): r[i] = f(y2[i]) return r return rhs > what did you run into that was problematic? For more complex model building the right hand side functions actually make use of fairly complex class hierarchies. That was the major stumbling block. But people also were using dictionaries and other non-numpy data structures and just generally idiomatic Python that is not always supported. Some of that stuff is inherently slow/bad design of course, but it still ended up killing the use of my solver for this project. They are now rewriting in C++, which is absolutely a great choice for their case (and probably would have been viable for us too if we had had more people with a C/C++ background in the team). > passing @numba.jit'd functions to scipy is in the Numba guide I wanted to use scipy.root from numba. Not the other way around. Now if all of the numerical Python community was standardized on numba, a lot of this would not be an issue. Scipys LowLevelCallable is a great step in the right direction. But fundamentally I don't see how you will ever get the different libraries to play together nicely in a performant way. It would require every API to expose numba jitable functions. Last I checked, the only functions you could call from within numba code were other numba functions and the handful of numpy features the numba authors implemented themselves (I remember waiting for dot and inv support). If I have an algorithm by a student implemented on a networkx graph as a data structure I can't just jit that. In Julia it automatically is.
- marmaduke 8y agoI see what you mean. I’ve done exactly that sort of thing in C with an array of function pointers, but I’m not sure it would work in Numba. The churn is exhausting but I see the merit of starting over and getting everything done in a fully fledged JITd language.
- Certhas 8y agoAs said in another part of the thread, we were in a very good spot to do so. There absolutely are many reasons to not jump in at this time.