Kinetic Monte Carlo simulations for the diffusion of molecules on a substrate
The molecule shown in the image (hexaiodo-substituted macrocycle cyclohexa-m-phenylene (CHP) ), when deposited on a noble metal substrate such as Cu(111) , Ag(111) or Au(111), at room temperature looses the I atoms and starts diffusing.
The relative probability of diffusion and of binding to a neighboring molecule determine the shape of the network that will be obtained. The experiments performed at Empa [ http://dx.doi.org/10.1021/ja107947z J. AM. CHEM. SOC. 2010, 132, 16669–16676 ] show that on a Cu substrate dendrites will form while on a Au substrate 2D networks will form.
During the execution the program shows snapshots of the positions of the molecules.
Molecules free to diffuse will be represented via blue dots. Molecules that “irreversibly” formed a bond with a neighboring molecule will be represented by red dots. A few snapshots of the system during time evloution are saved as images.
Execute the program
The program asks you for some input:
update graph each steps coverage temperature in K diffusion barrier binding barrier
update graph each steps 1000 coverage 0.3 number of steps 130000 temperature in K 300 diffusion barrier 0.1 binding barrier 0.1
Observe how events occur. Observe how time evolves. Did the job perform all the 300000 steps? Why? Observe the patterns obtained.
update graph each steps 1000 coverage 0.3 number of steps 300000 temperature in K 300 diffusion barrier 0.3 binding barrier 0.4
Do you notice differences in the way events occur? How is evolving time compared to the previous case? How does the final geometry differ from the previous case?
#### MAIN KMC LOOP t=0 #### at the beginning we have to check possible events for all molecules tobeupdated=[iu for iu in range(len(molecules))] possible_events= for i in range(nsteps): #### check possible events for selected set of molecules possible_events=possible_events+events(molecules,tobeupdated,nx,ny) if possible_events==: print "no more events possible" break #### compute total rate (can be imporved) R=total_rate(possible_events,rates) #### decide which event to apply selected_event=find_event(R,rates,possible_events) #### apply the event end update partially the list of events possible_events,tobeupdated=apply_event(molecules,selected_event,possible_events) #### update time rho2=np.random.random() dt=-np.log(rho2)/R t=t+dt
is very simple and reflects the basic steps of the KMC approach. On the contrary, the function needed to create the list of events is quite complex:
def events(m,selection,nx,ny): allevents= set= #### list of first neighbors (relative position) set.append([[0,0,1,0],[1,0,0,0],[0,1,1,0],[-1,1,1,0],[-1,0,0,0],[-1,0,1,0]]) set.append([[1,-1,-1,0],[1,0,0,0],[1,0,-1,0],[0,0,-1,0],[-1,0,0,0],[0,-1,-1,0]]) for i in selection: le= #### consider only molecules that are not binded . . .