Implementing something like tt.trapz() for numerical integration over a set domain?

It should work as the gradient of a integrate class is just the original function itself. So as long as it is some integrate approximation, you should be able to just replace the approximation.