Question on pair-style table

Dear LAMMPS users,

Hi,

I’m trying to use pairstyle table. In the manual there is some explanation on the file that we should use as the input of the potential. First, should I necessarily include force component as well as potential at each distance ? I tried to do the simulation without including force at each distance and it only give me warning as follow:

W****ARNING: 5800 of 5800 distance values in table with relative error
over 1e-006 to re-computed values (…/pair_table.cpp:490)

The formula that I have for the potential is in a way that I prefer not to calculate its derivative to calculate force because of its complexity. Moreover, since we use Verlet integration I doubt that the force is a requirement. Rather I guess it is only for the accuracy.

Another thing is that the simulation cannot change the box too much too reach to the desired pressure:

Step Lx Ly Lz Press Pxx Pyy Pzz Pxy Pxz Pyz PotEng Temp Density Volume
0 7.6090002 6.5895872 2.908 -382581.75 -366914.02 -367016.96 -413814.26 -14.134249 1.1799159e-010 1.4346748e-010 2453.3368 0 3.1955163 145.80762
1 7.6022536 6.5837428 2.905092 -382845.79 -367271.63 -367373.81 -413891.94 166.77058 -3.0889215e-011 1.6468195e-010 2451.0758 0 3.2043958 145.40358
2 7.5955024 6.5778945 2.902184 -383093.48 -367615.77 -367716.7 -413947.97 351.38823 -7.8587978e-012 3.1914386e-011 2448.8139 0 3.213312 145.00012
3 7.5921245 6.5749683 2.90073 -383211.09 -367782.67 -367882.86 -413967.73 445.15562 -1.5443877e-010 4.4368037e-011 2447.6824 0 3.2177845 144.79858
4 7.590435 6.5735048 2.900003 -383268.21 -367864.76 -367964.54 -413975.33 492.36974 1.1040666e-010 1.9763907e-010 2447.1166 0 3.2200244 144.69785
5 7.5904284 6.573499 2.9000002 -383268.43 -367865.08 -367964.86 -413975.36 492.5545 4.9710024e-011 -1.9245347e-010 2447.1144 0 3.2200331 144.69746
6 7.5904282 6.5734989 2.9000001 -383268.44 -367865.09 -367964.87 -413975.36 492.56028 2.3625709e-011 4.3044889e-011 2447.1143 0 3.2200334 144.69745
7 7.5904281 6.5734988 2.9 -383268.44 -367865.09 -367964.87 -413975.36 492.56316 8.6666019e-011 1.9089766e-010 2447.1143 0 3.2200335 144.69744
8 7.590428 6.5734987 2.9 -383268.44 -367865.09 -367964.88 -413975.36 492.56461 -8.9931363e-011 7.8944447e-011 2447.1143 0 3.2200336 144.69744
9 7.590428 6.5734987 2.9 -383268.44 -367865.09 -367964.88 -413975.36 492.56479 -3.7301732e-011 -5.2482691e-010 2447.1143 0 3.2200336 144.69744
10 7.590428 6.5734987 2.9 -383268.44 -367865.09 -367964.88 -413975.36 492.56488 1.9860963e-011 -7.5659898e-011 2447.1143 0 3.2200336 144.69744
11 7.590428 6.5734987 2.9 -383268.44 -367865.09 -367964.88 -413975.36 492.56488 1.9860963e-011 -7.5659898e-011 2447.1143 0 3.2200336 144.69744
Loop time of 0.0830579 on 1 procs for 11 steps with 14 atoms

Finally, nothing change during equilibration toward 0 pressure:

Step Lx Ly Lz Press Pxx Pyy Pzz Pxy Pxz Pyz PotEng Temp
0 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
1000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
2000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
3000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
4000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
5000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
6000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
7000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
8000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
9000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
10000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
11000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
12000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
13000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
14000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
15000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
16000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
17000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
18000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
19000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100
20000 7.590428 6.5734987 2.9 -382028.03 -366717.35 -367001.1 -412365.65 212.60573 271.53689 -534.98751 2447.1143 100

Any comment would be really helpful.
Thanks,
Ali

I would definitely be worried about that warning. It implies to me that best case your potential tabulation is not accurate enough and worst case it is complete bogus. Did you double-check it (by plotting or whatever)? Furthermore, please attach your input script, otherwise it is impossible to help you any further.

Hi Stefan,

Thanks for the response. Do you know if there is any way I can just use a tabulated file containing potential value without forces ? I’m asking because for the potential that I intend to use, doing the derivative is not so easy I think.

Thanks,

Ali

ali,

the basic rule to follow here is: when it is not explicitly stated in the LAMMPS documentation, that one can do something different than how it is described, this option is not available. a prerequisite of adding features to LAMMPS is that they need to be documented. this means, that when something is not documented, it is not available. the only case where this doesn’t apply strictly is when combining multiple commands, as people may not have been aware of later added commands and their incompatibility.

that said, i have to wonder. how can somebody state, that MD and specifically the velocity verlet algorithm doesn’t require forces? forces are explicitly mentioned in the verlet formulation and in any MD time integrator i am aware of. the story is different for MC, but there you don’t use time integration and thus forces are not usually required, only energy (which in turn is for the most part a diagnostic in MD). i am surprised to see this kind of gross oversight from somebody who has been so extremely nitpicky about what is given in the LAMMPS documentation and what not. whatever little credibility you had left with me is now gone.

when forces are not easy to compute explicitly, they can be easily approximated numerically. while this is not as good as an analytical expression (numerical differentiation adds noise and is best done in higher numerical precision than what is later used in the table to reduce noise), it is a hell of a lot better to do this from the analytical energy curve than expecting the MD code to do this from the already discretized tabulated energies.

axel.

Axel,

The reason that I asked the question was that why does lammps need Forces as well as potential. Because, I think since there exist a lot of approximation in the way we implement pairstyle table, you could include some approximation to calculate force as well. For many times, I tried to ask my questions very briefly on the form, and due to that I think several times you thought that my questions are simple, basic. Although it is a year that I’m working with LAMMPS and I really enjoy it, I have experience using other MD software such as NAMD and GROMACS before that time. Moreover, It is almost a long time that I searched the community before asking the question. However, sometimes I missed some previous comments unintentionally. But, It doesn’t mean that I do not spend time and respect you. About the Verlet algorithm, I didn’t mean that it is not used in formula. what I was saying and didn’t mention to make my Email short, is that I’m using MATLAB to build the tabulated file. and since the formula that I have for the potential is so big, it takes a while to use diff command. Rather than I can use some approximation to calculate. But, I was asking did you include such approximation in LAMMPS software or not ?

Finally about your last sentence, i.e. whatever little credibility you had left with me is now gone. I’m not happy to hear that, and I cannot blame you because your comments really helped me before. But, at least I expect you to read the previous Email that I sent one more time as well as this Email before expressing such comment.

Kind regards,

Ali