I want to integrate a control problem, i.e. an ODE of the form dx/dt = A.f(t), where
x(t) is a function in R^3, f(t) is a function in R^4 and A is a matrix 3x4. In my special case, f(t) = F'(t), i.e. a time derivative of a function F. Furthermore, F is 1-periodic. Hence, integrating the ODE over the interval [0, 1] should yield the starting position again. However, methods like solve_ivp from scipy.integrate do not respect this periodicity at all (I have tried all the possible methods like RK45, Radau, DOP853, LSODA).
Is there a special ODE solver that respects such periodicity to a high degree of precision?
All algorithms will suffer from precision loss anyway, so you will most likely never achieve exact periodicity. Nonetheless, you can also try to increase the precision of the integration by using the parameters "atol" and "rtol", which, roughly, will keep the error of the integration below (atol + rtol*y) at each time step. You can typically go as low as atol=rtol=1e-14.