Obtaining an accurate mechanical model of a soft deformable robot compatible with the computation time imposed by robotic applications is often considered an unattainable goal. This paper should invert this idea. The proposed methodology offers the possibility to dramatically reduce the size and the online computation time of a finite element model (FEM) of a soft robot. After a set of expensive offline simulations based on the whole model, we apply snapshot-proper orthogonal decomposition to sharply reduce the number of state variables of the soft-robot model. To keep the computational efficiency, hyperreduction is used to perform the integration on a reduced domain. The method allows to tune the error during the two main steps of complexity reduction. The method handles external loads (contact, friction, gravity, etc.) with precision as long as they are tested during the offline simulations. The method is validated on two very different examples of FEMs of soft robots and on one real soft robot. It enables acceleration factors of more than 100, while saving accuracy, in particular compared to coarsely meshed FEMs and provides a generic way to control soft robots.