Hi, I am trying to build my own atomistic simulator using kimpy. To begin with, I have been able to compute the quantities such as partial energy and partial forces using kimpy. However, I have done this only for a finite set of atoms. Now, I wish to extend my simulator for an infinite set of atoms, i.e., with periodic boundary conditions. For this, as I understand, my neighbour list function will have to be modified to check for neighbours in the neighbouring boxes.
My questions are:
What should my neighbour list data structure look like? Keeping in mind that it shouldn’t just have the neighbours that are inside the simulation box but also the ones that are outside in the neighbouring boxes. Also, it must be compatible with kim-api.
How are the quantities such as energy and forces computed by kim-api in the presence of periodic boundary conditions?
Is there anything else that I need to take care of? I will be grateful for any advice.
kimpy currently allows two types of neighbor list.
The first one is shipped together with kimpy. It is implemented in C++ and has a Python wrapper. An example usage is in the ASE KIM calculator.
The second approach allows you to use python objects. As you said, the data structure of the neighbor list object and the get_neigh function should conform to some predefined pattern to be compatible with the kim api. See this for an example, where the ASE raw neighbor list is wrapped as a neighbor list for kim.
In both cases, compute_args.set_callback registers the neighbor list object and get_neigh function; that’s where the data structure is checked.
The way PBCs are dealt with in kim-api is the same as in LAMMPS. Briefly, the atoms within the box are called contributing atoms, and their periodic images (called non-contributing atoms) outside the box needs to be provided by a simulator to satisfy the PBCs (for example, created via the create_padding function from the kimpy neighbor list). The total energy will be the sum of the energies of contributing atoms, and the force on an atom (either contributing or non-contributing) will be the negative derivative of this total energy w.r.t. the coords of the atom.
Note, the forces are partial forces, and the forces on non-contributing atoms need to be added back to their corresponding contributing atoms to get the total force on atoms in the box. See this kim-api doc page for more information.
As written in the kim-api documentation partialEnergy is the sum of partialParticleEnergies of the contributing atoms (partialParticle Energy of the non contributing atoms is taken to be zero).
Is this true for any potential? I have noticed that it isn’t. Consider the following two cases:
Case 1: When LJ potential is used, I find this to be true. So when i sum the partialParticleEnergy for only the contributing atoms, as expected, it comes out to be equal to partialEnergy (that i get from KIM).
Case 2: However, when Tersoff potential is used, I get non-zero partialParticleEnergy for non-contributing atoms. So when i sum the partialParticleEnergy of only the contributing atoms, the sum doesn’t come out to be equal to the partialEnergy. It should have according to the kim-api doc. Surprisingly, when I sum the partialParticleEnergy of all atoms (contributing and non-contributing), the sum comes out to be equal to the partialEnergy.
Could you clarify why this may be happening? Or am I missing something here?
The sum of the paritalParticleEnergy needs to be equal to the partialEnergy, and the partialParticleEnergy for a non-contributing particle should be zero. So, it seems there is a bug in the Tersoff model.
I’m tagging Dr. Tobias Brink, contributor of the KIM Tersoff model. @brink any idea?
Hello. It seems that the requirement of non-contributing atoms having zero energy was added or formalized after the Tersoff code was implemented. Indeed, it is a port from LAMMPS, where (if I remember correctly) the correct distribution of energies from ghost atoms to “real” atoms is done by the core code later and not by the potential. There does not seem to be a test in KIM to catch it, or is there?
It seems that the simulator would know best which ghost particle corresponds to which real particle (if at all), so simply redistributing the energy without changing anything else cannot be done by the model driver. That means I would have to look at the places where particle energies are assigned to make sure that they are only assigned to contributing particles. That would change the formulation of the potential, we’ll have to see how to do that correctly.
From the KIM side, it would be good to have some mandatory tests that check for these requirements, it should be simple enough. It is easy to fall into this trap, since the arrays for the particle energy (have to) contain entries for the non-contributing particles, it is just that they should not be used!
Actually, it seems we already fixed that a while ago! See the comment in pair_tersoff.cpp regarding “Non-symmetric assignment of bond energy”. Are you using the newest version of the driver? If yes, do you still have the problem? If yes, can you provide a simple test case that we could use to easily reproduce and debug the issue?
Actually, there was a test for this added to the KIM API’s utilities in this commit. I just ran this test on Tersoff_LAMMPS_Tersoff_1988_C__MO_579868029681_004 and it passes in both cases (a configuration with all contributing particles and a configuration with some non-contributing particles). This test confirms that:
The sum of partialParticleEnergy for contributing atoms adds up to partialEnergy
There are no non-contributing atoms with non-zero partialParticleEnergy
So, there should be no problems if the latest version is really being used.