feat: downsample frequencies in mode monitors (FXC-3351)

[dbochkov-flexcompute]

Addresses #2530.

Notebooks re-run examples:

• photonics: https://github.com/flexcompute/tidy3d-notebooks/blob/d6109a7182701e0e487e1cf4d35434f5063c7bca/RingResonator.ipynb (1001 mode solves -> 21 mode solves with no visible change)
• RF: https://github.com/flexcompute/tidy3d-notebooks/blob/d6109a7182701e0e487e1cf4d35434f5063c7bca/DifferentialStripline.ipynb

Still not sure about making it default in WavePorts

Greptile Overview

Greptile Summary

This PR implements frequency downsampling/interpolation for mode monitors and mode solvers to reduce computational cost in broadband simulations. The implementation introduces a new ModeInterpSpec class with three interpolation methods (linear, cubic, Chebyshev polynomial), allowing modes to be computed at fewer strategically-chosen frequencies and interpolated to the full frequency set.

Key Changes:

• New ModeInterpSpec class with UniformSampling, ChebSampling, and CustomSampling strategies for flexible frequency selection
• Integration throughout the mode computation pipeline: ModeSpec, AbstractModeMonitor, ModeSolver, and ModeSolverData
• Comprehensive interpolation implementation in ModeSolverData.interp() supporting all three methods with proper handling of complex-valued data
• Requires mode tracking to be enabled (mode_spec.sort_spec.track_freq) to maintain consistent mode ordering across frequencies
• Extensive test coverage with 965 lines of well-structured tests
• Extrapolation warnings when interpolating outside the computed frequency range

Notable Implementation Details:

• Chebyshev nodes use barycentric formula for optimal polynomial interpolation
• Complex data is interpolated by treating real and imaginary parts independently
• reduce_data flag allows storing only sampled frequencies to minimize storage costs
• Grid correction factors can be recalculated or interpolated based on user preference

Confidence Score: 4/5

• This PR is generally safe to merge with minor issues that should be addressed
• Score reflects well-designed core functionality with comprehensive testing, but deducted one point for: (1) missing upper bound validation on num_points in ModeInterpSpec that could exceed total frequencies, (2) f-string usage in warning message violating style guide rule 9b45957a-2d0d-4c24-a5d3-eef6f721cbfd, (3) potential edge case with single-frequency inputs not explicitly handled, and (4) extrapolation tolerance using multiplicative factor may behave unexpectedly near zero frequencies
• Pay closest attention to tidy3d/components/mode_spec.py (validation gaps) and tidy3d/components/monitor.py (style guide violation)

Important Files Changed

File Analysis

Filename	Score	Overview
tidy3d/components/mode_spec.py	4/5	Introduces ModeInterpSpec with three interpolation methods and sampling strategies; validation ensures num_points meets method requirements but lacks upper bound check and single-frequency edge case handling
tidy3d/components/monitor.py	4/5	Adds interp_spec field to AbstractModeMonitor with tracking requirement validation; f-string usage in warning violates style guide
tidy3d/components/mode/mode_solver.py	5/5	Implements core interpolation logic in ModeSolver; debug print statements have been removed (no longer present in latest commit)
tidy3d/components/data/monitor_data.py	4/5	Adds comprehensive interp() method to ModeSolverData supporting all interpolation methods; extrapolation tolerance uses multiplicative factor which may be problematic near zero frequencies
tests/test_components/test_mode_interp.py	5/5	Comprehensive test coverage (965 lines) for all interpolation functionality with well-structured test cases

Sequence Diagram

