How to choose initial state vectors for a resonant orbit numerical simulation? How do higher numbers affect orbital resonance? is a really interesting question, and I thought I would try to do a simple numerical simulation of various resonance ratios to see if some are inherently stable or inherently unstable.
I understand that numerical simulation is inexact and so it can introduce its own instabilities so attention will have to be paid to the choice of integrators and results taken with a large grain of salt.
In order to get something to happen within say a hundred or a thousand orbits the masses of the orbiting bodies have to be large enough to have a significant impact on each other1, and so the state vectors I choose to initialize the calculation will not necessarily be those of circular orbits each body neglecting the other.
Are there any guidelines or rules of thumbs for initializing resonant orbits?

1I'll assume a central force, i.e. the orbiting bodies will not affect the position of the body they orbit, and so there won't be any indirect coupling through it.
Related:

Just how &quot;locked&quot; are resonant-chains of exoplanets thought to be? (e.g. K2-138 and TOI-178) (currently unanswered)
Why does the exoplanetary system TOI-178 challenges current theories of planet formation?
I don't know if there are rules of thumb for initializing resonant orbital simulations.  I might suggest:

Setting initial state vectors' Z-Components to zero.  Keeping the whole simulation in the X-Y plane can help with debugging and visualization.

Choosing initial state vectors that are equivalent to orbits with non-zero eccentricity. Orbital resonances tend to amplify or dampen eccentricities, but the effects can be multiplicative. So if you choose circular orbits, you may see less of an effect.

Choosing initial state vectors that line up apogees and perigees with closest contacts.

Running multiple simulations with slight perturbations to initial state vectors to test stability.

Some interesting initial experiments might be to see if you could replicate the instability of the L1 and L2 points of an orbit.  Or perhaps put some little masses in the Kirkwood Gaps and see if they get expelled.


Notes:

I tried writing some of my own numerical simulations to simulate natural planetary flybys. I used an algorithm outlined by Voesenek, but I was getting round-off errors that shot my simulation. Instead of investing more time to fix my simulation, I ended up using Universe Sandbox 2, that I paid for and downloaded from Steam.

Universe Sandbox probably doesn't have the fidelity you need for resonance stability simulations.  Gallardo wrote a pretty slick long term orbital integration software package called ORBE, that is free to the public for educational purposes.  I wish I had read Tabaré Gallardo's 2017 Exploring the orbital evolution of planetary systems before writing my own simulator, since he uses a neat trick by Encke to account for the Sun's gravity analytically via 2-body problem equations while separately treating the inter-planetary gravitational accelerations numerically. This allows larger time steps and prevents some kinds of round-off error.  Gallardo's software is designed to test long term stellar system stability.  Also, it uses orbital elements rather than state vectors, which are a bit less intuitive to me. Gallardo's used of ORBE to run Kirkwood Gap stability simulations is totally AWESOME!!  Here is a capture of one of the figures from Gallardo's 2017: