Add periodic (minimum-image) distance support to apd_system#5
Add periodic (minimum-image) distance support to apd_system#5ligerzero-ai wants to merge 1 commit into
Conversation
Introduce an opt-in `periodic` flag on `apd_system`. When enabled, every seed-to-point distance uses the minimum-image convention across the rectilinear domain, so the anisotropic power diagram becomes periodic and grains may wrap across the domain boundary. - new `periodic=False` constructor flag (fully backward-compatible) - `_displacement(y, x)` helper computing the minimum-image displacement as a KeOps LazyTensor via `dy - L * round(dy / L)` - all four `D_ij` computations (assemble_apd, OT_dual_function, check_optimality, adjust_X) route through `_displacement` - new public `grain_of(points)` query returning the owning grain index for arbitrary points, periodic-aware - `adjust_X` (Lloyd step) wraps periodic centroids back into the box - CPU tests in tests/test_apds_periodic.py Co-Authored-By: Claude Opus 4.8 (1M context) <noreply@anthropic.com>
|
thanks, looks great, will possibly combine it or at least review against my own attempt at doing PBC (currently just privately). Claude Code should make it easy for both of us to review it if I do make some changes and turn into a bit of a "hybrid". My main thinking was to perhaps exploit LazyTensors here a bit and first go an argmin reduction along the "possible periodic images", it might make it much faster on GPUs to do it like that |
|
yes so I did a first quick test, it appears that your periodic BC approach is unfortunately not quite always correct (Claude Code initially hasn't realised that and tried to convince me that your approach is better and way faster computationally, so not surprised it didn't flag it up for you either). What you are doing is a nice quick approximation that works well if anisotropy is mild, but if a grain is really elongated, it can be very incorrect. And since the purpose of this library is optimising over the diagrams to make sure properties like grain size are prescribed, we cannot accept the approximation as valid. So this will need to be a hybrid of your code and my corrections. I can push directly to this branch if you are happy with it, @ligerzero-ai |
|
sure, that's fine by me |
|
just for curiosity's sake, what counts as mild anisotropy? just so I have a number in my head for future reference |
|
Each grain is defined by a weight, a centre point and an anisotropy matrix A, which also defines an equation of the ellipse/ellipsoid (see the ellipses in the logo of the library). A grain is isotropic if matrix A is identity, which also means that the ellipse is actually a circle / sphere. Your approach works exactly if every matrix is identity (but then the diagram is a power diagram, so all boundaries are straight lines) and becomes an approximation when matrix A is not identity. One way to measure whether we have mild anisotropy is to whether eigenvalues of the matrix A are close to 1. There is no precise way of saying at what point your approach gets it wrong (depends on number of grains, how close they are to the boundary etc), but it's very likely to be wrong when a grain is really elongated in one direction. |
Introduce an opt-in
periodicflag onapd_system. When enabled, every seed-to-point distance uses the minimum-image convention across the rectilinear domain, so the anisotropic power diagram becomes periodic and grains may wrap across the domain boundary.periodic=Falseconstructor flag (fully backward-compatible)_displacement(y, x)helper computing the minimum-image displacement as a KeOps LazyTensor viady - L * round(dy / L)D_ijcomputations (assemble_apd, OT_dual_function, check_optimality, adjust_X) route through_displacementgrain_of(points)query returning the owning grain index for arbitrary points, periodic-awareadjust_X(Lloyd step) wraps periodic centroids back into the box