sequenceDiagram
    participant User
    participant ModeSolver
    participant ModeInterpSpec
    participant ModeSolverData
    participant compute_modes

    User->>ModeSolver: Create with interp_spec
    ModeSolver->>ModeInterpSpec: Validate num_points
    ModeInterpSpec-->>ModeSolver: Validation complete
    
    User->>ModeSolver: Call data_raw property
    
    alt interp_spec is None
        ModeSolver->>compute_modes: Solve at all frequencies
        compute_modes-->>ModeSolver: Return raw data
        ModeSolver-->>User: ModeSolverData (all freqs)
    else interp_spec provided
        ModeSolver->>ModeInterpSpec: sampling_points(freqs)
        
        alt method is linear/cubic
            ModeInterpSpec->>ModeInterpSpec: Generate uniform points
        else method is cheb
            ModeInterpSpec->>ModeInterpSpec: Generate Chebyshev nodes
        end
        
        ModeInterpSpec-->>ModeSolver: Return freqs_reduced
        ModeSolver->>ModeSolver: Create mode_solver_reduced
        ModeSolver->>compute_modes: Solve at reduced freqs
        compute_modes-->>ModeSolver: Return data_reduced
        
        ModeSolver->>ModeSolverData: interp(freqs, method)
        
        alt method is cheb
            ModeSolverData->>ModeSolverData: Barycentric interpolation
        else method is linear/cubic
            ModeSolverData->>ModeSolverData: xarray interpolation
        end
        
        ModeSolverData-->>ModeSolver: Interpolated data
        ModeSolver-->>User: ModeSolverData (all freqs)
    end

Loading

[greptile-apps]

_{7 files reviewed, 9 comments}

_{Edit Code Review Agent Settings | Greptile}

[github-actions]

Diff Coverage

Diff: origin/develop...HEAD, staged and unstaged changes

• tidy3d/init.py (100%)
• tidy3d/components/data/dataset.py (97.1%): Missing lines 173
• tidy3d/components/data/monitor_data.py (100%)
• tidy3d/components/microwave/data/dataset.py (100%)
• tidy3d/components/microwave/data/monitor_data.py (100%)
• tidy3d/components/mode/mode_solver.py (89.1%): Missing lines 481,622,625,989,1979
• tidy3d/components/mode/simulation.py (100%)
• tidy3d/components/mode_spec.py (97.2%): Missing lines 201,298,685
• tidy3d/components/monitor.py (90.9%): Missing lines 826
• tidy3d/components/validators.py (100%)

Summary

• Total: 297 lines
• Missing: 10 lines
• Coverage: 96%

tidy3d/components/data/dataset.py

Lines 169-177

  169         num_freqs, num_modes = sort_inds_2d.shape
  170         modify_data = {}
  171         for key, data in self.data_arrs.items():
  172             if "mode_index" not in data.dims or "f" not in data.dims:
! 173                 continue
  174             dims_orig = data.dims
  175             f_coord = data.coords["f"]
  176             slices = []
  177             for ifreq in range(num_freqs):

---

tidy3d/components/mode/mode_solver.py

Lines 477-485

  477         return self.data
  478 
  479     def _freqs_for_group_index(self, freqs: FreqArray) -> FreqArray:
  480         """Get frequencies used to compute group index."""
! 481         return self.mode_spec._freqs_for_group_index(freqs=self.freqs)
  482 
  483     @cached_property
  484     def _sampling_freqs(self) -> FreqArray:
  485         """Get frequencies used to compute group index and interpolation."""

---

Lines 618-629

  618         # TODO: At a later time, we should ensure that `eps_spec` is automatically returned to
  619         # to compute the backward propagation mode solution using a mode solver
  620         # with direction "-".
  621         eps_spec = []
! 622         for _ in solver.freqs:
  623             eps_spec.append("tensorial_complex")
  624         # finite grid corrections
