Compatibality with C/C++ wrapped functions

I have a multidimensional probability distribution density of the following form (using Tex notation):
\rho \sim \exp(-(T + U) ), where T is, for example:
T = \frac{p_R^2}{2 \mu} + \frac{p_\theta^2}{2 \mu R^2} + \frac{p_\varphi^2}{2 \mu R^2 \sin^2 \theta}
(Hamiltonian of a particle in spherical coordinates).
So, in this example T is a 5-dimensional function of (R, p_R, \theta, p_\theta, p_\varphi) variables; \mu is fixed parameter. Potential function U is a function of variable R, but the problem is that it’s written in some compiled language (C/C++/Fortran). (I can wrap it to call from Python).
Does PyMC3 provide functionality to sample from described probability function? If so, could you provide a minimal example?
Kind thanks

You have to wrap the potential function U into a theano function, and then call it within with pm.Model():

There are quite a look of similar discussion here if you do a search of (C/C++/Fortran), eg:

You can also have a look at this discussion on Github:

As for sampling from the target potential function, you can write it as a pm.Potential('rho', -(T+U)).

Thank you for reply.
The link to github issue is unavailable…

oops, please find the updated link above