Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 1 addition & 1 deletion gplugins/tidy3d/materials.py
Original file line number Diff line number Diff line change
Expand Up @@ -9,7 +9,7 @@

material_name_to_tidy3d = {
"si": td.material_library["cSi"]["Li1993_293K"],
"sio2": td.material_library["SiO2"]["Horiba"],
"sio2": td.material_library["SiO2"]["Palik_Lossless"],
"sin": td.material_library["Si3N4"]["Luke2015PMLStable"],
}

Expand Down
48 changes: 42 additions & 6 deletions gplugins/tidy3d/modes.py
Original file line number Diff line number Diff line change
Expand Up @@ -92,6 +92,8 @@ class Waveguide(BaseModel, extra="forbid", arbitrary_types_allowed=True):
surface_k: absorption coefficient added to the core material
index on the top-surface layer.
bend_radius: radius to simulate circular bend.
target_neff: target effective index for the mode solver. Defaults
to the real part of the core refractive index if not specified.
num_modes: number of modes to compute.
group_index_step: if set to `True`, indicates that the group
index must also be calculated. If set to a positive float
Expand Down Expand Up @@ -147,6 +149,7 @@ class Waveguide(BaseModel, extra="forbid", arbitrary_types_allowed=True):
surface_thickness: float = 0.0
surface_k: float = 0.0
bend_radius: float | None = None
target_neff: float | None = None
num_modes: int = 2
group_index_step: bool | float = False
precision: Precision = "double"
Expand Down Expand Up @@ -178,6 +181,10 @@ def filepath(self) -> pathlib.Path | None:
h = hashlib.md5(named_args_string.encode()).hexdigest()[:16]
return cache_path / f"{self.__class__.__name__}_{h}.npz"

def _resolve_target_neff(self, n_core: complex) -> float:
"""Return target_neff if set, otherwise fall back to n_core.real."""
return self.target_neff if self.target_neff is not None else n_core.real

@property
def waveguide(self):
"""Tidy3D waveguide used by this instance."""
Expand Down Expand Up @@ -222,9 +229,11 @@ def waveguide(self):
else None
)

target_neff = self._resolve_target_neff(n_core)

mode_spec = td.ModeSpec(
num_modes=self.num_modes,
target_neff=n_core.real,
target_neff=target_neff,
bend_radius=self.bend_radius,
bend_axis=1,
num_pml=(12, 12) if self.bend_radius else (0, 0),
Expand Down Expand Up @@ -502,6 +511,8 @@ class WaveguideCoupler(Waveguide):
index must also be calculated. If set to a positive float
it defines the fractional frequency step used for the
numerical differentiation of the effective index.
target_neff: target effective index for the mode solver. Defaults
to the real part of the core refractive index if not specified.
precision: computation precision.
grid_resolution: wavelength resolution of the computation grid.
max_grid_scaling: grid scaling factor in cladding regions.
Expand Down Expand Up @@ -554,9 +565,11 @@ def waveguide(self):
else None
)

target_neff = self._resolve_target_neff(n_core)

