I suggest the following default strategy for calculating
The idea is to build up two partial Schreier trees Tx, Ty (w.r.t.
G.generators), rooted at x and y, respectively. One alternately
adds nodes to Tx and Ty until they have a node z, say, in common (if
one of the trees is completed and there is no such node then return
false). Given such a z, one uses Tx,Ty in the usual way to
get elements gx,gy such that z=action(x,gx)=action(y,gy),
and so the required element taking x to y is gx/gy.
If x and y are in the same G-orbit, then this approach should be
much faster than computing an orbit for x until y is found, and
much much faster than trying to compute a StabChain if the action is
Regards, Leonard Soicher.