! 625         grid_factors, relative_grid_distances = solver._grid_correction(
  626             simulation=solver.simulation,
  627             plane=solver.plane,
  628             mode_spec=solver.mode_spec,
  629             n_complex=n_complex,

---

Lines 985-993

  985                 # Extract coordinate values from one of the six field components
  986                 xyz_coords = solver.grid_snapped[field_name].to_list
  987                 x, y, z = (coord.copy() for coord in xyz_coords)
  988 
! 989             f = np.atleast_1d(self._sampling_freqs)
  990             mode_index = np.arange(self.mode_spec.num_modes)
  991 
  992             # Initialize output arrays
  993             shape = (x.size, y.size, z.size, len(f), mode_index.size)

---

Lines 1975-1983

  1975         """Check if there are media with a complex-valued epsilon in the plane of the mode.
  1976         A separate check is done inside the solver, which looks at the actual
  1977         eps and mu and uses a tolerance to determine whether to use real or complex fields, so
  1978         the actual behavior may differ from what's predicted by this property."""
! 1979         check_freqs = np.unique(
  1980             [
  1981                 np.amin(self._sampling_freqs),
  1982                 np.amax(self._sampling_freqs),
  1983                 np.mean(self._sampling_freqs),

---

tidy3d/components/mode_spec.py

Lines 197-205

  197     def _validate_freqs(cls, val):
  198         """Validate custom frequencies."""
  199         freqs_array = np.asarray(val)
  200         if freqs_array.size < 2:
! 201             raise ValidationError("Custom sampling requires at least 2 frequency points.")
  202         return val
  203 
  204     def sampling_points(self, freqs: FreqArray) -> FreqArray:
  205         """Return the custom frequency sampling points.

---

Lines 294-302

  294     def _validate_method_needs_points(cls, val, values):
  295         """Validate that the method has enough points."""
  296         sampling_spec = values.get("sampling_spec")
  297         if sampling_spec is None:
! 298             return val
  299 
  300         num_points = sampling_spec._num_points
  301         if val == "cubic" and num_points < 4:
  302             raise ValidationError(

---

Lines 681-689

  681         if track_freq is not None:
  682             return track_freq
  683         if sort_spec is not None:
  684             return sort_spec.track_freq
! 685         return None
  686 
  687     @pd.validator("interp_spec", always=True)
  688     @skip_if_fields_missing(["sort_spec", "track_freq"])
  689     def _interp_spec_needs_tracking(cls, val, values):

---

tidy3d/components/monitor.py

Lines 822-830

  822         """Size of monitor storage given the number of points after discretization."""
  823         amps_size = 3 * BYTES_COMPLEX * len(self._stored_freqs) * self.mode_spec.num_modes
  824         fields_size = 0
  825         if self.store_fields_direction is not None:
! 826             fields_size = (
  827                 6 * BYTES_COMPLEX * num_cells * len(self._stored_freqs) * self.mode_spec.num_modes
  828             )
  829             if self.mode_spec.precision == "double":
  830                 fields_size *= 2

---

[momchil-flex]

Nice! I didn't know this was coming!

I think we might want to enable it by default for ModeMonitor-s, or at least for all ports. But we could also sit on it for a bit, maybe let @lucas-flexcompute test this out as a default in photonforge first?

[caseyflex]

This is awesome @dbochkov-flexcompute !

[weiliangjin2021]

A very pleasant PR to review: both feature-rich, but not so many lines of code change!

As a future step, maybe it's useful to support adaptive sampling, as the mode profile probably changes more rapidly in high-frequency compared to low-frequency range. Or perhaps expose an API for custom sampling frequencies for now?

[momchil-flex]

Or perhaps expose an API for custom sampling frequencies for now?

Yes I was thinking about suggesting the same thing. Another application of custom sampling frequencies could be when you have a mode crossing, you could put more frequencies around it so that the tracking works.

[lucas-flexcompute]

Sorry if this is a dumb question, but I need to understand it first: is the ModeInterpSpec a superset of what ModeSource.num_freqs already does? How do they combine?

[dbochkov-flexcompute]

Sorry if this is a dumb question, but I need to understand it first: is the ModeInterpSpec a superset of what ModeSource.num_freqs already does? How do they combine?

There is definitely a big overlap between them, but ModeInterpSpec was not intended as a replacement for num_freqs. More like ModeInterpSpec is analogue of broadband source's num_freqs for monitors with more options. I think we could consider using ModeInterpSpec in sources as well and limit it there to only method="cheb" (and also maybe rename it to FreqInterpSpec as it would also apply to gaussian beam sources).

@momchil-flex @caseyflex this could actually make the decision whether to put ModeInterpSpec into ModeSpec or not easier. Because, if sources use ModeInterpSpec too, then no issue with them ignoring it in ModeSpec. However, for gaussian beams we still need to pass it separately. Btw, I don't think there is any benefit in using profile interpolation in incoming gaussian monitors, right?

[momchil-flex]

@momchil-flex @caseyflex this could actually make the decision whether to put ModeInterpSpec into ModeSpec or not easier. Because, if sources use ModeInterpSpec too, then no issue with them ignoring it in ModeSpec. However, for gaussian beams we still need to pass it separately.

So are you driving towards separating it also because of this argument? And I guess deprecating ModeSource.num_freqs for a FreqInterpSpec?

Btw, I don't think there is any benefit in using profile interpolation in incoming gaussian monitors, right?

I don't think so as the field computation is just analytical. Maybe a tiny benefit to interpolating.

[caseyflex]

@dbochkov-flexcompute I think it makes sense to incorporate it into ModeSpec and use it in mode sources as well

[dbochkov-flexcompute]

I think I will move ModeInterpSpec into ModeSpec, but not switch broadband sources to use it. This is because we have broadband gaussian beams and even plane wave, and it would also be weird to require method="cheb" only

[dbochkov-flexcompute]

Did a revision of this (turned out a bit more than expected):

• moved interp spec into mode spec
• added renormalization after interpolation
• slight reorg of basic dataset classes so that we can easily reorder and interpolate both main mode data and additional RF data attached to it (there is an option just to recalculate RF priorities though)
• swapped order of interpretors in data_raw so that group index and dispersion are calculated as previously

Still not sure whether we need to worry much about interpolating grid correction factor or better to recompute them.

Another question is whether we want to implement any sort of automatic fall back strategy "cheb" -> "cubic" -> "linear" -> "nearest" of we detect a mode discontinuity

[caseyflex]

I think it's better to recompute the grid correction factor, interpolating the complex exponential can lead to bad behavior (abs < 1)

I'm not sure about a fallback, but I am going to add method="poly" where the sampling points are closely concentrated around the central frequency (for more reliable tracking), and polynomial extrapolation is used

[momchil-flex]

One thing I remembered, you should add the ModeInterpSpec to the api docs. And while your at it could you also add the ModeSortSpec which I forgot in my PR?

[dbochkov-flexcompute]

I think it's better to recompute the grid correction factor, interpolating the complex exponential can lead to bad behavior (abs < 1)

Agree. I did have to add some new fields (grid_distances_primal and grid_distances_dual) into ModeSolverData though to enable that.

ModeSolverData.interpolated_data is working now too and I followed the following logic:

mode_solver_full_mnt = td.ModeSolverMonitor(..., reduce_data=False, name="full")
mode_solver_reduced_mnt = td.ModeSolverMonitor(..., reduce_data=True, name="reduced")
...
assert sim_data["full"] == sim_data["reduced"].interpolated_copy

[momchil-flex]

I think it's better to recompute the grid correction factor, interpolating the complex exponential can lead to bad behavior (abs < 1)

Agree. I did have to add some new fields (grid_distances_primal and grid_distances_dual) into ModeSolverData though to enable that.

Curious about this - the grid should be the same at all frequencies, no?

[dbochkov-flexcompute]

I think it's better to recompute the grid correction factor, interpolating the complex exponential can lead to bad behavior (abs < 1)

Agree. I did have to add some new fields (grid_distances_primal and grid_distances_dual) into ModeSolverData though to enable that.

Curious about this - the grid should be the same at all frequencies, no?

yeah, just after n_complex is interpolated onto new set of frequencies, we need to recompute grid_factor = |d2| / (|d1| + |d2|) * exp(1j * k * n * d1) + |d1| / (|d1| + |d2|) * exp(1j * k * n * d2) for each new frequency. Currently, grid factors are computed in ModeSolver where distances in normal direction d1 and d2 are known from the simulation. To allow this recalculation to happen in ModeSolverData.interp_in_freq() -> ModeSolverData we need to keep d1 and d2 in ModeSolverData as well

[caseyflex]

@dbochkov-flexcompute I left a few more comments, but I'd like to get this in soon, so I can clean up my EME PR. I guess the important thing is to get the API stabilized -- maybe we can do that and merge, and then I can change the internals while working on EME. To that end, I'd just consider moving reduce_data into ModeInterpSpec as the primary API change, then see if we can get this merged

[caseyflex]

This looks great to me, very flexible and extensible. thanks @dbochkov-flexcompute !

[dbochkov-flexcompute]

@greptile

[greptile-apps]

_{18 files reviewed, no comments}

_{Edit Code Review Agent Settings | Greptile}