-
Notifications
You must be signed in to change notification settings - Fork 66
fix(mode): ensure degenerate modes are orthogonal from mode solver #3057
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
base: develop
Are you sure you want to change the base?
fix(mode): ensure degenerate modes are orthogonal from mode solver #3057
Conversation
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
3 files reviewed, 2 comments
| scale_factor = np.sqrt(self_dot) | ||
| eigvecs[:, i] = vec / scale_factor | ||
| # 3. Calculate Inverse Square Root of Eigenvalues | ||
| lambda_inv_sqrt = np.diag(1.0 / np.sqrt(eigvals)) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
style: The eigenvalues from np.linalg.eig(S) could theoretically be complex (since S is a complex symmetric matrix, not Hermitian). Taking np.sqrt(eigvals) of complex eigenvalues will produce complex results, which is mathematically correct for this complex-symmetric decomposition. However, if any eigenvalue is zero or very close to zero, 1.0 / np.sqrt(eigvals) will cause a division by zero or produce very large values. Consider adding a check for near-zero eigenvalues similar to the eigenvector self-orthogonality check above.
Prompt To Fix With AI
This is a comment left during a code review.
Path: tidy3d/components/mode/solver.py
Line: 1257:1257
Comment:
**style:** The eigenvalues from `np.linalg.eig(S)` could theoretically be complex (since S is a complex symmetric matrix, not Hermitian). Taking `np.sqrt(eigvals)` of complex eigenvalues will produce complex results, which is mathematically correct for this complex-symmetric decomposition. However, if any eigenvalue is zero or very close to zero, `1.0 / np.sqrt(eigvals)` will cause a division by zero or produce very large values. Consider adding a check for near-zero eigenvalues similar to the eigenvector self-orthogonality check above.
How can I resolve this? If you propose a fix, please make it concise.| n = 1.5 | ||
| num_modes = 4 | ||
| mode_spec = td.ModeSpec(num_modes=num_modes) | ||
| medium = td.Medium(permittivity=n**2) # , conductivity=0.1) |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
style: The commented-out # , conductivity=0.1) should be removed to keep the code clean.
Context Used: Rule from dashboard - Remove commented-out or obsolete code; rely on version control for history. (source)
Prompt To Fix With AI
This is a comment left during a code review.
Path: tests/test_plugins/test_mode_solver.py
Line: 1517:1517
Comment:
**style:** The commented-out `# , conductivity=0.1)` should be removed to keep the code clean.
**Context Used:** Rule from `dashboard` - Remove commented-out or obsolete code; rely on version control for history. ([source](https://app.greptile.com/review/custom-context?memory=0ee31792-33ef-48e0-8328-48277c2ef438))
How can I resolve this? If you propose a fix, please make it concise.
caseyflex
left a comment
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Still testing with EME, sent you a slack message. Just a preliminary comment, will give more feedback tomorrow. This is great though, thanks @dmarek-flex !
| ) | ||
| assert len(degen_sets) == 1 | ||
|
|
||
| S = mode_data.outer_dot(mode_data, conjugate=True).isel(f=0).values |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
shouldn't this be conjugate=False?
Remaining considerations:
Greptile Overview
Greptile Summary
This PR fixes a mathematical issue with degenerate modes in the mode solver. When eigenmodes have identical or nearly identical eigenvalues (degenerate modes), the eigensolver doesn't guarantee orthogonality between them. This PR ensures degenerate modes respect the bi-orthogonality condition by:
TOL_DEGENERATEconstant for identifying degenerate modes_identify_degenerate_modes()to group modes with similar effective indices_dot()and_outer_dot()for computing bi-orthogonal overlap integrals_normalize_modes()to normalize modes to unit self-overlap_make_orthogonal_basis_for_degenerate_modes()using a complex-symmetric eigendecomposition approachThe implementation correctly handles the complex-symmetric (not Hermitian) nature of the bi-orthogonal inner product, using unconjugated dot products for eigenvector normalization.
Confidence Score: 4/5
tidy3d/components/mode/solver.py- the_make_orthogonal_basis_for_degenerate_modesmethod could benefit from a check for near-zero eigenvalues.Important Files Changed
File Analysis
test_degenerate_mode_processingthat verifies degenerate modes are bi-orthogonal by checking off-diagonal overlap values are below threshold.TOL_DEGENERATEconstant and implemented degenerate mode orthogonalization via_identify_degenerate_modes,_dot,_outer_dot,_normalize_modes, and_make_orthogonal_basis_for_degenerate_modesmethods. The implementation handles complex symmetric matrices correctly.Sequence Diagram
sequenceDiagram participant User participant ModeSolver participant EigSolver participant compute_modes User->>ModeSolver: data_raw ModeSolver->>compute_modes: compute_modes(...) compute_modes->>EigSolver: solver_em() EigSolver-->>compute_modes: E, H, neff, keff Note over compute_modes: Transform back to original axes compute_modes->>EigSolver: _normalize_modes(E, H, dl_f, dl_b, all_modes) Note over EigSolver: Normalize each mode to unit self-overlap EigSolver-->>compute_modes: E, H (normalized) compute_modes->>EigSolver: _identify_degenerate_modes(n_complex, TOL_DEGENERATE) Note over EigSolver: Group modes with similar n_complex EigSolver-->>compute_modes: degenerate_sets: [{0,1}, ...] compute_modes->>EigSolver: _make_orthogonal_basis_for_degenerate_modes(...) loop For each degenerate group EigSolver->>EigSolver: _outer_dot() → compute overlap matrix S EigSolver->>EigSolver: np.linalg.eig(S) → diagonalize EigSolver->>EigSolver: Renormalize eigenvectors (v.T @ v = 1) EigSolver->>EigSolver: W = V @ Λ^(-1/2) @ V.T EigSolver->>EigSolver: Transform E, H by W end EigSolver-->>compute_modes: E, H (orthogonalized) compute_modes-->>ModeSolver: fields, n_complex, eps_spec ModeSolver-->>User: ModeSolverData