Gillespie Algorithm for Turing patterns on networks

Hello everyone,

I tried to implement the Gillespie algorithm for the study of Turing patterns from a three-species model on a network structure. However, the simulation results do not match my expectations. First, a pattern also occurs if diffusion is switched off (epsilon =0). Second, the resulting pattern / node differentiation seems to be purely random and does not show any similarity to the deterministic (temporally oscillating) pattern for large system sizes O. I don't see were my implementation is wrong....

The (deterministic) model was taken from here. Maybe I made a mistake with the derivation of the reaction rates? Or maybe something is wrong with the for-loops? The code is given in the attached file.

I would be very grateful for any help!

Alex