Open frobnitzem opened 7 years ago
If you have a simple formula then I can try to implement this. I was always a little confused about this though. Is the slab correction still applied in the z direction, or in some skewed direction normal to one of the box sides?
i think that because LAMMPS requires that that the unit cell vector a of a triclinic cell is collinear with x, and the a-b plane in the xy plane, the z-direction is the only meaningful direction for a slab correction. also, in the orthogonal case, a fixed box size in z direction is required, so that would logically mean for triclinic cells, that only changes to the xy tilt factor would be allowed.
I'm creating a pull request that just uses the code as-is to apply the correction in the z-direction.
I believe this is correct whether or not there is tilt in xz and yz. Consider the cases of a slab within a cell that has a c vector that of [1,1,10] and one that is [0,0,10]. The first c vector simply means that the c-image of the slab is out of register with the lower one. In both cases, the correction is based on treating the image as a 2D sheet of charge, so the xz and yz tilts do not matter.
there may be tilt in xz or yz, but it must not change, as that would change the distance between slabs and that is assumed to be fixed.
my previous remarks about the tilt factors are a red herring. due to the way how LAMMPS defines a tilted box, the z-direction spacing is not affected. however, the changes required to support this are quite invasive, since the effect of slab_volfactor
has to be taken into account for the regular kspace calculation.
the suggested change of transferring the coordinates from lamda to cartesian coordinates works fine, but the rest of the computation, needs to take the change in volume into account.
According to the docs, http://lammps.sandia.gov/doc/kspace_modify.html
... which is supported by the checks at ewald.cpp:109 "Cannot (yet) use Ewald with triclinic box and slab correction" and pppm.cpp:192 "Cannot (yet) use PPPM with triclinic box and slab correction"
From my reading of Yeh and Berkowitz, subtracting the z-dipole should not be fundamentally more difficult in non-orthogonal box shapes, and is indeed implemented straightforwardly already on pppm.cpp:2944. The only issue I noticed when trying it out was that the coordinates, x[i][2], refer to scaled versions and need to be transformed back into Cartesian when the dipole is summed on line 2954.
Rather than do that, the simple formula in Ballenger, JCP 140, 161102 makes it really tempting to implement corrections for arbitrary box shapes by accepting a tensor, 'J' and computing the whole dipole.