some kind of `wrap` function for PBC's?

[rkurchin]

Would it make sense to have this be part of the interface? If so, we would need to be able to dispatch on boundary condition types, which might require rejiggering the type parameters of AbstractSystem...

[Gregstrq]

Could you pleas elaborate your idea a bit more? What do you mean by wrap function?

[rkurchin]

Sorry yeah that wasn't very clear. I'll just do a concrete simplified example. Suppose you have periodic boundary conditions in 1D on a "box" of length 2. If you had a coordinate of a particle at point 3, it should "wrap" back to the point 1 (essentially just a modular arithmetic operation). This is a pretty standard thing for packages that deal with periodic BC's (e.g. Python's ASE, or the Xtals.jl package), but it's also very possible it's not actually important enough to belong in a high-level interface like this, hence I opened the issue to discuss :)

[eahenle]

From my (admittedly biased) perspective, periodic boundaries are very important. Can't handle crystalline materials without them. Molecular dynamics simulations usually use them, too.

The math is, as you say, not complicated. In Xtals we keep positions in fractional coordinates so the math simple: anything less than 0 gets 1 added to it, and anything greater than 1 gets 1 taken away. (+1.5, 0.5, -1) maps to (0.5, 0.5, 0).

[mfherbst]

I think if you have PBC and you add an absolute coordinate that is "larger than the border" it should always wrap. That is the point as @eahenle has described it. If you don't want that you first "disable" periodicity and then add the position.

Naturally one could also have a function to do these subtleties (like disable wrapping and favour extension of the system), but I think this belongs to another package using this interface.

[rkurchin]

Just noting that this is probably quite related to #5

[Gregstrq]

I think if you have PBC and you add an absolute coordinate that is "larger than the border" it should always wrap. That is the point as @eahenle has described it. If you don't want that you first "disable" periodicity and then add the position.

Technically, I imagine certain subtleties associated with fast implementations. For example, it may be better first to update all the positions and then to wrap them all back inside the boundary. @jgreener64, how is it usually done in MD simulations?

By the way, how the coordinates are wrapped back if the bounding box is skewed (not rectangular)? It seems that the only approach is to convert coordinates to fractional.

[jgreener64]

@jgreener64, how is it usually done in MD simulations?

It depends on the software, but usually after one or all coordinates are updated as you say.

I think there might be too many subtleties to support this in a way that pleases everyone. We also haven't made setters part of the interface, I think with decent justification, so it is unclear how we would modify the coordinates of a system in a general way.

[rkurchin]

These are all fair points. I really opened this issue as a "thought that popped into my head" and have no strong feelings about whether it should be implemented here or not. Certainly given that we're not doing setters, at the most if we had it it would just return the wrapped coordinates rather than modifying them, but at that point I'm not sure there's all that much use to having it in the interface...

[mfherbst]

I agree with @rkurchin and @jgreener64. I suppose we close this for now and revisit once we move towards setters.

[cortner]

It feels all these methods to manipulate systems could become part of AtomsBuilder

[rkurchin]

Agreed – I'm fine with this being closed and having an issue over there instead as a reminder

[cortner]

I'm planning to move over the rest of JuLIP functinality to AtomsBuilder over the coming week. If you file an issue there and describe how you think wrap should work then I'll try to match that.