Four leading theories for high-energy rearranging collisions of nuclei with atomic targets are critically reviewed. The focus is on the central aspects of charge exchange, ionisation and detachment in scatterings of fast protons with atomic hydrogen, helium and negative hydrogen ions. The methods used include the impulse approximation (IA), reformulated impulse approximation (RIA), exact boundary corrected second Born (CB2 or B2B) approximation and continuum distorted wave (CDW) approximation. Regarding three-body problems, the most comprehensive computations to date are carried out on both differential (dQ/dΩ) and total (Q) cross sections for proton–hydrogen charge exchange at impact energies ranging from 25 keV to 7.5 MeV. Throughout this energy range, for which the measured results exist, the quantitative agreement between the RIA and the available experimental data is found to be systematically excellent. Moreover, the RIA consistently outperforms the IA, CB2 and CDW approximations. As to four-body problems, the analysis is centered on single electron capture and transfer ionisation as well as on one-electron detachment involving two-electron atoms and ions. For single charge exchange and transfer ionisation in proton–helium collisions, emergence of dynamic electron correlations is evidenced due to their progressively rising importance with increasing incident energies. Static electron correlations are found to be crucial for accurate quantitative predictions on electron detachment in collisions of protons with negative hydrogen ions. Above all, the major role in detachment is attributed to consistency between perturbation potentials and the corresponding scattering wavefunctions. Usefulness of perturbative distorted-wave collision theories critically depends upon a judicious intertwining of a powerful set of analytical tools from mathematical physics with highly accurate and efficient computational methods. The presently employed analytical methods deal with two- and three-center bound-free atomic form factors evaluated with the Cauchy complex contour integration technique followed by the Feynman-Dalitz-Lewis integrals. The final numerical task is reduced to evaluation of integrals of dimensions ranging from one to thirteen (1D–13D). For higher dimensions (6D–13D) that are encountered in this work, stochastic method are employed and especially the adaptive and iterative exact Monte Carlo code VEGAS proves to be remarkably useful. For lower dimensions (1D–5D), deterministic quadrature rules are advantageously used and particularly the plurivariate fast Padé transform (FPT) is established as a well suited and robust method for benchmark computations. The FPT remarkably accelerates slowly converging Riemann partial sums of varying length from the trapezoidal quadrature rule by using the Padé approximant. This method is operationally implemented through the Wynn recursive epsilon algorithm with the benefits of having a stable, efficient and low storage computational method of unprecedented accuracy.