Saturation limits/controlling output scipy's odeint using sympy/pydy

Question

I'm using PyDy to create the system of ODE's for three-body system as such:

right_hand_side = pydy.codegen.generate_ode_function(
       forcing_vector,
       model.coordinate_symbols,
       model.speed_symbols,
       model.constant_symbols,
       mass_matrix=mass_matrix)

Which is then integrated:

y = odeint(right_hand_side,
           x0_parameters.values(),
           t,
           args=(symbolic_constant_parameters,))

However, I need to constrain some of the values output by odeint; How can I go about this?

I understand that If I manually created the system of ODE's that I could create the saturation logic there, but I'm not entirely sure how to do that when I'm using the PyDy generate ODE function.

Not sure if its pertinent info, but the forcing vector and mass matrix come from setting up kinematic equations and feeding these into Sympy mechanic's KaneMethod.