So it turns out my interpretation of the problem was incorrect. I’ve used a debugger to monitor the functions in PairTable. Am I really the only person who has had a problem with this before?
First, here are parts of my input files, as a refresher to what I am doing.
3-d 2 metal surfaces coming in contact and then sliding
dimension 3
boundary p p p
units metal
atom_style atomic
read_data data
#pair potentials
pair_style table spline 299
pair_coeff 1 1 somemore_pots AlCo
pair_write 1 1 50000 r 1.16 16 alco AlCo
and here is “somemore_pots”
AlCo
N 299
1 1.11127 17.4427 52.4443
2 1.16419 14.8643 45.5055
3 1.21711 12.6264 39.5015
4 1.27002 10.6838 34.2826
5 1.32294 8.99832 29.7286
6 1.37586 7.53732 25.7515
7 1.42878 6.27278 22.2678
8 1.4817 5.1805 19.2131
…
I don’t think the “read_data data” is relevant to the issue.
The number “299” in the file, can be increased to say 10000, and this fixes all the problems, for a reason that later became apparent. Note, however, that there are only 299 entries in my potential file.
First, here is the plot of my potential file (col. 2 vs 3) against the file that I wrote using the pair_write command.
http://shell.cas.usf.edu/~kmclaugh/spline/spline.pdf
and zoomed in to the problematic region
http://shell.cas.usf.edu/~kmclaugh/spline/spline2.pdf
And to demonstrate that this is not just a spline problem, but also a problem which arises using “pair_style linear”
http://shell.cas.usf.edu/~kmclaugh/spline/linear.pdf
and zoomed in
http://shell.cas.usf.edu/~kmclaugh/spline/linear2.pdf
…
So here are my findings using the debugger.
When spline is called, everything works as expected. yp1, ypn are the correct negative forces read from the file, y[] and x[] both match the file perfectly. y2[] is created and appears to be reasonable.
Now, splint is called from compute table, but it doesn’t use the expected values for x each time it is called. In fact, the only reason I can imagine it being called at this point (pre-simulation) would be for the pair write command, but this is not the case as I could then match the values for x being used with the points on my pair_write graphs (let’s call this the output potential). So anyhow, I made a list of the x’s being used the first 6 calls:
1.11127
1.47885
1.771751
2.02267
2.2451
2.4485
At first, I couldn’t make any sense of these numbers, but then looking closer at spline.pdf and spline2.pdf, I noticed that these are the points which are in agreement between my input potentials and my output potentials. In fact, if I use the same input file, and a different potential file, I see the same behavior. My input and output potentials agree primarily at the listed points.
I tried replacing spline with linear, and made the linear graphs above, which appears to suffer the same fate. My potentials (input and output), agree at the above specified points which are then connected with lines. Though I haven’t confirmed the following, my current thought is that the spline routine is being called once to estimate the values at some given points (listed above). Lammps then creates a table using these x values versus the splint calculated y values, then using a second pass of splint to connect these new dots, upon request of the simulation, or “pair_write”.
Additionally, I have tried specifying the inner and outer cutoff of my potential file, but that leads to even more unexpected garbage, which is even uncorrected by using a large n.
At this point, the easiest way to correct the problem is to use a very large value of n (though this increases seek times when the simulation requests a potential). I am going to try to modify the code to correct this behavior, though I fear my intervention may undermine some other usage of pair_style table.
It would be good to hear back from anyone with experience, thoughts, ideas or advise on this matter.
Thanks,
Keith McLaughlin