I would like to study the thermodynamics of the Fe-C alloy system using the eam/fs potential of Hepburn and Ackland. Specifically, I need the chemical potential as a function of composition for any selected temperature. Since Carbon sits at interstitial sites of the BCC Fe lattice I need the fix gcmc command of lammps. I won’t describe the several tests I did with the insertion of C because, well, it’s boring. Instead let’s focus on the most recent simulations.
The test run I did consisted of a starting configuration of 2000 Fe atoms, that is, no carbon. The initial configuration was equilibrated and the temperature is T=1774. The input file, which you will notice is quite simple, is:
Essentially, you are running NVT MD on a metallic crystal, and at the same time you are stochastically displacing atoms according to Metropolis Monte Carlo. You observe the kinetic energy (temperature) to rise by about 10% over 100 timesteps. I think the problem is with the center-of-mass motion of the system. All those random displacements produce random changes in total momentum. Over many MC moves (5000/atom), the total momentum becomes large, as does the kinetic energy. This will result in an increase in the reported temperature, even if the true temperature (with total momentum subtracted) is not increasing. I boldly predict that if you change this line:
I don’t think we we ever got to the bottom of your problem. We recently made some changes to fix gcmc that might have fixed it. Also, I just added two gcmc examples, LJ and CO2. They are probably on github now and will be in the next release of LAMMPS. Maybe they will give you some ideas.