In BifurcationKit, I tested a lot the shooting paradigm. My “general” conclusion is
- you need multiple shooting (unless the orbit is super attracting)
- you need an implicit solver with relatively high tolerances, I often use
Rodas5P
You can see that in the tutorials. This is especially true when you do something fancy like computing curves of bifurcations of periodic orbits.
It would be nice to try orthogonal collocation on this, I will see if I have some time.