mode_spec = td.ModeSpec(
num_modes=self.num_modes,
target_neff=n_core.real,
target_neff=target_neff,
bend_radius=self.bend_radius,
bend_axis=1,
num_pml=(12, 12) if self.bend_radius else (0, 0),
Expand Down Expand Up @@ -863,7 +876,10 @@ def sweep_mode_area(waveguide: Waveguide, **sweep_kwargs) -> np.ndarray:


def sweep_bend_mismatch(
waveguide: Waveguide, bend_radii: tuple[float, ...]
waveguide: Waveguide,
bend_radii: tuple[float, ...],
track_modes: bool = False,
modes_to_track: Sequence[int] = (0,),
) -> np.ndarray:
"""Overlap integral squared for the bend mode mismatch loss.

Expand All @@ -873,7 +889,22 @@ def sweep_bend_mismatch(
Args:
waveguide: base waveguide geometry.
bend_radii: radii values to sweep.
track_modes: if True, for each radius select the bend mode with
the best overlap for each tracked straight mode.
modes_to_track: straight mode indices to track. Required when
track_modes is True.
"""
if track_modes:
if len(modes_to_track) == 0:
raise ValueError("modes_to_track must be provided when track_modes is True")
if waveguide.num_modes < max(modes_to_track) + 1:
raise ValueError(
f"num_modes ({waveguide.num_modes}) must be >= "
f"{max(modes_to_track) + 1} to track modes {modes_to_track}"
)
if waveguide.num_modes < 2:
raise ValueError("Track modes requires num_modes >= 2")
Comment on lines +897 to +906
Copy link
Copy Markdown
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

suggestion (bug_risk): Revisit squaring strategy to avoid double-squaring and align with the documented quantity.

In the track_modes branch you compute best = [np.max(np.abs(overlap[:, m]) ** 2) ...] and later return np.abs(results) ** 2, so the output is effectively |overlap|**4. This matches the non-tracking branch, but the docstring says it returns the “overlap integral squared”. Please decide whether the intended quantity is |overlap|**2 or |overlap|**4 and update both branches consistently. If only the squared magnitude is intended, remove the inner ** 2 (and possibly the final ** 2) to avoid redundant work and mismatch with the documentation.

Suggested implementation:

        # Take the maximum absolute overlap for each tracked straight mode.
        # The squaring to obtain |overlap|**2 is applied once at the end of the function.
        best = [np.max(np.abs(overlap[:, m]), axis=1) for m in modes_to_track]
        # Take the maximum absolute overlap across all modes; squaring is done once at the end.
        results = np.max(np.abs(overlap), axis=1)
    # Return the overlap integral squared, i.e. |overlap|**2, in line with the docstring.
    return np.abs(results) ** 2

These edits assume:

  1. The tracked-modes branch contains a line like:
    best = [np.max(np.abs(overlap[:, m]) ** 2, axis=1) for m in modes_to_track]
  2. The non-tracking branch contains a line like:
    results = np.max(np.abs(overlap) ** 2, axis=1)
  3. The function currently ends with return np.abs(results) ** 2.

If the actual code uses slightly different variable names or shapes (e.g. axis=0 vs axis=1, or builds results from best before returning), you should:

  • Remove any inner ** 2 that is applied directly to np.abs(overlap...) in both the tracking and non-tracking branches.
  • Keep a single squaring step right before returning, so the final returned quantity is |overlap|**2, matching the docstring “Overlap integral squared”.
  • Ensure intermediate comments (if any) do not claim |overlap|**4.

Comment on lines +905 to +906
Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

medium

The restriction waveguide.num_modes < 2 for mode tracking seems unnecessary. Even if a waveguide has only one mode, a user might want to use the track_modes logic for consistency in their scripts, and the underlying overlap calculation overlap[:, m] works correctly for a single mode (where m=0).


kwargs = dict(waveguide)
kwargs.pop("bend_radius")
straight = Waveguide(**kwargs)
Expand All @@ -882,9 +913,14 @@ def sweep_bend_mismatch(
for radius in tqdm(bend_radii):
bend = Waveguide(bend_radius=radius, **kwargs)
overlap = bend.overlap(straight)
results.append(
np.diagonal(overlap) ** 2 if straight.num_modes > 1 else overlap**2
)

if track_modes:
best = [np.max(np.abs(overlap[:, m]) ** 2) for m in modes_to_track]
Copy link
Copy Markdown

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

high

There is a potential logic issue with double squaring when track_modes is enabled.

Line 918 calculates the power overlap as np.abs(overlap) ** 2. Then, line 925 squares the entire results list again: return np.abs(results) ** 2. This results in returning the overlap to the fourth power ($|\eta|^4$).

While the docstring mentions the loss is squared for a round trip, the naming 'Overlap integral squared' usually implies $|\eta|^2$. If $|\eta|^4$ is indeed the intended return value for this function (to represent the round-trip mismatch), then the implementation is consistent with the else branch. However, if the intention was to return the power overlap $|\eta|^2$, then line 918 should not be squared, or line 925 should be adjusted.

results.append(best)
else:
results.append(
np.diagonal(overlap) ** 2 if straight.num_modes > 1 else overlap**2
)

return np.abs(results) ** 2

Expand Down
2 changes: 1 addition & 1 deletion gplugins/tidy3d/tests/test_modes_coupler.py
Original file line number Diff line number Diff line change
Expand Up @@ -18,7 +18,7 @@ def test_neff() -> None:
clad_material="sio2",
)
n_eff = wg.n_eff[0].real
n_eff_ref = 2.5743837634515767
n_eff_ref = 2.6076116125083564
assert np.isclose(n_eff, n_eff_ref, rtol=0.01), n_eff


Expand Down
Loading