We propose a new algorithm for simulating noncommutative phi-four theory on the fuzzy sphere based on, i) coupling the scalar field to a U(1) gauge field, in such a way that in the commutative limit N\longrightarrow \infty, the two modes decouple and we are left with pure scalar phi-four on the sphere, and ii) diagonalizing the scalar field by means of a U(N) unitary matrix, and then integrating out the unitary group from the partition function. The number of degrees of freedom in the scalar sector reduces, therefore, from N^2 to the N eigenvalues of the scalar field, whereas the dynamics of the U(1) gauge field, is given by D=3 Yang-Mills matrix model with a Myers term. As an application, the phase diagram, including the triple point, of noncommutative phi-four theory on the fuzzy sphere, is reconstructed with small values of N up to N=10, and large numbers of statistics.