What are the dependency requirements? Does it have to be 100% pure C++ of some version? Can you use libraries like Eigen or Boost, or is strict no-dependency C++?
If this has to be a 100% pure C++ project and you only need to be able to fit this one model and can put some time/compute into tuning the tuning parameters, then I’d suggest directly coding the multinomial version of HMC.* It’s a rejection-free algorithm that also solves the nasty harmonics problem of vanilla HMC. It’s not as efficient as NUTS, but it’s not that far off once you have a reasonable step size and max number of steps.
It it doesn’t have to be pure C++, then you can do what @jessegrabowski suggests and embed the C++ definition in a foreign function interface through ctypes. It’s maybe a bit easier than writing multinomial HMC but not much.
I don’t know any existing C++ implementations of something like NUTS that are detached from probabilistic programming frameworks and don’t have dependencies into things
* randomly split total number of steps L into F + B, where F are forward steps and B backward steps, then select a point proportional to its negative exponentiated energy.