Determining the trajectories of multiple acoustically and hydrodynamically interacting as well as colliding particles is one of the challenges in numerical acoustophoresis. Although the acoustic forces between multiple small spherical particles can be obtained analytically, previous research did not address the particle-particle contacts in a rigorous way. This article extends existing methods by presenting an algorithm on displacement level which models the hard contacts using set-valued force laws, hence allowing for the first time the computation of a first approximation of complete trajectories of multiple hydrodynamically and acoustically interacting particles. This work uses a semi-analytical method to determine the acoustic forces, which is accurate up to the dipole contributions of the multipole expansion. The hydrodynamic interactions are modeled using the resistance and mobility functions of the Stokes' flow. In previous experimental work particles have been reported to interact acoustically, ultimately forming stacked lines near the pressure nodes of a standing wave. This phenomenon is examined experimentally and numerically, the simulation shows good agreement with the experimental results. To demonstrate the capabilities of the method, the rotation of a particle clump in two orthogonal waves is simulated. The presented method allows further insight in self-assembly applications and acoustic particle manipulation.