I want to write a LAMMPS script to calculate defect formation energy where the defects (both vacancies and interstitials) are located within/around the grain boundary.
More specifically, first, all atoms within 20 A of the grain boundary plane are identified as potential sites. Then, for vacancy formation energy simulations, an atom at a particular site, say i, is removed and the simulation cell is relaxed through an energy minimization.
Since, the number of atoms in this slab around the grain boundary will be many and we need to remove (add) one atom at a time for simulating a vacancy (interstitial), we need to implement a script which can select one of these sites at a time and remove it.
I can use the region command to select the atoms in the slab around the grain boundary. But how to store the atom ids in an array variable which I can call later iteratively in a loop and do the rest of the calculation ?
If you do have any idea about the implementation of the above using a LAMMPS script, can you please share it ? It would be helpful to have an example input file.
Thanks and kind regards,
I have done the same simulations and you can read the manuscript in a few days after it is published.
I used my self-written cpp code to iterate the loop (deleting or creating one atom in a pre-defined lattice site) and just used the lammps as MD engine to perform cell minimization and relaxation every time. If you need my cpp code or further help, just mail me.
It is possible to write a LAMMPS loop over all atoms, and then act on only those that are in the desired group. The action could be to set the atom type to that of a noninteracting atom, calculate the relaxed energy and then set it back. Something similar could be done for interstitials, by populating the interstitial lattice with noninteracting atoms.
But all of this is error-prone and inefficient. Better to write a driver script to handle adding and removing atoms. The Python interface is a nice way to do this, or you can do it by writing a script that edits a LAMMPS data file. It sounds like that is what Afshin did.
Thanks for your response. Yes, you are correct. I had a discussion with Afshin yesterday. He has developed a driver code in c++ for this purpose. I am waiting to see his script.