Currently, AFAICT, dual averaging is kind of hard-coded into NUTS, and the mass matrix adaptation simply depends on the type of the provided metric. Would be nice if there were a different way of specifying the used adaptation method, which "of course" should be orthogonal to the used metric.
Currently, this involved too much hacky boilerplate code for my tastes, see e.g. #473 (comment)