This paper presents a novel computational approach for SIMP-based Topology Optimisation (TO) of hyperelastic materials at large strains. During the TO process for structures subjected to very large deformations, and especially in the presence of intermediate density regions, the standard Newtonsolver (or its arc length variant) have been reported not to converge (refer to References [15,27,33]). In this paper, the new TO stabilisation technique proposed in [1] in the context of level-set TO, initially devised to alleviate numerical instabilities inherent to level-set TO, is extended for the TO by means of the SIMP method. The success of the methodology rests on the combination of two distinct key ingredients. First, the nonlinear equilibrium equations of motion for intermediate TO design stages are solved in a non-exact albeit consistent incrementally linearised fashion by splitting the design load into a number of load increments. Second, the resulting linearised tangent elasticity tensor is locally stabilised (regularised) in order to prevent its loss of positive definiteness and, thus, avoid the loss of convexity of the discrete tangent operator. This solution strategy is shown to be extremely robust in the context of density-based TO, where the constitutive law of the underlying evolving solid structure is a mixture of solid and void constituents, the latter classically defined by means of a fictitious strain energy. The robustness and applicability of this TO methodological approach are thoroughly demonstrated through an ample spectrum of challenging numerical examples, ranging from benchmark two-dimensional (plane stress) examples to larger scale three-dimensional applications. Crucially, the performance of all the final designs has been tested at a post-processing stage without adding any source of artificial stiffness. Specifically, an arc-length Newton-Raphson method has been employed in conjunction with a ratio of the material parameters for void and solid regions of 10 −12 .