Using cppyy for julia/C++ interoperability

@jlperla: I think that simply ‘avoiding python’ is defeatist. Look at it as an opportunity instead. To get back to my b-spline library, the python wrapper I did with cppyy adds a fair amount of functionality, namely the pythonization layer. With this layer (and cppyy in the back) I can, at run time, pass template arguments to the construction of objects, and cppyy instantiates templates it hasn’t yet ready and the JIT creates code which is as fast as precompiled binary. The python layer is genuinely useful - if it weren’t there I had to write this part in julia to produce the same functionality, which would require an ‘intimate’ wrap with deep ties to the C++ code.

The different layers are quite well decoupled: the (julia) caller can set up arrays and parameters and then calls into the python layer. The python layer inspects the parametrization and builds code if necessary, then it delegates to the C++ layer which eventually does the number-crunching. It seems like a good way to split up the work, and the wrap of the already functional python module in julia can take over the python functionality pretty much 1:1, because in julia you wouldn’t really want a different interface.

So for the concrete example, and to answer @jstrube’s post, I’ll simply go ahead and see if it works - beyond my initial proof of concept. All the yes-butting might be silenced by a pudding proof. Succeeding, this is certainly much more exciting than just yet another b-spline module. If using cppyy doesn’t work for every C++/python module, then the code and/or tools can be worked over until they do.

Note that you will need my bspline module and your python3 must have cppyy installed for this to work. Here we go:

julia> using PyCall

julia> bspline = pyimport("bspline")
PyObject <module 'bspline' from '/home/kfj/.local/lib/python3.8/site-packages/bspline/__init__.py'>

julia> src=bspline.make_array((1000,1000),"double",3)
PyObject <bspline.array object at 0x7f9f66d3c460>

julia> trg=bspline.make_array((1000,1000),"double",3)
PyObject <bspline.array object at 0x7f9f66d3fb50>

julia> bspl=bspline.bspline("double", 3, 2, [1000,1000])
PyObject <bspline.bspline object at 0x7f9f66d3fbb0>

julia> set!(src, (5,5), (1.0,2.0,3.0))
(1.0, 2.0, 3.0)

julia> bspl.prefilter(src)

julia> bspline.index_transform(bspl,trg)

julia> get(trg,(5,5))
(0.9999999999999997, 1.9999999999999993, 2.9999999999999996)

This is a demonstration that the higher-level functions in the bspline module can be accessed from julia. I’ve resorted to creating the source and target arrays with the bspline module for this example - the arrays are normal numpy arrays (you can get this ‘aspect’ by accessing src.as_np, for example). If one were to ‘julianize’ the module, suitable array-to-array conversion would be nice-to-have, especially when it comes to multi-channel data.

So I create two arrays src and trg, each with a million tuples of three doubles. To make sure that the spline ‘works’ I assign (1,2,3) to the tuple at coordinates (5,5), which would be (6,6) in julia. I also create a bspline object, passing the number of channels (3), the dimensionality (2) and the shape (1000,1000). Note that I might have passed more arguments to, e.g. to get a spline of different degree (up to 45). The magic here is that the data type, number of channels and the dimensionality, which are template arguments in C++ (and therefore have to be known at compile-time) are run-time arguments to the python module, and are used to create suitable code on-the-fly if necessary.

Next I call ‘prefilter’, passing the ‘src’ array. This takes the values in src as knot point values for the spline, performs the prefiltering, and deposits the spline coefficients inside the b-spline object bspl, where they are held during the object’s scope, and can be accessed by referring to, e.g., bspl.core.as_np, which exposes them as a numpy array.

Now I call index_transform, passing the spline as a functor (the spline’s built-in default evaluator will be used) and trg, the target array. index_transform passes discrete indices into the target array to the functor passed as it’s first argument and writes the result to the target at that very index. This function is not a common feature, but it’s an ideal strting point for stencil code, where you fill a target array by means of a functor which is invoked with the target coordinates and takes it from there (an inverse transformation).

Finally, I expect that the index_transform should have deposited the evaluations of the spline in the target array, so I expect to find the tuple (1,2,3) at trg[6.6], which is indeed the case, minus a bit of error which has accrued along the way. I refrain from showing the result for the other 999999 tuples, they should all be near-zero :wink:

Beyond the mere proof of viability, this little demonstration should reveal that the process is quite fast - after all the arrays involved are reasonably large (you may have to repeat the process to see the true speed). The python layer has added benefits like accessing the module’s documentation, try:

py"help"(bspline)

Which is not what you would get “for free” when wrapping the C++ code - I hope this shows that having the python layer is indeed an added benefit. What you will notice is that the array handling is a bit clumsy, but this is due to my inadequacy using julia - NumPy arrays and julia arrays are compatible at a deep level (via python’s array interface), the trouble is to get the parametrization right, hence my ‘chickening out’ and using the bspline module to create suitable arrays. As these arrays are transparently accessible by julia, this is not really a problem either.

To conclude this lengthy post, and as an answer to @jstrube, I think that adding a little julia wrapper to ‘grok’ the python module would be much less work than producing an equivalent by wrapping the C++ code. How about giving that a shot? I’d also like to invite mkitti (I’d like to invite you via an at reference, but julialang.org forbids me to mention more than two other users because I’m new) , who has kindly offered to help wrapping my code for julia to consider the path via python - it may be a bit circuitous and pull in a lot of collateral code, but the result is extremely versatile and fast, as I hope to have demonstrated above - and it’s almost there. Wrapping the C++ code would - as it stands - not give us run-time template instantiation, all we could hope for is to wrap already instantiated templates with CxxWrap with additional intrusion into the C++ code - or waiting for Cxx to maybe become available. And with a julia wrap of the python module, a rich body of code would be immediately available to all other julia code requiring cardinal b-splines. I don’t know of another package with this scope.

1 Like