A data-driven surrogate framework to accelerate particle-resolved modelling of quasi-dilute suspensions of rigid, non-spherical particles in Stokes flow is introduced. A regularized-Stokeslet boundary element method (BEM) is implemented to compute hydrodynamic responses in canonical linear flows, focusing on the particle stresslet and angular velocity for spheroids, and additionally the chiral thrust for helicoidal particles. For spheroids, the BEM solver is validated against available analytical benchmarks (Faxén-type relations for the stresslet and Jeffery's theory for rotation), and parameter choices for surface discretization and regularization are selected through systematic convergence studies. For helicoidal particles, where no analytical solutions exist, accuracy is quantified via Richardson-style self-convergence, complemented by tests of linearity, frame objectivity, and chirality-dependent symmetries. The resulting datasets are used to train a neural-operator surrogate that maps local flow descriptors and particle configuration to the corresponding stresslet, rotation, and thrust at negligible evaluation cost. Across independent test sets spanning random orientations and flow types, the surrogate achieves median relative errors below 1% for the deviatoric stresslet (95th percentile below 3%) and comparable accuracy for angular velocity and thrust. The combination of validated BEM generation and fast inference provides a practical route to coupling complex particle shapes into mesoscale solvers such as the force-coupling method, enabling large-ensemble studies of microstructure and suspension rheology.