The construction of numerical schemes for the Kuramoto model is challenging
due to the structural properties of the system which are essential in order to
capture the correct physical behavior, like the description of stationary
states and phase transitions. Additional difficulties are represented by the
high dimensionality of the problem in presence of multiple frequencies. In this
paper, we develop numerical methods which are capable to preserve these
structural properties of the Kuramoto equation in the presence of diffusion and
to solve efficiently the multiple frequencies case. The novel schemes are then
used to numerically investigate the phase transitions in the case of identical
and non identical oscillators
