add a TerminalSpec for RF specific mode information #2444
Conversation
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
|
Still very rough, and I am not particularly happy about needing to define all of these different path spec types. But right now it is the best way I found to get around the circular dependency issue I am having. Basically, I cannot import path integral code into @weiliangjin2021 @dbochkov-flexcompute Any ideas to solve this? Otherwise you can retrieve impedance info by doing something as simple as:
|
|
Nah, having path specs seems to be the best. I'll just add some factory class to convert path specs to the correct path integral version. |
dmarek-flex
force-pushed
the
dmarek/add_terminal_mode_monitor
branch
2 times, most recently
from
44634db to
e12f1a7
Compare
|
Ready for review finally! The main feature this PR provided is automatically computing impedance and adding it to the
to some object that calculates and returns mode data. The simulation will not be affected if used in a Now, to get around all of the circular dependence stuff, I had to split the path integral classes in the microwave plugin into the specification part, and the part that actually computes the integral. At some later point we will most likely move the everything out of the microwave plugin, since these path integrals are needed for many things in RF. If anyone has suggestions on how to better get around these issues, please let me know! The main issue was the need to import the monitor data module in the path integral modules. Finally, the auto generated paths are bounding boxes around each conductor, which will work for 99% of transmission lines, but there might be improvements needed in the future to auto generate the line string type path integrals for handling more complex structures. |
|
@yuanshen-flexcompute It would be nice to test out this automatic version of impedance calculation on some of your test cases to make sure I did not miss anything. |
|
Code Coverage SummaryDiff against developResults for commit: c5f03f3 Minimum allowed coverage is ♻️ This comment has been updated with latest results |
|
Changed Files CoverageResults for commit: c5f03f3 Minimum allowed coverage is ♻️ This comment has been updated with latest results |
dmarek-flex
force-pushed
the
dmarek/add_terminal_mode_monitor
branch
from
3781b06 to
c5f03f3
Compare
Comment on lines
2016
to
2020
| from tidy3d.components.microwave.path_integral_factory import ( | ||
| make_current_integral, | ||
| make_voltage_integral, | ||
| ) | ||
| from tidy3d.plugins.microwave.impedance_calculator import ImpedanceCalculator |
There was a problem hiding this comment.
any reason not to import these at the top level?
There was a problem hiding this comment.
circular dependency 😢 Many of those popped up in this PR
There was a problem hiding this comment.
Yes I think this is part of why we need to do the RF refactor soon on the lines of what we chatted back on Friday too.
dmarek-flex
force-pushed
the
dmarek/add_terminal_mode_monitor
branch
3 times, most recently
from
3f3195e to
39318d6
Compare
Diff CoverageDiff: origin/develop...HEAD, staged and unstaged changes
Summary
tidy3d/components/data/monitor_data.pyLines 2040-2049 2040 info["wg TE fraction"] = self.pol_fraction_waveguide["te"]
2041 info["wg TM fraction"] = self.pol_fraction_waveguide["tm"]
2042
2043 if self.Z0 is not None:
! 2044 info["Re(Z0)"] = self.Z0.real
! 2045 info["Im(Z0)"] = self.Z0.imag
2046
2047 return xr.Dataset(data_vars=info)
2048
2049 def to_dataframe(self) -> DataFrame:tidy3d/components/geometry/utils.pyLines 71-90 71 Flat list of non-empty geometries matching the specified types.
72 """
73 # Handle single Shapely object by wrapping it in a list
74 if isinstance(geoms, Shapely):
! 75 geoms = [geoms]
76
77 flat = []
78 for geom in geoms:
79 if geom.is_empty:
! 80 continue
81 if isinstance(geom, keep_types):
82 flat.append(geom)
! 83 elif isinstance(geom, (MultiPolygon, MultiLineString, MultiPoint, GeometryCollection)):
! 84 flat.extend(flatten_shapely_geometries(geom.geoms, keep_types))
! 85 elif isinstance(geom, BaseGeometry) and hasattr(geom, "geoms"):
! 86 flat.extend(flatten_shapely_geometries(geom.geoms, keep_types))
87 return flat
88
89
90 def merging_geometries_on_plane(Lines 509-517 509 else: # SnapType.Contract
510 min_upper_bound_idx += snap_margin
511 max_upper_bound_idx -= snap_margin
512 if max_upper_bound_idx < min_upper_bound_idx:
! 513 raise SetupError("The supplied 'snap_buffer' is too large for this contraction.")
514 min_snap = get_upper_bound(interval_min, coords, min_upper_bound_idx, rel_tol=rtol)
515 max_snap = get_lower_bound(interval_max, coords, max_upper_bound_idx, rel_tol=rtol)
516 return (min_snap, max_snap)tidy3d/components/microwave/path_integral_factory.pyLines 117-130 117 raise SetupError(
118 f"Failed to auto-generate path specification for impedance calculation in monitor '{monitor.name}'."
119 ) from e
120
! 121 try:
! 122 v_integral = make_voltage_integral(v_spec)
! 123 i_integral = make_current_integral(i_spec)
! 124 except Exception as e:
! 125 raise SetupError(
126 f"Failed to construct path integrals from the terminal specification for monitor '{monitor.name}'. "
127 "Please create a github issue so that the problem can be investigated."
128 ) from e
! 129 return (v_integral, i_integral)tidy3d/components/microwave/path_spec.pyLines 103-111 103 """Axis for performing integration."""
104 for index, value in enumerate(self.size):
105 if value != 0:
106 return index
! 107 raise Tidy3dError("Failed to identify axis.")
108
109 def _vertices_2D(self, axis: Axis) -> tuple[Coordinate2D, Coordinate2D]:
110 """Returns the two vertices of this path in the plane defined by ``axis``."""
111 min = self.bounds[0]Lines 161-183 161 -------
162 VoltageIntegralAxisAlignedSpec
163 The created path integral for computing voltage between the two terminals.
164 """
! 165 axis_positions = Geometry.parse_two_xyz_kwargs(x=x, y=y, z=z)
166 # Calculate center and size of the future box
! 167 midpoint = (plus_terminal + minus_terminal) / 2
! 168 length = np.abs(plus_terminal - minus_terminal)
! 169 center = [midpoint, midpoint, midpoint]
! 170 size = [length, length, length]
! 171 for axis, position in axis_positions:
! 172 size[axis] = 0
! 173 center[axis] = position
174
! 175 direction = "+"
! 176 if plus_terminal < minus_terminal:
! 177 direction = "-"
178
! 179 return VoltageIntegralAxisAlignedSpec(
180 center=center,
181 size=size,
182 extrapolate_to_endpoints=extrapolate_to_endpoints,
183 snap_path_to_grid=snap_path_to_grid,Lines 265-273 265 """Axis normal to loop"""
266 for index, value in enumerate(self.size):
267 if value == 0:
268 return index
! 269 raise Tidy3dError("Failed to identify axis.")
270
271 def _to_path_integral_specs(
272 self, h_horizontal=None, h_vertical=None
273 ) -> tuple[AxisAlignedPathIntegralSpec, ...]:Lines 464-478 464
465 @staticmethod
466 def _compute_dl_component(coord_array: xr.DataArray, closed_contour=False) -> np.array:
467 """Computes the differential length element along the integration path."""
! 468 dl = np.gradient(coord_array)
! 469 if closed_contour:
470 # If the contour is closed, we can use central difference on the starting/end point
471 # which will be more accurate than the default forward/backward choice in np.gradient
! 472 grad_end = np.gradient([coord_array[-2], coord_array[0], coord_array[1]])
! 473 dl[0] = dl[-1] = grad_end[1]
! 474 return dl
475
476 @classmethod
477 def from_circular_path(
478 cls, center: Coordinate, radius: float, num_points: int, normal_axis: Axis, clockwise: boolLines 499-530 499 :class:`.CustomPathIntegral2DSpec`
500 A path integral defined on a circular path.
501 """
502
! 503 def generate_circle_coordinates(radius: float, num_points: int, clockwise: bool):
504 """Helper for generating x,y vertices around a circle in the local coordinate frame."""
! 505 sign = 1.0
! 506 if clockwise:
! 507 sign = -1.0
! 508 angles = np.linspace(0, sign * 2 * np.pi, num_points, endpoint=True)
! 509 xt = radius * np.cos(angles)
! 510 yt = radius * np.sin(angles)
! 511 return (xt, yt)
512
513 # Get transverse axes
! 514 normal_center, trans_center = Geometry.pop_axis(center, normal_axis)
515
516 # These x,y coordinates in the local coordinate frame
! 517 if normal_axis == 1:
518 # Handle special case when y is the axis that is popped
! 519 clockwise = not clockwise
! 520 xt, yt = generate_circle_coordinates(radius, num_points, clockwise)
! 521 xt += trans_center[0]
! 522 yt += trans_center[1]
! 523 circle_vertices = np.column_stack((xt, yt))
524 # Close the contour exactly
! 525 circle_vertices[-1, :] = circle_vertices[0, :]
! 526 return cls(axis=normal_axis, position=normal_center, vertices=circle_vertices)
527
528 @cached_property
529 def is_closed_contour(self) -> bool:
530 """Returns ``true`` when the first vertex equals the last vertex."""Lines 562-578 562
563 @cached_property
564 def sign(self) -> Direction:
565 """Uses the ordering of the vertices to determine the direction of the current flow."""
! 566 linestr = shapely.LineString(coordinates=self.vertices)
! 567 is_ccw = shapely.is_ccw(linestr)
568 # Invert statement when the vertices are given as (x, z)
! 569 if self.axis == 1:
! 570 is_ccw = not is_ccw
! 571 if is_ccw:
! 572 return "+"
573 else:
! 574 return "-"
575
576
577 class CustomVoltageIntegral2DSpec(CustomPathIntegral2DSpec):
578 """Class for specfying the computation of voltage between two points defined by a custom path.Lines 697-705 697 linestr = shapely.LineString(coordinates=self.vertices)
698 is_ccw = shapely.is_ccw(linestr)
699 # Invert statement when the vertices are given as (x, z)
700 if self.axis == 1:
! 701 is_ccw = not is_ccw
702 if is_ccw:
703 return "+"
704 else:
705 return "-"tidy3d/components/mode/mode_solver.pyLines 484-492 484
485 mode_solver_data = self._filter_components(mode_solver_data)
486 # Calculate and add the characteristic impedance
487 if self.mode_spec.terminal_spec is not None:
! 488 mode_solver_data = self._add_characteristic_impedance(mode_solver_data)
489 return mode_solver_data
490
491 @cached_property
492 def bend_axis_3d(self) -> Axis:Lines 1296-1312 1296 def _add_characteristic_impedance(self, mode_solver_data: ModeSolverData):
1297 """Calculate and add characteristic impedance to ``mode_solver_data`` with the path specifications.
1298 If they were not supplied by the user, then create a specification automatically.
1299 """
! 1300 v_integral, i_integral = make_path_integrals(
1301 self.mode_spec.terminal_spec, self.to_monitor(name=MODE_MONITOR_NAME), self.simulation
1302 )
! 1303 impedance_calc = ImpedanceCalculator(
1304 voltage_integral=v_integral, current_integral=i_integral
1305 )
1306 # Need to operate on the full symmetry expanded fields
! 1307 Z0 = impedance_calc.compute_impedance(mode_solver_data.symmetry_expanded_copy)
! 1308 return mode_solver_data.updated_copy(Z0=Z0)
1309
1310 @cached_property
1311 def data(self) -> ModeSolverData:
1312 """:class:`.ModeSolverData` containing the field and effective index data.tidy3d/plugins/microwave/custom_path_integrals.pyLines 229-279 229 """Current integral comprising one or more disjoint paths"""
230
231 def compute_current(self, em_field: MonitorDataTypes) -> IntegralResultTypes:
232 """Compute current flowing in loop defined by the outer edge of a rectangle."""
! 233 if isinstance(em_field, FieldTimeData) and self.sum_spec == "split":
! 234 raise DataError(
235 "Only frequency domain field data is supported when using the 'split' sum_spec. "
236 "Either switch the sum_spec to 'sum' or supply frequency domain data."
237 )
238
! 239 from tidy3d.components.microwave.path_integral_factory import make_current_integral
240
! 241 current_integrals = [make_current_integral(path_spec) for path_spec in self.path_specs]
242
243 # Initialize arrays with first current term
! 244 first_term = current_integrals[0].compute_current(em_field)
! 245 current_in_phase = xr.zeros_like(first_term)
! 246 current_out_phase = xr.zeros_like(first_term)
247
248 # Get reference phase from first non-zero current
! 249 phase_reference = None
! 250 for path in current_integrals:
! 251 term = path.compute_current(em_field)
! 252 if np.any(abs(term) > 0):
! 253 phase_reference = np.angle(term)
! 254 break
! 255 if phase_reference is None:
! 256 raise DataError("Cannot complete calculation of current. No non-zero current found.")
257
258 # Accumulate currents based on phase comparison
! 259 for path in current_integrals:
! 260 term = path.compute_current(em_field)
! 261 if np.all(abs(term) == 0):
! 262 continue
263
264 # Compare phase to reference
! 265 phase_diff = np.angle(term) - phase_reference
266 # Wrap phase difference to [-pi, pi]
! 267 phase_diff = np.mod(phase_diff + np.pi, 2 * np.pi) - np.pi
268
269 # Check phase consistency across frequencies
! 270 freq_axis = term.get_axis_num("f")
! 271 all_in_phase = np.all(abs(phase_diff) <= np.pi / 2, axis=freq_axis)
! 272 all_out_of_phase = np.all(abs(phase_diff) > np.pi / 2, axis=freq_axis)
! 273 consistent_phase = np.logical_or(all_in_phase, all_out_of_phase)
! 274 if not np.all(consistent_phase) and self.sum_spec == "split":
! 275 log.warning(
276 "Phase alignment of current is not consistent across frequencies. "
277 "The current calculation may be inaccurate."
278 )Lines 277-305 277 "The current calculation may be inaccurate."
278 )
279
280 # Add to in-phase or out-of-phase current
! 281 is_in_phase = abs(phase_diff) <= np.pi / 2
! 282 current_in_phase += xr.where(is_in_phase, term, 0)
! 283 current_out_phase += xr.where(~is_in_phase, term, 0)
284
! 285 if self.sum_spec == "sum":
! 286 return CurrentIntegralAxisAligned._set_data_array_attributes(
287 current_in_phase + current_out_phase
288 )
289
290 # For split mode, return the larger magnitude current
! 291 freq_axis = current_in_phase.get_axis_num("f")
! 292 in_all_larger = np.all(abs(current_in_phase) >= abs(current_out_phase), axis=freq_axis)
! 293 in_all_smaller = np.all(abs(current_in_phase) < abs(current_out_phase), axis=freq_axis)
! 294 consistent_max_current = np.logical_or(in_all_larger, in_all_smaller)
! 295 if not np.all(consistent_max_current) and self.sum_spec == "split":
! 296 log.warning(
297 "There is not a consistently larger current across frequencies between the in "
298 "and out of phase components. The current calculation may be inaccurate."
299 )
300
! 301 current = xr.where(
302 abs(current_in_phase) >= abs(current_out_phase), current_in_phase, current_out_phase
303 )
! 304 return CurrentIntegralAxisAligned._set_data_array_attributes(current) |
adding support for colocated fields fixing bug for more complicated shapes adding tests, fixing symmetric conditions, fixing for 2D structures fix simulation validator for terminal spec refactor by splitting the path integral specification away from the integral computation, now the path specification handles pretty much everything, except for the final integral computation and results preparation rename auto patch spec to generator type name, fixed tests add test for sim validation and improved doc strings fix regression fix python 3.9 tests rebase on ruff changes reorg of ModeData, the impedance is now calculated by the ModeSolver class avoiding the need for a duplicate terminal spec
dmarek-flex
force-pushed
the
dmarek/add_terminal_mode_monitor
branch
from
39318d6 to
469d8ae
Compare
|
Planning on changing the name to |
That does make more sense. Are there any other impedance specs that might be introduced in future? Just wondering if it would also make sense to expand it to |
Yea I think |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Add this suggestion to a batch that can be applied as a single commit.
This suggestion is invalid because no changes were made to the code.
Suggestions cannot be applied while the pull request is closed.
Suggestions cannot be applied while viewing a subset of changes.
Only one suggestion per line can be applied in a batch.
Add this suggestion to a batch that can be applied as a single commit.
Applying suggestions on deleted lines is not supported.
You must change the existing code in this line in order to create a valid suggestion.
Outdated suggestions cannot be applied.
This suggestion has been applied or marked resolved.
Suggestions cannot be applied from pending reviews.
Suggestions cannot be applied on multi-line comments.
Suggestions cannot be applied while the pull request is queued to merge.
Suggestion cannot be applied right now. Please check back later.
Calculating impedance automatically:
CompositeCurrentIntegral) made up of multiple current integrals. The total current flowing in a transmission line is net 0 because of current flowing in and out . To isolate these two terms, we split the currents based on their phase (using an initial reference phase). Then we choose the maximum of the two at the end. We choose the max, because it is possible that some current is flowing back in the PEC boundary or ground plane, which are not totally enclosed by the modal plane and their contribution will be missing.Need at least one small backend change
Todo on: