Periodic microstructures and bug fixes

[rimoli]

This branch adds periodic microstructures in 2D and 3D, fixes a round of bugs found while going through the package, and gets the CI running again. The full list of changes is in CHANGELOG.rst under "Unreleased".

Periodic microstructures

• <periodic> in <domain> selects the periodic axes (True, or e.g. xy); the same option is an argument of cli.run, SeedList.position and PolyMesh.from_seeds.
• Seeds are placed with their periodic images, and the Laguerre tessellation is periodic: cells crossing a face are cut, and their pieces tile the domain.
• The triangular, tetrahedral and raster meshes have matching nodes on opposite faces. The quality and size settings (mesh_min_angle, mesh_max_volume, per-phase max_volume, mesh_max_edge_length) act as on non-periodic meshes.
• The Abaqus output has one node set per periodic face, in matching order, ready for periodic boundary conditions. The text files store the pairs of nodes and facets.
• A new setting, periodic_margin, keeps seeds from leaving thin slivers on the faces. It is off by default, and the examples use auto. edge_opt now also thickens thin pieces of cells at the faces and opens narrow corners there.
• New examples, each with a page in the docs: periodic_2D.xml, periodic_3D.xml, pbx_2D.xml, pbx_3D.xml, pbx_interface_2D.xml, pbx_interface_3D.xml and periodic_tiling.py.

Bug fixes and changes

The fixes are listed under "Fixed" in the changelog. They include reproducible seeding, cdf size distributions, TetGen volume constraints, the orientation of some ellipsoid breakdowns, cells clipped by circular and elliptical domains, several hangs, and the raster mesh and TriMesh.write outputs.

Some changes affect existing inputs, listed under "Changed":

• mesh_max_edge_length now also refines the grain boundaries of 3D meshes.
• The facets of triangular and tetrahedral meshes are sorted, so repeated runs give identical files. Triangle and TetGen list edges in an order that varies between runs.

Tests and docs

• 222 tests, against 13 on master. They include geometric checks of periodic tessellations and meshes against tiled non-periodic references.
• The documentation builds with Sphinx without warnings, including the new example pages.

CI

CI runs again and passes on Linux, macOS and Windows. The fixes on this branch:

• The import order (isort) in the files I touched.
• The test jobs could not start. The tox==3.14.0 pin in requirements.txt forced an old pluggy, with which pytest itself fails to import. This would also have hit master on its next push.
• pyvoro-mmalahe on PyPI bundles Voro++ 0.4.6, whose radical tessellation can return a cell uncut when the compiler fuses multiply-adds. On macOS, pip compiles pyvoro-mmalahe from source, and Apple clang fuses them by default on Apple silicon. The cells then overlap, and Triangle segfaulted on the result in the macOS jobs. For example, built that way, four discs of radius 0.1 at the quarter points of a unit square give four cells of area 1.0 instead of 0.25. Built without fused multiply-adds, as in the Linux wheels of pyvoro-mmalahe, the same case gives 0.25. Voro++ made these tests robust upstream in February 2026. As we discussed, I published a fork, pyvoro-rimoli (repository), which bundles the current Voro++ and has wheels for Linux, macOS and Windows, and MicroStructPy now depends on it. With it, the whole test suite passes, and our regression cases give bitwise the same meshes as the build we used during development.
• The Linux jobs were killed for a different reason. One of the new tests built a reference tessellation that needed about 10 GB of memory, and the runners ran out of it and were shut down. The reference now uses only the copies of the seeds near the domain, with the same values and a third of the memory.
• Read the Docs used Python 3.8, which pyvoro-rimoli does not support. It now uses Python 3.10 and installs requirements.txt, like the documentation check.

One check is still red. Read the Docs now installs everything, but the build runs every example, and with the new periodic examples it goes past the 15-minute limit of the free plan. The documentation checks on GitHub build the same pages and pass.

Trying it

pyvoro-rimoli installs the same pyvoro module as pyvoro-mmalahe, so uninstall the old one first:

pip uninstall pyvoro-mmalahe
pip install -e .
microstructpy --demo=periodic_2D.xml

🤖 Generated with Claude Code

[kip-hart]

I'm getting ready to post feedback on this PR. Nothing major. Currently we're over the docs builder's (readthedocs) limit of 15 minutes. I'm testing if tweaking the parameters for the basalt circle example can get us under that limit. That testing is in #115 .

[kip-hart]

I pulled the branch and ran it instead of reading the diff alone. Python 3.12.3,
a fresh venv built from this branch's requirements.txt, package installed with
pip install -e ..

The suite passes. 217 passed in 169s, no failures, which matches what the
description claims. The jump in coverage is in real tests rather than in
parametrization: master has 14 test functions across tests/, this branch has
191, and those 191 collect as the 217 cases.

CI claims that check out

The tox==3.14.0 diagnosis is right, and I reproduced it in the order the
workflow actually installs things. Installing flake8 pytest pytest-cov coveralls first and then requirements.txt pulls pluggy back down to 0.13.1,
and pytest then fails to import at all:

ImportError: cannot import name 'HookimplOpts' from 'pluggy'

One correction for the description. It says the old pluggy is one "that the
current pytest-cov cannot load". It is pytest itself that cannot load, not
pytest-cov. Same root cause, same fix.

The pyvoro fork is what it says it is. I diffed both sdists from PyPI.
pyvoro/__init__.py, pyvoro/vpp.cpp and pyvoro/vpp.h are byte identical to
pyvoro-mmalahe 1.3.4, so nothing changed in the Python facing layer. The bundled
Voro++ really is updated, 2396 lines of difference across the 27 sources in
src/.

CI's blocking lint passes. flake8 src tests setup.py --exclude=__init__.py --select=E9,F63,F7,F82 returns 0. The five remaining cosmetic warnings are all
in docs_banner.py and logo.py, which this branch does not touch.

What I found

Four defects I confirmed by running code, in the order I would fix them:

1. The 2D orientation error statistics are wrong by up to 360 degrees, because
   the angles are subtracted without wrapping. A fit whose true worst error is
   0.8 degrees reports mae=359.50 and R^2=-2.997. Details in the inline
   comment on verification.py.
2. <periodic>0</periodic> raises TypeError instead of turning periodicity
   off. periodic_axes has no integer branch.
3. Six geometry classes lost hashability, because __eq__ was added without
   __hash__.
4. One whole domain Triangle build per non converging input is computed and
   thrown away.

One more that is not a defect but is a cost, in the inline comment on
calc_rtol. The new fit lowers the tolerance on high cv inputs, which rejects
more seed placements. On basalt_circle.xml, 1300 seeds at cv 4.27 and no
explicit rtol, positioning goes from 6.1 s to 35.8 s with everything but the
tolerance held fixed. That accounts for 29.7 s of the 32.5 s that example got
slower overall. Every input with a long tailed size distribution and a default
rtol pays it.

Two more that are real but narrow. The two np.allclose(..., atol=tol) calls in
the periodic 2D boundary code leave rtol=1e-5 in place, which defeats the
tol = 1e-9 * scale the code sets deliberately. I worked out what it takes to
bite: a false match needs 1e-5 * |mins| >= lengths, so an offset to size ratio
of 1e5 or more. At exactly that ratio the periodic build does fail, with
RuntimeError: Triangle did not keep the input points., but I could not
attribute that failure to the tolerance leak rather than to double precision at
those coordinates. I am flagging the leak as latent rather than as a
demonstrated wrong answer.

Two things I want your input on

A note on the pyvoro reproducer, which does not change the conclusion. Four
discs of radius 0.1 at the quarter points of a unit square give the correct
[0.25, 0.25, 0.25, 0.25] under pyvoro-mmalahe here, not the four cells of area
1.0 the description reports. I ran a differential fuzz past that too, 200 2D
cases per periodicity setting at unequal radii plus 120 3D cases, identical
seeds against both builds, and the outputs are bitwise identical.

That is a statement about this machine and not a challenge to the swap. Two
things explain it. pyvoro-mmalahe publishes wheels only for manylinux x86_64,
so macOS and Windows have always compiled it from the sdist, with clang and
MSVC. And the upstream fix is a numerical robustness change: rad_option.hh
drops sqrt from the two radical cutting predicates and computes r_mul as a
difference of squares, while config.hh replaces the 1e-11 tolerance and the
whole marginal case buffer with 10*DBL_EPSILON and an extra search stack. In
strict IEEE double the old and new predicates agree on four million random
draws I threw at them. They diverge once a compiler contracts expressions or
keeps intermediates in wider precision, which is exactly how this ends up
showing on one platform and not another.

So the swap looks justified to me. The ask is smaller: say in the description
which platform the four disc case was observed on, since it reads as a flat
statement about the library and someone on Linux will try it and see 0.25.

