Add support for dft_force with mirror symmetry

[oskooi]

Fixes #129.

Context

#129 reports field instability when combining DFT force monitors with symmetry objects. The root cause is that fields::add_dft_force() in src/stress.cpp unconditionally calls S.reduce(where_) to apply symmetry reduction to the force region volumes. In contrast, fields::add_dft_flux() in src/dft.cpp has a use_symmetry parameter (defaulting to true) that allows bypassing S.reduce().

The S.reduce() operation is problematic for force calculations because the stress tensor involves pairing different field components (E in force direction vs. E in normal direction) as off-diagonal terms. When symmetry transforms the volume and modifies the component weights, the paired DFT chunk lists (offdiag1, offdiag2) can end up with mismatched structures, producing incorrect/unstable force values.

Fix

Step 1: Add use_symmetry parameter to add_dft_force in C++ (3 files)

src/meep.hpp (lines 2165-2174): Add bool use_symmetry = true parameter to all three add_dft_force overloads, matching the pattern of add_dft_flux (lines 2045-2057).

src/stress.cpp (lines 153-157): Accept the new use_symmetry parameter and conditionally apply S.reduce(): volume_list *where = use_symmetry ? S.reduce(where_): new volume_list(where_); This mirrors dft.cpp line 607.

Step 2: Pass use_symmetry=False from Python when symmetries are present (1 file)

python/simulation.py (lines 3399-3404): Modify _add_force to pass use_symmetry=False to add_dft_force when the simulation has symmetry objects. This is done by passing an extra argument through _add_fluxish_stuff, which already supports *args forwarding (line 3818: add_dft_stuff(vol_list, freq, decimation_factor, *args)).

Step 3: Add unit test (1 file)

1. Runs a force calculation on a geometry with mirror symmetry using symmetries=[mp.Mirror(mp.Y)].
2. Runs the same geometry without symmetry.
3. Asserts both produce the same force values (within tolerance).

[codecov-commenter]

⚠️ Please install the to ensure uploads and comments are reliably processed by Codecov.

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 73.91%. Comparing base (f29a8c7) to head (d7cb31a).
⚠️ Report is 111 commits behind head on master.
❗ Your organization needs to install the Codecov GitHub app to enable full functionality.

Additional details and impacted files

@@            Coverage Diff             @@
##           master    #3161      +/-   ##
==========================================
+ Coverage   73.81%   73.91%   +0.09%     
==========================================
  Files          18       18              
  Lines        5423     5455      +32     
==========================================
+ Hits         4003     4032      +29     
- Misses       1420     1423       +3     

Files with missing lines	Coverage Δ
python/simulation.py	77.57% <100.00%> (+0.16%)	⬆️

... and 1 file with indirect coverage changes

🚀 New features to boost your workflow:

• ❄️ Test Analytics: Detect flaky tests, report on failures, and find test suite problems.

[stevengj]

This "fix" seems wrong. The volume_list *where = S.reduce(where_); line in stress.cpp knows about the force direction that is being computed (via where->c), and hence takes this into account in the symmetry reduction, exactly as we do for flux. Basically there is a where->weight that is update as necessary, and is passed as the extra_weight parameter in add_dft, so that it get included (once) in each stress term.

It could of course be that there is some other bug in how this procedure is being carried out, but the solution is not to remove it entirely, since the symmetry reduction is an important optimization.

[stevengj]

It is true, however, that the diagonal and off-diagonal stress-tensor terms should be transformed separately. And the problem is essentially that the where->c transformation is only looking at one of the components, and not both.

(Diagonal terms basically don't pick up any transformation, whereas off-diagonal terms pick up two transformations.)

The right fix here is maybe to add a second component to the volume_list, defaulting to Dielectric (or something that transforms as a scalar), which we can then set appropriately for these transformations.

(For the diagonal terms, we can just set where->c to Dielectric as well.)

[stevengj]

An even easier fix is to optionally let S.reduce look at the normal direction in addition to where->c in order to decide the phase.