I want to do cascade irradiation on In2O3, so based on what I found in references I should apply zbl potential for short-range distances and Buckingham and coulombic potential for long-range.
since I tried splicing and it did not work somehow, I started to use pair write for each of them separately and merge files together.
What I did: first use pair-style hybrid zbl 0.1 1
and pair- write ( start with 1000000 and a couple of times reduce it even to 10 numbers ) for 0.1 to 1 distances
then use pair-style hybrid/overlay buck/coul/long 10 10
and again pair write ( start with 1000000 and a couple of times reduce it even to 10 numbers ) for 1 to 10 distances
then I prepared a merged file from 0.1 to 10 and used it as a pair coeff and my pair style was table linear ewald …
so the problem I have is I faced a lot of warnings when LAMMPS is trying to calculate force and check with my potentials files, force values are inconsistency with -dE/
So, this inconsistency warning per se may not be necessarily a problem (see pair_style table and use of RSQ). It could be just that the precise values of energy and force LAMMPS is using are slightly different than yours as it creates new data by interpolation. As it could be a matter of difference in very-far-away decimals, it may not be relevant - you can use the pair_write to print the tables of your simulation with the “merged” tables and compare the tables you are inputing versus the tables lammps are outputting for some insights. Furthermore, you can also save a trajectory and do a rerun (see rerun command — LAMMPS documentation) with the two sets of tables (the ones you input and the ones LAMMPS is using) to get a glimpse if this would cause a significant difference in the values of the observables as well as in the instantaneous fluctuations. These are some tests that came to mind now, but maybe there are more proper ways to evaluate if the difference in decimals is significant.
One thing that I think will be a problem though is the fact that you are trying to tabulate the coul/long energy and forces: I am not very familiar with long range solvers (ppm, ewald and etc), but you should not have in your tables the effect the coul/long would have for distances larger than the cutoff, so you will be losing information. I think that the idea of the coul/long is that it computes the Coulombic energies/forces within the cutoff region using the standard expression (+ possible modifications due to the long range solver) and also considers this interaction to exist in some fancy fashion for pairs sitting at distances > than the cutoff. This will for sure not be in your table for example.
For proper energy conservation and reasonable behavior, you want potential functions to be continuous and with a continuous derivative. This will not be the case for just slapping chunks from two different potential functions together. Usually, one uses a switching function to smoothly switch from one potential function to the other. The python scripts in the tools/tabulate folder can be used for implementing this, but they have to be used with care and a good understanding of the underlying physics.
Your explanation of how pair styles for long-range Coulomb works is incorrect. With a lattice solver (Ewald or PPPM) you take the Coulomb function and add and subtract a damping function (erfc() based) and then split the summation of the charge interaction into a real space and a reciprocal space part. The former adds the damping, the latter subtracts it. This turns a single conditionally converging sum into two well converging sums. The damping function is chose in such a way that the damped Coulomb is negligible at the cutoff distance. It is thus very much possible to tabulate this function, but additional care is needed. The charge information has to be included in the pair_write command and each atom with a different charge needs a different atom type and tables have to be created for each pair of atom types, which may be quite a few. Also, the table needs to indicated that it requires the Ewald or PPPM kspace style to complement and thus complete the Coulomb forces for the atoms beyond the cutoff.
Ahh… I was reasoning within analytical functions: I mean, if she was summing two analytical functions and their derivatives it should be equivalent to summoning a pair style given by the sum of these same two analytical functions. But I guess that since ultimately the tables corresponds to piecewise functions, summing them would be problematic… I have never thought of that actually
Ah, I didnt know it was possible to use a table with a long range solver. The specifics of what the solver does with the original Coulombic function before the cutoff is still a mistery to me though I would need to read more carefully a book to learn it properly
Wait, in fact the problem you are pointing out is that when we declare the two different potential functions (even if analytical) in a given pair_style (like, for example, the lj/cut/coul/long), there is a function (like the equivalent of a scaling factor) mediating how much each of them contributes at a given pair distance (assuming the potentials are all pairwise)?
I see. Good to know that this occurs for some combinations of potentials though - I will add it to my “LAMMPS notes” together with the other stuff you taught me
Thanks for your response.
I tried to run a simulation using a pair-style table and RSQ and added pair-write to see the potential lammps used. Then I plotted potential and force from my input file and also the file pair-write command generated. Although warnings were related to the force, potential plots are different and force plots are more matched. I do not know what causes this issue.
considering I am using pair-style table linear ( or even splined) 1000
and in pair-coeff I used merged files so based on that LAMMPS should read the file, considering the same potentials, and calculate forces so I should see matching for potential files ( merged and the one LAMMPS gave me ) and see differences in forces, am I right?
So, if I understand correctly, you are trying to built a tabulated potential that corresponds to a pair_style that cherishes many different potentials simultaneously (Buckingham, coulombic and this zbl thing) for each given pair of atom types. To do this, you initially created a table for the zbl and then another table for the buck/coul/long using the pair_write command. You then summed the values of energy and forces (3rd and 4th column) of these two tables aiming to generate a unique potential table.* You did this to get the potential for every single pair of atom types and are now using these tables to do the simulation you want to do.
*(or you summed the energies (3rd columns) and did some derivation to get the new values of forces)
Notably, tabulated potentials are non analytical. They correspond to piecewise functions that will be optimally fit according to the (pair distances, energy) and/or (pair distances, forces) points you have.
From my understanding of this now you will encounter some problems in your approach (regardless of this warning you are getting):
From what Axel said, some pair_styles that implicate in more than one potential function for a same atom pair contemplate a switching function (meaning that the individual potentials are not fully “activated”/contributing for the values of energy/force at all pair distances - basically, the full potential is not always a simple sum of the individual ones). I have no idea how you are going to get rid of this problem if your original pair_style fall within this rule
I am not sure (still ruminating) about summing the potential energies/forces of different potentials in attempt to make a unique potential/force table corresponding to their sum. I mean, if you were summing analytical functions and their derivatives, there should be no problem from basic derivation rules. But if they are table stuff (non analytical), it is in the realm of possibilities that your final potential will not cherish, within each bin, spline functions equivalent to the summation of those originally fit in the bins of the two individual potentials. This could maybe cause differences in what you are expecting. In fact, I am lagging really hard right now but I am not sure it is okay even if you could reinforce the summation of the piecewise functions precisely since I am not sure this would lead to continuous force profiles (if you derive the potential) (maybe this is what Axel meant)
Also, I just noticed you are using this /overlay option and I have no idea what it does (or if it poses problems on its own for the tabulating), so maybe you will have extra problems because of this? Better to check
Once you fix all of that (and also that you are also sure that you are building the table in a way that indeed encapsulates the set of potentials you want), you can determine whether or not this error is a matter of precision in the values coming from the LAMMPS interpolation when using your tables versus the values at the bins or not.
I am not doing any summation. what I want is to apply zbl potential from 0 to 2 and then switch to buck/coul/long from 2 to 10 angstrom. So what I did was apply ZBL potential and use pair-write for distances below 2 A and then separately apply buck/coul/long with cutoff 10 and used pair-write for 2 to 10 A distances. Then I merged these two files and had a potential file from 0 to 10 A distance. Considering this do you think I am right to tabulate the potential I want?
I think that’s possibly worse though, because then you will have a huge discontinuity of everything in the region where the two files meet. (I am guessing the shape of the zbl potential is verry different than the other ones in the range where data meet?). I mean, for sure lammps will find a way to unify this data in the interpolation, but it will be a very weird/abrupt change probably. Unless you are okay with this? I mean, I don’t know how exactly your potentials are supposed to look like (I have 0 experience with that you are trying to do).
EDIT: What if you carefully tailor a switching function like the one Axel meant in a common range where you allow the two functions to exist?
I would need to read more carefully a book to learn it properly