The other item is calc_rtol. See the inline comment for the numbers.

On the new dependency

pyvoro-rimoli 1.4.0 is a single release uploaded 2026-09-30, with no
attestations or Trusted Publishing provenance in its PyPI metadata. Both
setup.py and requirements.txt point at it, so every pip install microstructpy and every CI job pulls it. The author_email in the metadata is
still Joe Jordan's, inherited from upstream, which is worth fixing on your side
so the package reports its actual maintainer.

The wheels cover cp39 through cp314 and requires_python is >=3.9, so the
3.9 entry in the test matrix is fine.

Documentation

The docs hold up. The new example pages follow the same section structure as
the existing ones, the docstring convention is unchanged, and coverage of it
went up rather than down: Returns: goes from 44 to 95 and Raises: from 1 to
10 across the library.

Three small things, all in inline comments. The first Added bullet of the
changelog runs 30 lines against a median of 4 and buries the gmsh limitation at
the end of it. The added prose joins clauses with semicolons 54 times where the
existing docs never do. And labelled appears four times in a repo that is
otherwise consistently US spelled.

The remaining red check is the Read the Docs timeout, and it is now handled
separately in #115, so nothing is needed here for it. The short version is that
the build runs every example and cannot finish inside the fifteen minutes Read
the Docs allows, so #115 moves that work to the documentation workflow and has
Read the Docs restore the figures instead. That took it from a 901 s timeout to
96 s.

Two things came out of that work which belong to this branch rather than to
#115. The numpy churn in the inline comment on .readthedocs.yaml, and an
invalid XML declaration in elliptical_grains.xml, which declared
version="0.5" where the specification allows only 1.x. The second is older
than this branch and I fixed it in #115, since it was failing the documentation
check here intermittently: Expat 2.6.1 accepts it and newer Expat rejects it, so
it passed locally and failed on some runners.

[kip-hart]

I went through the fork before commenting here. The Python facing layer is
byte identical to pyvoro-mmalahe 1.3.4, and the Voro++ update is real and
does what the changelog says, so the substance of this is fine by me.

For the record of why, since the four disc case in the description does not
fail on Linux: the fix is in rad_option.hh, which drops sqrt from
r_ctest and r_scale_check and computes r_mul as (r-mx)*(r+mx) instead
of r*r - mx*mx, and in config.hh, which replaces the hard coded 1e-11
tolerance and the marginal case buffer with 10*DBL_EPSILON and an extra
search stack. cell.cc carries 1282 substantive changed lines on top of
that. The old and new predicates agree in strict IEEE double across four
million draws, so the divergence needs a compiler that contracts expressions
or widens intermediates. That fits a failure on the macOS and Windows
runners, which is where this always built from source, given pyvoro-mmalahe
ships wheels for manylinux x86_64 only.

Three things I would still like.

Name the platform in the description for the four disc case. As written it
reads as a property of the library, and anyone on Linux will try it and get
0.25.

Keep the bitwise regression check you mention. Thirteen years of upstream
changes to the cell cutting core is a lot of moving parts under our meshes,
and a stored reference mesh is the cheapest guard against the next Voro++
bump moving results quietly.

Turn on Trusted Publishing in the fork's release workflow so the wheels carry
PyPI attestations. It is cheap to set up and this is now a dependency of a
published package. The author_email in the metadata is also still
joe.jordan@imperial.ac.uk, inherited from upstream, so the package does not
name its actual maintainer.

[rimoli]

Thanks for going through the fork. Three answers.

Platform: you were right, and I confirmed your explanation. The 0.4.6 sources built with Apple clang's defaults on Apple silicon give the four cells of area 1.0. The same sources built with -ffp-contract=off give 0.25 and pass all the tests of the fork. The PR description, the changelog and the fork's README now say this.

Regression guard: pyvoro-rimoli 1.4.1 adds reference tessellations to its tests: the volumes and neighbors of the cells of four fixed configurations, 2D and 3D, bounded and periodic. They are compared to a relative 1e-10, since a bitwise comparison would not survive different compilers, which is the very issue here. The old build with fused multiply-adds fails them, and the same sources built without them pass. The bitwise regression cases of MicroStructPy stay too, and I reran them for this round.

Trusted publishing: 1.4.0 was already published with it, with attestations. PyPI's integrity API returns the provenance of each file, for example https://pypi.org/integrity/pyvoro-rimoli/1.4.0/pyvoro_rimoli-1.4.0.tar.gz/provenance, with publisher rimoli/pyvoro, workflow wheels.yml and environment pypi. The JSON metadata does not show it. The maintainer email is in 1.4.1.

