A possible way to simplify (in steps) the potential surface

V = V(r,θ,φ), exact, useless
V = V(r) + V(θ,φ) + δV(r,θ,φ), separation + coupling, not very representative
V = V(θ,φ) + ½k[r-r₀(θ,φ)]², closest to reality, note that k, if treated as a variable of (θ,φ), has a variance less than 6%.

A possible approach to the third potential: V = ½k[r-r₀(θ,φ)]² + ε*V(θ,φ). The first part resembles a study of the Monodromy problem (of course it is much more complicated, because now we have 3 DOF, however locally we can approximate r₀(θ,φ)≈r₀(φ), for example, and from there we can further separate the phase space into tiles and continue the study.