Dear All

Recently I tried to reproduce the result of a paper about two honeycomb bilayer (PRB.85.184412 http://prb.aps.org/forward/PRB/v85/i18/e184412) . The paper used simple MC simulations ( in fact authors used exchange MC technique http://dx.doi.org/10.1143/JPSJ.65.1604) with classical Heisenberg model to derive susceptibility of this system. Authors used three different exchange parameters: 1-J1 as NN interaction in honeycomb surfaces 2-J2 as NNN interactions in honeycomb surfaces 3-Jc as face-to face intra-bilayer interactions.

I used exactly the same parameters of the paper (I use minus sign for exchange). My setup is as follows:

****************************************** ****************************************** Exchange parameters setting: J1=30.2 K, J2/J1=0.1 and Jc/J1=0.1.

Landee g-factor: g=2

steps for thermalization: 50000

steps for sampling: 1000000

lattice size: L=24 ****************************************** ******************************************

All of my parameters are the same as the paper except that I used more steps for sampling and for thermalization because the system is frustrated. The following plots are mine and from the paper ( I derived the data of fig(8,c) of the paper by using G3DATA):

The ALPS result (mine) [image: Inline image 1] The PRB.85.184412 result (please see fig 8_c of the paper):

[image: Inline image 2]

*As you can see the plot of mine has a maximum around 40 K while the plot of PRB.85.184412 has a maximum around 60-70 K and also the plot shows a minimum around 20 K.

W**hat do you think? what is my fault? Is it a bug or something else?*

If somebody like to repeat the calculation, these are my input files:

the python script:

#================================= import pyalps import matplotlib.pyplot as plt import pyalps.plot

#prepare the input parameters parms = [] for t in [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 22, 24, 26, 28, 30, 32, 34, 36, 38, 40, 43, 46, 49, 52, 55, 58, 61, 64, 67, 70, 75, 80, 85, 90, 95, 100, 110, 120, 130, 140, 150, 160, 170, 180, 190, 200, 210, 220, 230, 240, 250, 260, 270, 280, 290, 300, 310, 320, 330, 340, 350, 360, 370, 380, 390, 400 ]: parms.append( { 'LATTICE_LIBRARY' : "graphene_bilayer.xml", 'LATTICE' : "graphene bilayer lattice", 'S' : 1.5, 'T' : t, 'J1' : -30.2, 'J2' : -3.02, 'J3' : -3.02, 'g' : 2.00, 'THERMALIZATION' : 50000, 'SWEEPS' : 1000000, 'UPDATE' : "local", 'MODEL' : "Heisenberg", 'L' : 24 } )

#write the input file and run the simulation input_file = pyalps.writeInputFiles('parm2b',parms) pyalps.runApplication('spinmc',input_file,Tmin=300, MPI=4, mpirun='mpirun')

#get the list of result files result_files = pyalps.getResultFiles(prefix='parm2b') print "Loading results from the files: ", result_files

#print the observables stored in those files: print "The files contain the following mesurements:", print pyalps.loadObservableList(result_files)

#load a selection of measurements: data = pyalps.loadMeasurements(result_files,['Susceptibility'])

#make a plot for the magnetization: collect Magnetziation as function of T plotdata = pyalps.collectXY(data,'T','Susceptibility')

# convert the data to text file for plotting using another tool print "The results in plain text are:" print pyalps.plot.convertToText(plotdata)

print "Saving into file parm2b.txt" f = open ('parm2b.txt','w') f.write(pyalps.plot.convertToText(plotdata)) f.close()

#================================

the lattice file (xml format)

===================================== <LATTICES>

<LATTICE name="hexagonal lattice" dimension="3"> <PARAMETER name="a" default="1"/> <PARAMETER name="c" default="1"/> <BASIS> <VECTOR>a 0 0</VECTOR> <VECTOR>a/2 a*sqrt(3)/2 0</VECTOR> <VECTOR>0 0 c</VECTOR> </BASIS> <RECIPROCALBASIS> <VECTOR>2*pi/a -2*pi/a/sqrt(3) 0 </VECTOR> <VECTOR>0 4*pi/a/sqrt(3) 0 </VECTOR> <VECTOR>0 0 2*pi/c</VECTOR> </RECIPROCALBASIS> </LATTICE>

<UNITCELL name="graphene bilayer" dimension="3" vertices="4"> <VERTEX><COORDINATE> 0.000000000000 0.0000000000 0 </COORDINATE></VERTEX> <VERTEX><COORDINATE> 0.333333333333 0.3333333333 0</COORDINATE></VERTEX> <VERTEX><COORDINATE> 0.000000000000 0.0000000000 0.5 </COORDINATE></VERTEX> <VERTEX><COORDINATE> 0.333333333333 0.3333333333 0.5</COORDINATE></VERTEX>

<EDGE type="1"><SOURCE vertex="1"/><TARGET vertex="2" offset=" 0 0 0"/></EDGE> <EDGE type="1"><SOURCE vertex="3"/><TARGET vertex="4" offset=" 0 0 0"/></EDGE> <EDGE type="1"><SOURCE vertex="2"/><TARGET vertex="1" offset=" 1 0 0"/></EDGE> <EDGE type="1"><SOURCE vertex="2"/><TARGET vertex="1" offset=" 0 1 0"/></EDGE> <EDGE type="1"><SOURCE vertex="4"/><TARGET vertex="3" offset=" 1 0 0"/></EDGE> <EDGE type="1"><SOURCE vertex="4"/><TARGET vertex="3" offset=" 0 1 0"/></EDGE>

<EDGE type="2"><SOURCE vertex="1"/><TARGET vertex="3" offset=" 0 0 0"/></EDGE> <EDGE type="2"><SOURCE vertex="2"/><TARGET vertex="4" offset=" 0 0 0"/></EDGE>

<EDGE type="3"><SOURCE vertex="1"/><TARGET vertex="1" offset=" 1 0 0"/></EDGE> <EDGE type="3"><SOURCE vertex="1"/><TARGET vertex="1" offset=" 0 1 0"/></EDGE> <EDGE type="3"><SOURCE vertex="1"/><TARGET vertex="1" offset=" 1 -1 0"/></EDGE> <EDGE type="3"><SOURCE vertex="2"/><TARGET vertex="2" offset=" 1 0 0"/></EDGE> <EDGE type="3"><SOURCE vertex="2"/><TARGET vertex="2" offset=" 0 1 0"/></EDGE> <EDGE type="3"><SOURCE vertex="2"/><TARGET vertex="2" offset=" 1 -1 0"/></EDGE> <EDGE type="3"><SOURCE vertex="3"/><TARGET vertex="3" offset=" 1 0 0"/></EDGE> <EDGE type="3"><SOURCE vertex="3"/><TARGET vertex="3" offset=" 0 1 0"/></EDGE> <EDGE type="3"><SOURCE vertex="3"/><TARGET vertex="3" offset=" 1 -1 0"/></EDGE> <EDGE type="3"><SOURCE vertex="4"/><TARGET vertex="4" offset=" 1 0 0"/></EDGE> <EDGE type="3"><SOURCE vertex="4"/><TARGET vertex="4" offset=" 0 1 0"/></EDGE> <EDGE type="3"><SOURCE vertex="4"/><TARGET vertex="4" offset=" 1 -1 0"/></EDGE>

</UNITCELL>

<LATTICEGRAPH name = "graphene bilayer lattice"> <FINITELATTICE> <LATTICE ref="hexagonal lattice"/> <PARAMETER name="W" default="L"/> <PARAMETER name="H" default="1"/> <EXTENT dimension="1" size="L"/> <EXTENT dimension="2" size="W"/> <EXTENT dimension="3" size="H"/> <BOUNDARY type="periodic"/> </FINITELATTICE> <UNITCELL ref="graphene bilayer"/> </LATTICEGRAPH>

</LATTICES> =======================================================================

Thanks in advanced

Mojtaba Alaei