Function inside a kernel

[pierrick-giffard]

I noticed it was not possible to use functions inside kernels, even though it would be really useful in some cases. Does anyone has any idea to implement this?

[CKehl]

Hi. General functions inside a Kernel makes little sense, as the Kernel is specific to be called for each particle, and translated buy the dedicated Just-In-Time compiler, which only allows a subset of (mostly arithmetic) instructions. in #719 though, I introduced a parameter called 'postIterationFunctions' to ParticleSet.execute(), that takes a python array of (non-kernel) function handles, which are called after each kernel iteration and which can be of any sort and flavour you have in Python. If you like, please have a look on #719 , look on the 'Changes' tab, browse to the 'particleset.py' and see if that new functionality is what you want to have.

Regards, Christian

[erikvansebille]

I think what @pierrick-giffard meant is functions inside the Kernel itself to simplify the kernel code. For example, I think that @pierrick-giffard needs to calculate the spatial gradient of a field on a location (using e.g. central differencing), on every step of the RK4 scheme, and for many fields. Being able to create a function for this gradient-calculation would make the code much cleaner and easier to debug

By the way, I now remember that @delandmeterp has previously implemented a way to provide your own header files in C; see for example https://github.com/OceanParcels/parcels/blob/master/tests/test_kernel_language.py#L254-L288. Not for the faint-at-heart, but perhaps this could work for you too, @pierrick-giffard?

[willirath]

Let me quickly chime in to add the remark that as a far target, using Numba with ScipyParticles might be the way to go. :) . (Also sorry for the close / reopen.)

[pierrick-giffard]

@erikvansebille, for the moment I managed to implement a code not to dirty using for loops within the kernel.

[erikvansebille]

Let me quickly chime in to add the remark that as a far target, using Numba with ScipyParticles might be the way to go. :) . (Also sorry for the close / reopen.)

Yes, I agree @willirath that this is something that we should certainly explore. But for now on the far horizon!

for the moment I managed to implement a code not to dirty using for loops within the kernel.

Wow, I thought that for-loops were not supported by the Parcels code converter. Can you share a snippet of your code? I'd be keen to see how you did it!

[pierrick-giffard]

That's true, I hadn't tested it yet but it doesn't work with for loops...

[pierrick-giffard]

I think I can deal with it but I will get more troubles later because I'd like to use the random.vonmisesvariate function...

[erikvansebille]

Could you raise a separate Issue to ask for a random.vonmisesvariate functionality? And ideally also provide a c-implementation that we can use, see for example the normalvariate and expovariate implementation in https://github.com/OceanParcels/parcels/blob/master/parcels/include/random.h

[erikvansebille]

Closing for now, as no urgent need to implement functions within functions yet