We give a new algorithm for computing the robustness of magic-a measure of the utility of quantum states as a computational resource. Our work is motivated by the magic state model of fault-tolerant quantum computation. In this model, all unitaries belong to the Clifford group. Non-Clifford operations are effected by injecting nonstabiliser states, which are referred to as magic states in this context. The robustness of magic measures the complexity of simulating such a circuit using a classical Monte Carlo algorithm. It is closely related to the degree negativity that slows down Monte Carlo simulations through the infamous sign problem. Surprisingly, the robustness of magic is submultiplicative. This implies that the classical simulation overhead scales subexponentially with the number of injected magic states-better than a naive analysis would suggest. However, determining the robustness of n copies of a magic state is difficult, as its definition involves a convex optimisation problem in a 4 n -dimensional space. In this paper, we make use of inherent symmetries to reduce the problem to n dimensions. The total run-time of our algorithm, while still exponential in n, is super-polynomially faster than previously published methods. We provide a computer implementation and give the robustness of up to 10 copies of the most commonly used magic states. Guided by the exact results, we find a finite hierarchy of approximate solutions where each level can be evaluated in polynomial time and yields rigorous upper bounds to the robustness. Technically, we use symmetries of the stabiliser polytope to connect the robustness of magic to the geometry of a low-dimensional convex polytope generated by certain signed quantum weight enumerators. As a by-product, we characterised the automorphism group of the stabiliser polytope, and, more generally, of projections onto complex projective 3-designs. R reg (ρ) := lim n→∞ R(ρ ⊗n ) 1/n .Unfortunately, computing R(ρ ⊗n ) seems to be a difficult task. For ρ being a singlequbit state, the tensor power ρ ⊗n lives in an 4 n -dimensional space, and the sum over the s i in the definition (2) of the RoM has to range over the 2 O(n 2 ) stabiliser states defined for n-qubit systems. Any direct implementation of the optimisation problem (2) will thus quickly became computationally intractable-and, indeed, Howard and Campbell [23] could carry it out only up to n = 5.The starting point of this work is the observation that there is a large symmetry group shared by ρ ⊗n and the stabiliser polytope. Thus, we formulate the optimisation