[kip-hart]

Thanks for confirming that with trusted publishing. My comment was an AI hallucination that I didn't check myself.

[kip-hart]

I think the last thing is to bump the version pin for pyvoro-rimoli to 1.4.1

[rimoli]

Thanks for the careful review, Kip, and for running everything rather than reading the diff. All of it is addressed in eight commits on develop, with a reply on each inline comment:

• df1095c: <periodic> 0 </periodic> and 1 work.
• c428e83: the geometry classes are hashable again, and Circle and Sphere now too.
• fe7d152: rtol=0 in the tolerance checks of the 2D periodic mesher (your two sites, plus a third).
• 35b7478: no unused Triangle build after the last pass.
• 29eb126: the 2D orientation statistics compare angles modulo the symmetry of the shape.
• 543684e: numpy<2.0 in docs/requirements.txt.
• f8e2e41: the changelog entry split up, sentences instead of semicolons, "labeled", and your pytest correction.
• 3697c37: rtol='fit' back to master's coefficients, which match Fig. 10 of the paper. Eqs. (14) and (15) as printed do not (details in the inline reply).

Checks: with all the code fixes but the last, the six bitwise regression cases are identical. The last one changes only the case that uses rtol='fit', as expected. The full suite (222 tests) and a nitpicky docs build pass locally.

On the pyvoro side, pyvoro-rimoli 1.4.1 is out, with the maintainer email and the reference tessellations. The four-disc case is Apple silicon, and your explanation with fused multiply-adds is right. The description now says so. Trusted publishing was on from the start, as explained in the inline reply.

One note for #115: these commits change the library and the periodic_3D.xml example, so the example digest changes, and the docs cache will need publishing for the new head.

[kip-hart]

max_volume is ignored in 2D whenever the caller does not pass phases.

These defaults set the key to infinity, so phase.get('max_volume', max_volume) at line 1485 finds it present and never falls back to the
max_volume argument. Every region is then given an area constraint of
infinity. On master that did not matter, because the global constraint was
passed to Triangle underneath it:

master   meshpy.triangle.build(info, volume_constraints=True,
                               max_volume=max_volume, ...)
this PR  meshpy.triangle.build(info, volume_constraints=True,
                               max_volume=None, ...)

_build_2d dropped the global in 24fa757 and relies on the per-region
values alone. That is right for the periodic mesher, which needs different
constraints for the cells and their copies, and a single global value would
override them. The non-periodic path goes through the same helper, and for
it the infinity in these defaults became load bearing.

Measured on this branch, same polygonal mesh, only phases differing:

phases=None, max_volume=0.4      elements=208   largest element 2.288
phases=[{max_volume: 0.4}]       elements=252   largest element 0.385

The cap is 0.4 and the largest element is 2.288, so it is not being applied
at all. With the suggestion below the first line gives 252 elements and a
largest element of 0.385, and the full suite still passes at 222.

This is visible in the documentation. docs_banner.py calls

TriMesh.from_polymesh(pmesh, min_angle=12, max_edge_length=0.2,
                      max_volume=0.4)

with no phases, and the third panel of banner.png now has a coarse
ellipse where the published one is finely triangulated.

Dropping the key lets the argument through, which is what every reader of a
.get with a default expects:

Suggested change

	return [{'material_type': 'solid', 'max_volume': float('inf')}
	for _ in range(n_phases)]
	return [{'material_type': 'solid'} for _ in range(n_phases)]

Two things worth deciding rather than inheriting.

The same .get appears at lines 2036 and 2892, so this changes the fallback
in _facet_sizes and in the periodic region building as well. Only the
phases=None path is affected, and the new behavior is the intended one in
all three, but it is more than one call site.

An explicit {'max_volume': float('inf')} means unconstrained on this
branch, where on master the global still applied to it. The suggestion keeps
that, which I think is the better reading. It is a change from master either
way and belongs in the changelog.

The reason 222 tests did not catch this is that nothing exercises
phases=None with a finite max_volume. docs_banner.py is the only
caller that takes the path, and nothing asserts on what it draws. A test
that meshes one cell that way and checks the largest element against the cap
would close it.

[kip-hart]

Sorry I didn't comment on this before, I only noticed the drift after seeing the new docs site.