Hi all,
Recently I want to add the correct local normals to the Kolmogorov-Crespi (KC) potential based on the present version of KC potential (pair_kolmogorov_crespi_z.cpp) in USER-MISC package in LAMMPS.
The present version of KC potential in LAMMPS assumes that all normals are along the z-axis, which leads this potential only valid for flat surfaces (e.g., not valid for carbon nanotubes). So I want to add the correct local normals to this KC potential.
When implementing the code, I met a problem. Assuming there are two layers of graphene with the atoms in the bottom layer setting to type 1 and the atoms in the top layer setting to type 2. This present KC potential is invoked in the input file like follows:
pair_style hybrid/overlay kolmogorov/crespi/z 20.0
pair_coeff * * none
pair_coeff 1 2 kolmogorov/crespi/z CC.KC C C
As a result, only the neighbor list between types 1 and 2 is built and there is no neighbor list built for the atoms within one layer (i.e., types 1 1 and types 2 2). Therefore, to calculate the local normals, I need to build a neighbor list for the atoms within one layer. I tried to build this neighbor list by introducing another cutoff but I failed.
Does anyone has an idea for how to build this neighbor list? Thank you in advance!
Hi all,
Recently I want to add the correct local normals to the Kolmogorov-Crespi
(KC) potential based on the present version of KC potential
(pair_kolmogorov_crespi_z.cpp) in USER-MISC package in LAMMPS.
The present version of KC potential in LAMMPS assumes that all normals are
along the z-axis, which leads this potential only valid for flat surfaces
(e.g., not valid for carbon nanotubes). So I want to add the correct local
normals to this KC potential.
When implementing the code, I met a problem. Assuming there are two layers
of graphene with the atoms in the bottom layer setting to type 1 and the
atoms in the top layer setting to type 2. This present KC potential is
invoked in the input file like follows:
pair_style hybrid/overlay kolmogorov/crespi/z 20.0
pair_coeff * * none
pair_coeff 1 2 kolmogorov/crespi/z CC.KC C C
As a result, only the neighbor list between types 1 and 2 is built and
there is no neighbor list built for the atoms within one layer (i.e., types
1 1 and types 2 2). Therefore, to calculate the local normals, I need to
build a neighbor list for the atoms within one layer. I tried to build this
neighbor list by introducing another cutoff but I failed.
Does anyone has an idea for how to build this neighbor list? Thank you in
advance!
​please replace:
​pair_coeff 1 2 kolmogorov/crespi/z CC.KC C C
with:
pair_coeff 1*2 1*2 kolmogorov/crespi/z CC.KC C C
and the neighbor list handed to the kolmogorov/crespi/z calculation will
contain all pairs you are asking for.
for details of the syntax, please see the documentation of the pair_coeff
command.
Dear Axel,
Thanks! I think it will be OK if there is only KC potential. However, people usually use REBO+KC potential to simulate few layers of graphene. The command is like follows, as is in the documentation of pair_style kolmogorov/crespi/z:
pair_style hybrid/overlay rebo kolmogorov/crespi/z 14.0
pair_coeff * * rebo CH.airebo C C
pair_coeff 1 2 kolmogorov/crespi/z CC.KC C C
In this case, I can not use “pair_coeff 12 12 kolmogorov/crespi/z CC.KC C C” instead since it will change the interaction.
I know that the neighbor list is already built for types 1 1 in pair_airebo.cpp. My question is how to invoke this neighbor list in the new pair style I’m implementing?
Do you have any suggestion?
Thank you very much!
Dear Axel,
Thanks! I think it will be OK if there is only KC potential. However,
people usually use REBO+KC potential to simulate few layers of graphene.
The command is like follows, as is in the documentation of pair_style
kolmogorov/crespi/z:
pair_style hybrid/overlay rebo kolmogorov/crespi/z 14.0
pair_coeff * * rebo CH.airebo C C
pair_coeff 1 2 kolmogorov/crespi/z CC.KC C C
In this case, I can not use "pair_coeff 1*2 1*2 kolmogorov/crespi/z
CC.KC C C" instead since it will change the interaction.
I know that the neighbor list is already built for types 1 1 in
pair_airebo.cpp. My question is how to invoke this neighbor list in the new
pair style I'm implementing?
​you cannot. the neighbor list handed to the airebo pair style contains
*all* pairs, anyway. you would have to create new a new pair style for the
sole purpose of generating a suitable neighbor list request.
​
Do you have any suggestion?
​my suggestion is the only way. since you are already modifying the pair
style, it is up to you ​to modify it in such a way, that you skip over the
unwanted pairs.
Dear Axel,
Thanks for your information!
Maybe you are right, it’s easier to modify the pair style to skip over the unwanted pairs. I will try to modify it in such a way to see if it works.
Dear Axel,
I just checked that the codes calculate the normals correctly by using “pair_coeff 12 12 kolmogorov/crespi/z CC.KC C C”.
About your suggestion to skip over the unwanted pairs, I tried to set the setflag[i][j] = 0 when i=j. It really skipped the pairs for 1 1 and 2 2, but the neighbor list also did not build. How can I build the neighbor list for types 1 1 and 2 2 but to avoid calculating the forces between them? Can I just set the force to be zero (i.e., f[i][j]=0) when i=j?
Could you give me more hints? Thank you very much!
Dear Axel,
I just checked that the codes calculate the normals correctly by using
"pair_coeff 1*2 1*2 kolmogorov/crespi/z CC.KC C C".
About your suggestion to skip over the unwanted pairs, I tried to set the
setflag[i][j] = 0 when i=j. It really skipped the pairs for 1 1 and 2 2, but
the neighbor list also did not build. How can I build the neighbor list for
types 1 1 and 2 2 but to avoid calculating the forces between them? Can I
just set the force to be zero (i.e., f[i][j]=0) when i=j?
Could you give me more hints? Thank you very much!
i said "skip the pairs" not "skip the pairs of types". the skipping
*must* be done inside the ::compute() method. what you are doing is a
hack and will not do what you expect it should do. the clean and
rational approach would be to have a neighborlist with all necessary
pairs and then, build the required sublists from them or do nested
loops and skip over unwanted pairs.
overall, i have to warn you that my patience is very short with people
that a) do not follow my advice and b) continue asking for something
that i have recommended against or that is not possible. worst of all,
it appears to me, you have not really thought things through, and are
trying to do quick hacks. i strongly dislike that kind of programming
and will not provide advice along those lines.
Hi all,
Thanks for Axel’s suggestion, I think I find a way to calculate the local normals and skip the unwanted pairs. Firstly, I build a full neighbor list for all pairs and then build a sublist from the full list, just as what is done in pair_airebo.cpp. Secondly, I use the same switching functions to skip the unwanted pairs and I think it’s a good way to do this because KC potential is also designed for graphite system. To check the validity, I manually set the normals all along z axis (i.e., for flat graphene surface), and the binding energy and corrugation potential are extremely close to the results given by pair_style kolmogorov/crespi/z, which means one can use the original parameters safely.
When calculating the local normals, I loop over all owned atoms (i.e., from 0 to inum-1), for each owned atom i, its sublist has three neighbor atoms and I use three coordinates to calculate the normals. While when calculating the forces in the compute(), I use the full list and skip the unwanted pairs. Here I meet a problem. When looping over all owned atoms (i.e., from 0 to inum-1) and its neighbors (from 0 to jnum-1), the neighbor atoms of i can be a ghost atom, that means index of neighbor atom j > nlocal. In KC potential, I need local normals n[i] and n[j] to calculate the force. I already calculate n[i], however, when I try to calculate n[j], I found that not all ghost atoms has three neighbors, so the local normal n[j] can not be calculated correctly. I try to map the ghost atom back to local index (so they have three neighbors as it should be) by using atom->map(atom->tag[j]), but it can not get the correct index for some ghost atoms. I also tried to loop over owned + ghost atoms, i.e., from 0 to inum+gnum-1, but meet the same problem: some ghost atoms don’t have 3 nearest neighbors.
Does any one know how to get the correct local index of a ghost atom or has a better idea to calculate the normals for the ghost atoms?
Thank you in advance!
Hi all,
Thanks for Axel's suggestion, I think I find a way to calculate the local
normals and skip the unwanted pairs. Firstly, I build a full neighbor list
for all pairs and then build a sublist from the full list, just as what is
done in pair_airebo.cpp. Secondly, I use the same switching functions to
skip the unwanted pairs and I think it's a good way to do this because KC
potential is also designed for graphite system. To check the validity, I
manually set the normals all along z axis (i.e., for flat graphene
surface), and the binding energy and corrugation potential are extremely
close to the results given by pair_style kolmogorov/crespi/z, which means
one can use the original parameters safely.
When calculating the local normals, I loop over all owned atoms (i.e.,
from 0 to inum-1), for each owned atom i, its sublist has three neighbor
atoms and I use three coordinates to calculate the normals. While when
calculating the forces in the compute(), I use the full list and skip the
unwanted pairs. Here I meet a problem. When looping over all owned atoms
(i.e., from 0 to inum-1) and its neighbors (from 0 to jnum-1), the neighbor
atoms of i can be a ghost atom, that means index of neighbor atom j >
nlocal. In KC potential, I need local normals n[i] and n[j] to calculate
the force. I already calculate n[i], however, when I try to calculate n[j],
I found that not all ghost atoms has three neighbors, so the local normal
n[j] can not be calculated correctly. I try to map the ghost atom back to
local index (so they have three neighbors as it should be) by using
atom->map(atom->tag[j]), but it can not get the correct index for some
ghost atoms. I also tried to loop over owned + ghost atoms, i.e., from 0 to
inum+gnum-1, but meet the same problem: some ghost atoms don't have 3
nearest neighbors.
Does any one know how to get the correct local index of a ghost atom or
has a better idea to calculate the normals for the ghost atoms?
​this sounds like you need to modify the init_one() method in your pair
style to increase the cutoff distance retu​rned for a specific pair of
types to include an additional buffer distance to have more neighbors of
neighbors (including ghosts) included.
Hi Axel,
Recently I’m reading the code of pair_airebo.cpp and I have a similar question about the number of neighbors of ghost atoms. Please see details at https://sourceforge.net/p/lammps/mailman/message/35981788/. Do you think it is also caused by the small cutoff in the AIREBO potential?
Hi Axel,
As you suggested, I increased the cutoff from 2 Angstrom to 20 Angstrom to build the sublist for calculating the normals, but I still meet the same problem. Although the number of neighbors of a ghost atom increased up to 70, the nearest neighbors (within one carbon-carbon length) still less than 3 (only the nearest neighbors is meaningful to calculate the local normals).
After looking Huang’s question (https://sourceforge.net/p/lammps/mailman/message/35981788/), I’m confused about the meaning of the number of neighbors of a ghost atoms because of my limited understanding. Does it have the same physical meaning as a owned atom?