Smooth a MaterialGrid object's own boundary, not just its level set - #3279
Smooth a MaterialGrid object's own boundary, not just its level set#3279jin-castle wants to merge 5 commits into
Conversation
|
This would be a nice feature to have since the workaround of defining a 3d material grid to smooth out the normal boundaries of a 2d design region with finite thickness could cause you to run out of memory when doing parallel meep calculations since you would have to load a large 3d material grid onto each process |
| matgrid_volavg *mgva = (matgrid_volavg *)mgva_; | ||
| get_uproj_w(mgva, x[0], u_proj, w); | ||
| return w * ((1 - u_proj) / mgva->eps1 + u_proj / mgva->eps2); | ||
| return w / ((1 - u_proj) * mgva->eps1 + u_proj * mgva->eps2); |
There was a problem hiding this comment.
This isn't related to the core PR goal (smooth MG boundaries), right? The native cpp MG smoothing code indeed needs to be revamped, but I think we should do that separately.
There was a problem hiding this comment.
Right on both counts — separate defect, separate PR. I've split it.
Stays here (124554c9, 4967bd35, and the two follow-ups): the
get_front_object() / eff_chi1inv_matrix() change that lets a grid object's
own boundary reach the two-material average, its TestObjectBoundarySmoothing
cases, and the symmetry-invariance and anisotropic-media coverage.
Moves out (ccbb655a, 28cb1512): the matgrid_inveps_func integrand and
the "skip the quadrature when no interface crosses the voxel" shortcut, which
only becomes safe once that integrand is continuous. I'll open them as a
follow-up against master.
One test crosses the line and I want to flag it rather than have it look like a
silent revert: test_interior_level_set_is_still_smoothed goes back to its
original form in this PR, and the rewrite of it — the one that distinguishes a
projected (β=∞) ramp, which keeps a real interface in the interior, from an
unprojected (β=0) ramp, which does not — travels with the follow-up, since it
is asserting the behaviour that only exists there.
The split is exact: the two branches applied in sequence reproduce this PR's
tree byte-for-byte, so nothing is being quietly dropped or rewritten in the
process. On a clean MPI-enabled double-precision build of the reduced branch,
test_subpixel_3d.py is 11/11 and test_material_grid.py 3/3 — the same
subpixel count as before the split, because the one test whose assertion
depended on the follow-up reverts along with it.
For the follow-up PR's benefit rather than this one's, the short version of why
that one line is not cosmetic: matgrid_inveps_func integrates
(1-u)/eps1 + u/eps2, the inverse of a binary mixture, which is the right
quantity only at u ∈ {0,1}. The medium the surrounding code assembles at an
interior node is the graded one, (1-u)*eps1 + u*eps2. The gap between them
keeps the anisotropic correction finite as |∇u| → 0, so the assembled epsilon
is discontinuous in the design weights at the |∇u| < 1e-8 gate — which is why
it surfaces as a gradient error rather than a small bias. For weights projected
upstream (SSP at β=∞, where 64 % of nodes sit exactly at 0 or 1) the directional
FD error goes 6.6 % → 0.094 % at dp=1e-3 and 47.5 % → 0.003 % at dp=3e-3;
conic-filtered (smooth) weights go 1.8 % → 0.033 %; the do_averaging=False
control is 0.000 % before and after. That belongs in its own thread with its own
tests, and it sounds like it overlaps whatever revamp you have in mind, so I'd
rather it be discussed there than bolted onto this one.
|
|
||
| double fill = box_overlap_with_object(pixel, *o, tol, maxeval); | ||
|
|
||
| if (front_is_matgrid) { |
There was a problem hiding this comment.
I'm not sure this part will properly get reflected in the
get_front_object() rejected any pixel touching a variable material, so a material-grid object's boundary never reached the analytic two-material average; every such pixel fell through to fallback_chi1inv_row(). That fallback takes its normal from matgrid_grad(), i.e. from grad(u) alone, so it only sees interfaces that live in the design level set. An (n,n,1) grid has no z structure, and the top and bottom faces of a 3d design slab are not in the level set at all -- they are the geometric object's boundary. Those faces therefore stayed staircased regardless of do_averaging. They do not exist in 2d, which is why only 3d runs paid for it. - Narrow the rejection to what actually breaks the average. mat_behind must still be constant over the pixel, so a variable material there is rejected as before. As the material of the *front* object a material grid is admissible: the object's boundary is an ordinary level set, and the caller localizes the grid before averaging. A user material function keeps the old treatment -- its fallback derives the normal from normal_vector() and integrates chi1p1(), both of which sample the actual material and therefore already see the object boundary. Only the material-grid fallback is blind to it. - Localize the grid inside the pixel before averaging (inside_sample_point()). The pixel center can sit on the far side of the boundary, where a material grid is evaluated on its mirrored extension, so step to the centroid of the part that is inside: (1-fill)*h/2 along the inward normal for a planar cut. libctl does not fix the normal's orientation, so try both and keep the candidate that is actually inside the object. - Average only the boundary. Pixels with fill >= 1 - 1e-6 are handed back to the fallback; averaging them against a filling fraction of 1 would return the unsmoothed material and silently disable level-set smoothing across the whole design-region interior. fill <= 1e-6 takes the no-average path. - Defer the metal test for a front-side grid until after localization. is_metal() reads the grid's `medium`, which before get_material_pt() still holds the previous pixel's interpolation. mat_behind is never a grid, so its test stays where it was. test_subpixel_3d.py gains TestObjectBoundarySmoothing. Filling `weights` with a constant removes the in-plane level set, leaving the block's six faces as the only material interfaces, so the grid must be smoothed exactly as an ordinary Medium block of the same epsilon -- the z profiles agree to 1e-8. A second test guards that assertion against passing for the wrong reason: before this change the do_averaging on/off difference was identically zero. A third drives a ramped grid to confirm the interior level set is still smoothed. This does not address the grad(u) ~ 0 regime of NanoComp#2757. Interior pixels still go through the fallback's |grad(u)| > 1e-8 gate, which applies the mixing rule without testing whether an interface lies inside the pixel.
Three regimes the smoothing path had no coverage for: - symmetries must not change the smoothed epsilon assembly (structure only -- the adjoint gradient under symmetries has its own independent, pre-existing defect). Measured invariant to 3.6e-15 with and without Mirror(Y), at beta = inf and beta = 0. - a uniform grid of anisotropic media gets the same analytic boundary average as a plain block of the interpolated anisotropic epsilon (identical to the last bit -- same code path), and the on/off guard shows the average is actually engaging. - in the interior, the level-set fallback deliberately declines anisotropic pixels (its 1d normal average is scalar), so do_averaging must be a finite no-op there rather than a half-applied average.
test_gradient_matches_fd_with_smoothing holds fd/adjoint to 1% and the adjoint gradient currently depends on the MPI chunk division -- the pre-existing defect NanoComp#3274 fixes. This branch moved the smoothing operating point, and the chunk error that used to hide inside the tolerance now reads 1.3% at np=2 and 7.3% at np=3 (it is np-dependent, the signature of NanoComp#3274, and passes at every rank count with NanoComp#3274 merged in -- verified on a build of this branch merged with its head). Keep the serial assertion at 1%; bound the parallel one at 0.1, still far below the 70-100% disagreement the test guards against, and tighten it back when NanoComp#3274 is in.
6a1c35c to
bf794b9
Compare
jin-castle
changed the title
4 participants
do_averaging=Truepromises subpixel smoothing for aMaterialGrid, but amaterial grid's own object boundary never gets it. In 3d the top and bottom
faces of a design slab — the interfaces that set the vertical confinement of
every strip-waveguide design — stay staircased at every resolution, no matter
what
do_averagingis set to.One commit for the defect, a second that adds symmetry-invariance and
anisotropic-media coverage this path never had.
The material-grid fallback is the only smoothing path that cannot see an object boundary
There are two paths.
eff_chi1inv_matrix()callsget_front_object()and, ifthat succeeds, averages two materials across the object's boundary using its
normal and the pixel's filling fraction. If it fails,
fallback_chi1inv_row()integrates along a level set instead.
get_front_object()rejected any pixel touching a variable material:So a material-grid object's own boundary never reached the first path. And the
fallback takes its normal from
∇ualone for a grid while a user materialfunction gets it from
normal_vector()and integrateschi1p1()— both ofwhich sample the actual material and therefore already see the object boundary.
Only the material-grid branch is blind to it.
For an
(n, n, 1)grid∇uhas no z component at all. The slab's z faces arenot in the design level set; they are the object's boundary. Nothing in either
path was smoothing them. Those faces do not exist in 2d, and every adjoint test
upstream is 2d or cylindrical, which is why this survived.
Why the front material may be a grid and the back material may not
The average needs one material that is constant over the pixel to act as the
background, and one boundary with a normal and a filling fraction.
mat_behindmust therefore stay constant — a variable material there is still rejected.
mat_frontneed not be: the boundary being averaged is the object's own levelset, and the grid can be evaluated at a point before averaging.
Which point matters. The pixel center can sit on the far side of the boundary,
and that is exactly where a
MaterialGridis evaluated on its mirroredextension. For a planar cut of a pixel of size
hwith inside fractionfill,the centroid of the inside part lies
(1-fill)*h/2from the center along theinward normal. libctl does not fix the normal's orientation, so both candidates
are tried and the one
point_in_fixed_objectp()accepts is kept.Interior pixels must go back to the fallback. A pixel with
fill = 1has noboundary in it; averaging it would return the unsmoothed material and silently
disable level-set smoothing across the whole design region. Those are handed
back.
Changes
get_front_object()to what actually breaks theaverage, and move it after the front/behind assignment so it can distinguish
the two. A material grid is admissible as the front object's own material
(not as
default_material, which has no boundary to average) whendo_averagingis set;front_is_matgridtells the caller to localize it.inside_sample_point()— evaluate the grid at the centroid of the insidepart of the pixel, per the derivation above.
fill >= 1 - 1e-6returnsfallback = true;fill <= 1e-6takes the no-average path.is_metal()reads the grid'smedium, which beforeget_material_pt()still holds the previous pixel's interpolation.
mat_behindis never a grid,so its test stays where it was.
Tests
TestObjectBoundarySmoothinginpython/tests/test_subpixel_3d.py. Fillingweightswith a constant removes the in-plane level set entirely, leaving theblock's six faces as the only material interfaces — so whatever smoothing the
grid gets there must be the smoothing an ordinary
Mediumblock of the sameepsilon gets.
1e-8;do_averagingon/off difference exceeds 0.5, which guards the firstassertion against passing because both sides are staircased — before this
commit that difference was identically zero;
that interior pixels are not swallowed by the
fill = 1boundary path.A second commit adds coverage for two regimes this path never had tests for:
symmetriesinvariance — the smoothed epsilon assembly is identical withand without
mp.Mirror(mp.Y)for a mirror-symmetric design (measured3.6e-15, at
beta=infandbeta=0). Structure only: the adjoint gradientunder
symmetrieshas its own independent, pre-existing defect that nothinghere touches.
analytic boundary average as a plain block of the interpolated anisotropic
epsilon (identical to the last bit — same code path), with an on/off guard
showing the average engages; and in the interior, where the level-set
fallback deliberately declines anisotropic pixels (its 1d normal average is
scalar),
do_averagingis pinned as a finite no-op rather than ahalf-applied average.
One CI note: the MPI jobs initially tripped on
test_gradient_matches_fd_with_smoothing— an existing test from #3263 thatholds fd/adjoint to 1%. The excess is the chunk-dependence of the adjoint
gradient that #3274 fixes. After splitting this boundary-only change from the
graded-harmonic work, the two-rank case differs by 10.29% in both precisions;
it vanishes on this branch merged with #3274's head. Until #3274 lands, the
parallel assertion is bounded at 0.11 — still far below the 70–100%
disagreement the test guards against — and the serial assertion keeps its 1%.
Also checked against the pending #3277, since both touch
meepgeom.cpp: themerge is conflict-free and both sides' suites pass on the combined build,
including #3277's
test_adjoint_symmetric_grids(3/3) andtest_adjoint_adjacent_grids(5/5). One composition result worth recording: aMaterialGridoverlapped with a mirrored copy of itself (U_MEAN, the #1984pattern) assembles the same smoothed epsilon as a single grid carrying the
explicitly averaged weights, to 2.7e-15.
What would falsify this
If the uniform-grid z profile does not match the equivalent block, the
localization point or the filling fraction is wrong. If the interior ramp's
on/off difference collapses, interior pixels are being averaged against
fill = 1and level-set smoothing is gone. Both are asserted, so the tests arethe falsifier.
Not addressed
#1965 (edge artifacts that vanish when the grid is
default_materialrather than a
Block) may have the same root as this commit, but I have notconfirmed it.
cc @smartalecH @stevengj