Example of Kinetic Monte Carlo simulation of ZGB model
Download latest version here: test_zgb.py
and ZGB
model
Original ZGB paper: https://journals.aps.org/prl/abstract/10.1103/PhysRevLett.56.2553
from susmost import load_lattice_task
from susmost import mc
from susmost.mc.metropolis import kmc_step, kmc_rate
import numpy as np
lt = load_lattice_task('models/ZGB')
lt.set_property('O_cnt', {"O":1}) # counter of O coverage
lt.set_property('C_cnt', {"C":1}) # counter of CO coverage
assert lt.states_count == 3 # 3 lattice site states
m = mc.make_metropolis(lt, 20, [1.0], 1.0) # lattice size is 20x20; temperature and k_B are irrelevant here so set to 1.0
empty_idx = 0
C_idx = lt.find_state_index('name','C')
O_idx = lt.find_state_index('name','O')
yCO = 0.54 # critical values y1 = 0.389, y2 = 0.525
rCO = 4*yCO/(1+3*yCO) # reaction rate is specified separately for each oxygen dimer orientation, so it should be adjsuted
m.add_singe_site_step(empty_idx, kmc_step(kmc_rate(rCO, 0, 0), C_idx, -1, 0)) # CO adsorption, 0-th counter
for ei in m.regr.nn_edges_range(): # nearest neighbors in all directions
m.add_pair_step(empty_idx, ei, empty_idx, kmc_step(kmc_rate(1.0 - rCO,0,0), O_idx, O_idx, 1)) # O2 adsorption, counter index == 1
m.add_pair_step(O_idx, ei, C_idx, kmc_step(kmc_rate(1000.0,0,0), empty_idx, empty_idx, 2)) # O + CO reaction, 2nd counter
m.add_pair_step(C_idx, ei, O_idx, kmc_step(kmc_rate(1000.0,0,0), empty_idx, empty_idx, 2)) # CO + O reaction, 2nd counter
m.setup_kmc(1000.0) # r0 == 1000.0
m.params_log_period = 100
O_cnt_idx = m.param_idx('O_cnt')
C_cnt_idx = m.param_idx('C_cnt')
for i in range(1000):
m.run(1000000)
m.save_as_xyz("zgb.xyz", comment='E')
# print mean O coverage, CO coverage, acceptance rate, 2nd counter (number of CO + O reactions
print (i, np.mean(m.params_log[:,O_cnt_idx]), np.mean(m.params_log[:,C_cnt_idx]), np.mean(m.params_log[:,-4]), m.kmc_counters[2].accept)
Visualization of zgb.xyz with y_CO = 0.530: