The Monte Carlo technique is used to model surface diffusion on energetically heterogeneous surfaces. Simulations are performed for the bivariate barrier model with a concentration c(H) of high barriers distributed at random. The variation of the chemical and tracer diffusion coefficients with particles coverage and with c(H) is discussed both for interacting and noninteracting particles. Two different regimes are found for c(H) below and above the critical concentration for bond percolation. Above percolation the diffusion coefficients are found to vary as c(H)(-1/T), a result which is not predicted by effective medium approximation (EMA).