diff --git a/docs/api/plugins/index.rst b/docs/api/plugins/index.rst index e25a596631..a70209857d 100644 --- a/docs/api/plugins/index.rst +++ b/docs/api/plugins/index.rst @@ -12,6 +12,7 @@ Plugins adjoint waveguide design + microwave .. include:: /api/plugins/mode_solver.rst @@ -22,3 +23,4 @@ Plugins .. include:: /api/plugins/adjoint.rst .. include:: /api/plugins/waveguide.rst .. include:: /api/plugins/design.rst +.. include:: /api/plugins/microwave.rst \ No newline at end of file diff --git a/docs/api/plugins/microwave.rst b/docs/api/plugins/microwave.rst new file mode 100644 index 0000000000..dc963f691f --- /dev/null +++ b/docs/api/plugins/microwave.rst @@ -0,0 +1,13 @@ +.. currentmodule:: tidy3d + +Microwave +---------------------------- + +.. autosummary:: + :toctree: ../_autosummary/ + :template: module.rst + + tidy3d.plugins.microwave.AxisAlignedPathIntegral + tidy3d.plugins.microwave.VoltageIntegralAxisAligned + tidy3d.plugins.microwave.CurrentIntegralAxisAligned + tidy3d.plugins.microwave.ImpedanceCalculator \ No newline at end of file diff --git a/tests/test_plugins/test_microwave.py b/tests/test_plugins/test_microwave.py index 9ddf8b6174..bad8ff8941 100644 --- a/tests/test_plugins/test_microwave.py +++ b/tests/test_plugins/test_microwave.py @@ -4,7 +4,11 @@ import tidy3d as td from tidy3d import FieldData from tidy3d.constants import ETA_0 -from tidy3d.plugins.microwave import VoltageIntegralAA, CurrentIntegralAA, ImpedanceCalculator +from tidy3d.plugins.microwave import ( + VoltageIntegralAxisAligned, + CurrentIntegralAxisAligned, + ImpedanceCalculator, +) import pydantic.v1 as pydantic from tidy3d.exceptions import DataError from ..utils import run_emulated @@ -17,7 +21,7 @@ FSTOP = 1.5e9 F0 = (FSTART + FSTOP) / 2 FWIDTH = FSTOP - FSTART -FS = np.linspace(FSTART, FSTOP, 5) +FS = np.linspace(FSTART, FSTOP, 3) FIELD_MONITOR = td.FieldMonitor( size=MON_SIZE, fields=FIELDS, name="strip_field", freqs=FS, colocate=False ) @@ -33,8 +37,13 @@ td.FieldMonitor( center=(0, 0, 0), size=(1, 1, 1), freqs=FS, fields=["Ex", "Hx"], name="ExHx" ), - td.ModeMonitor( - center=(0, 0, 0), size=(1, 1, 0), freqs=FS, mode_spec=td.ModeSpec(), name="mode" + td.FieldTimeMonitor(center=(0, 0, 0), size=(1, 1, 0), colocate=False, name="field_time"), + td.ModeSolverMonitor( + center=(0, 0, 0), + size=(1, 1, 0), + freqs=FS, + mode_spec=td.ModeSpec(num_modes=2), + name="mode", ), ], sources=[ @@ -44,9 +53,11 @@ source_time=td.GaussianPulse(freq0=F0, fwidth=FWIDTH), ) ], - run_time=2e-12, + run_time=5e-16, ) +SIM_Z_DATA = run_emulated(SIM_Z) + """ Generate the data arrays for testing path integral computations """ @@ -105,13 +116,13 @@ def test_voltage_integral_axes(axis): size = [0, 0, 0] size[axis] = length center = [0, 0, 0] - voltage_integral = VoltageIntegralAA( + voltage_integral = VoltageIntegralAxisAligned( center=center, size=size, + sign="+", ) - sim = SIM_Z - sim_data = run_emulated(sim) - _ = voltage_integral.compute_voltage(sim_data["field"].field_components) + + _ = voltage_integral.compute_voltage(SIM_Z_DATA["field"]) @pytest.mark.parametrize("axis", [0, 1, 2]) @@ -120,13 +131,12 @@ def test_current_integral_axes(axis): size = [length, length, length] size[axis] = 0.0 center = [0, 0, 0] - current_integral = CurrentIntegralAA( + current_integral = CurrentIntegralAxisAligned( center=center, size=size, + sign="+", ) - sim = SIM_Z - sim_data = run_emulated(sim) - _ = current_integral.compute_current(sim_data["field"].field_components) + _ = current_integral.compute_current(SIM_Z_DATA["field"]) def test_voltage_integral_toggles(): @@ -134,16 +144,14 @@ def test_voltage_integral_toggles(): size = [0, 0, 0] size[0] = length center = [0, 0, 0] - voltage_integral = VoltageIntegralAA( + voltage_integral = VoltageIntegralAxisAligned( center=center, size=size, extrapolate_to_endpoints=True, snap_path_to_grid=True, sign="-", ) - sim = SIM_Z - sim_data = run_emulated(sim) - _ = voltage_integral.compute_voltage(sim_data["field"].field_components) + _ = voltage_integral.compute_voltage(SIM_Z_DATA["field"]) def test_current_integral_toggles(): @@ -151,16 +159,14 @@ def test_current_integral_toggles(): size = [length, length, length] size[0] = 0.0 center = [0, 0, 0] - current_integral = CurrentIntegralAA( + current_integral = CurrentIntegralAxisAligned( center=center, size=size, extrapolate_to_endpoints=True, snap_contour_to_grid=True, sign="-", ) - sim = SIM_Z - sim_data = run_emulated(sim) - _ = current_integral.compute_current(sim_data["field"].field_components) + _ = current_integral.compute_current(SIM_Z_DATA["field"]) def test_voltage_missing_fields(): @@ -168,14 +174,14 @@ def test_voltage_missing_fields(): size = [0, 0, 0] size[1] = length center = [0, 0, 0] - voltage_integral = VoltageIntegralAA( + voltage_integral = VoltageIntegralAxisAligned( center=center, size=size, + sign="+", ) - sim = SIM_Z - sim_data = run_emulated(sim) + with pytest.raises(DataError): - _ = voltage_integral.compute_voltage(sim_data["ExHx"].field_components) + _ = voltage_integral.compute_voltage(SIM_Z_DATA["ExHx"]) def test_current_missing_fields(): @@ -183,14 +189,42 @@ def test_current_missing_fields(): size = [length, length, length] size[0] = 0.0 center = [0, 0, 0] - current_integral = CurrentIntegralAA( + current_integral = CurrentIntegralAxisAligned( center=center, size=size, + sign="+", ) - sim = SIM_Z - sim_data = run_emulated(sim) + with pytest.raises(DataError): - _ = current_integral.compute_current(sim_data["ExHx"].field_components) + _ = current_integral.compute_current(SIM_Z_DATA["ExHx"]) + + +def test_time_monitor_voltage_integral(): + length = 0.5 + size = [0, 0, 0] + size[1] = length + center = [0, 0, 0] + voltage_integral = VoltageIntegralAxisAligned( + center=center, + size=size, + sign="+", + ) + + voltage_integral.compute_voltage(SIM_Z_DATA["field_time"]) + + +def test_mode_solver_monitor_voltage_integral(): + length = 0.5 + size = [0, 0, 0] + size[1] = length + center = [0, 0, 0] + voltage_integral = VoltageIntegralAxisAligned( + center=center, + size=size, + sign="+", + ) + + voltage_integral.compute_voltage(SIM_Z_DATA["mode"]) def test_tiny_voltage_path(): @@ -198,10 +232,11 @@ def test_tiny_voltage_path(): size = [0, 0, 0] size[1] = length center = [0, 0, 0] - voltage_integral = VoltageIntegralAA(center=center, size=size, extrapolate_to_endpoints=True) - sim = SIM_Z - sim_data = run_emulated(sim) - _ = voltage_integral.compute_voltage(sim_data["field"].field_components) + voltage_integral = VoltageIntegralAxisAligned( + center=center, size=size, sign="+", extrapolate_to_endpoints=True + ) + + _ = voltage_integral.compute_voltage(SIM_Z_DATA["field"]) def test_impedance_calculator(): @@ -209,16 +244,42 @@ def test_impedance_calculator(): _ = ImpedanceCalculator(voltage_integral=None, current_integral=None) +def test_impedance_calculator_on_time_data(): + # Setup path integrals + length = 0.5 + size = [0, length, 0] + size[1] = length + center = [0, 0, 0] + voltage_integral = VoltageIntegralAxisAligned( + center=center, size=size, sign="+", extrapolate_to_endpoints=True + ) + + size = [length, length, 0] + current_integral = CurrentIntegralAxisAligned(center=center, size=size, sign="+") + + # Compute impedance using the tool + Z_calc = ImpedanceCalculator( + voltage_integral=voltage_integral, current_integral=current_integral + ) + _ = Z_calc.compute_impedance(SIM_Z_DATA["field_time"]) + Z_calc = ImpedanceCalculator(voltage_integral=voltage_integral, current_integral=None) + _ = Z_calc.compute_impedance(SIM_Z_DATA["field_time"]) + Z_calc = ImpedanceCalculator(voltage_integral=None, current_integral=current_integral) + _ = Z_calc.compute_impedance(SIM_Z_DATA["field_time"]) + + def test_impedance_accuracy(): field_data = make_field_data() # Setup path integrals size = [0, STRIP_HEIGHT / 2, 0] center = [0, -STRIP_HEIGHT / 4, 0] - voltage_integral = VoltageIntegralAA(center=center, size=size, extrapolate_to_endpoints=True) + voltage_integral = VoltageIntegralAxisAligned( + center=center, size=size, sign="+", extrapolate_to_endpoints=True + ) size = [STRIP_WIDTH * 1.25, STRIP_HEIGHT / 2, 0] center = [0, 0, 0] - current_integral = CurrentIntegralAA(center=center, size=size) + current_integral = CurrentIntegralAxisAligned(center=center, size=size, sign="+") def impedance_of_stripline(width, height): # Assuming no fringing fields, is the same as a parallel plate diff --git a/tests/utils.py b/tests/utils.py index 1fe9932f5e..311cf4f120 100644 --- a/tests/utils.py +++ b/tests/utils.py @@ -10,6 +10,7 @@ from tidy3d.log import _get_level_int from tidy3d.web import BatchData from tidy3d.components.base import Tidy3dBaseModel +from tidy3d import ModeIndexDataArray """ utilities shared between all tests """ np.random.seed(4) @@ -78,7 +79,7 @@ def cartesian_to_unstructured( xyz = [array.x, array.y, array.z] lens = [len(coord) for coord in xyz] - num_len_zero = sum(l == 1 for l in lens) + num_len_zero = sum(length == 1 for length in lens) if num_len_zero == 1: normal_axis = lens.index(1) @@ -245,7 +246,7 @@ def cartesian_to_unstructured( def make_spatial_data( size, bounds, - lims=[0, 1], + lims=(0, 1), seed_data=None, unstructured=False, perturbation=0.1, @@ -875,6 +876,72 @@ def make_field_data(monitor: td.FieldMonitor) -> td.FieldData: **field_cmps, ) + def make_field_time_data(monitor: td.FieldTimeMonitor) -> td.FieldTimeData: + """make a random FieldTimeData from a FieldTimeMonitor.""" + field_cmps = {} + coords = {} + grid = simulation.discretize_monitor(monitor) + tmesh = simulation.tmesh + for field_name in monitor.fields: + spatial_coords_dict = grid[field_name].dict() + + for axis, dim in enumerate("xyz"): + if monitor.size[axis] == 0: + coords[dim] = [monitor.center[axis]] + else: + coords[dim] = np.array(spatial_coords_dict[dim]) + + (idx_begin, idx_end) = monitor.time_inds(tmesh) + tcoords = tmesh[idx_begin:idx_end] + coords["t"] = tcoords + field_cmps[field_name] = make_data( + coords=coords, data_array_type=td.ScalarFieldTimeDataArray, is_complex=False + ) + + return td.FieldTimeData( + monitor=monitor, + symmetry=(0, 0, 0), + symmetry_center=simulation.center, + grid_expanded=grid, + **field_cmps, + ) + + def make_mode_solver_data(monitor: td.ModeSolverMonitor) -> td.ModeSolverData: + """make a random ModeSolverData from a ModeSolverMonitor.""" + field_cmps = {} + coords = {} + grid = simulation.discretize_monitor(monitor) + index_coords = {} + index_coords["f"] = list(monitor.freqs) + index_coords["mode_index"] = np.arange(monitor.mode_spec.num_modes) + index_data_shape = (len(index_coords["f"]), len(index_coords["mode_index"])) + index_data = ModeIndexDataArray( + (1 + 1j) * np.random.random(index_data_shape), coords=index_coords + ) + for field_name in ["Ex", "Ey", "Ez", "Hx", "Hy", "Hz"]: + spatial_coords_dict = grid[field_name].dict() + + for axis, dim in enumerate("xyz"): + if monitor.size[axis] == 0: + coords[dim] = [monitor.center[axis]] + else: + coords[dim] = np.array(spatial_coords_dict[dim]) + + coords["f"] = list(monitor.freqs) + coords["mode_index"] = index_coords["mode_index"] + field_cmps[field_name] = make_data( + coords=coords, data_array_type=td.ScalarModeFieldDataArray, is_complex=True + ) + + return td.ModeSolverData( + monitor=monitor, + symmetry=(0, 0, 0), + symmetry_center=simulation.center, + grid_expanded=grid, + n_complex=index_data, + **field_cmps, + ) + def make_eps_data(monitor: td.PermittivityMonitor) -> td.PermittivityData: """make a random PermittivityData from a PermittivityMonitor.""" field_mnt = td.FieldMonitor(**monitor.dict(exclude={"type", "fields"})) @@ -920,6 +987,8 @@ def make_mode_data(monitor: td.ModeMonitor) -> td.ModeData: MONITOR_MAKER_MAP = { td.FieldMonitor: make_field_data, + td.FieldTimeMonitor: make_field_time_data, + td.ModeSolverMonitor: make_mode_solver_data, td.ModeMonitor: make_mode_data, td.PermittivityMonitor: make_eps_data, td.DiffractionMonitor: make_diff_data, diff --git a/tidy3d/plugins/microwave/__init__.py b/tidy3d/plugins/microwave/__init__.py index a7a6808ad2..5dbb3a1ae6 100644 --- a/tidy3d/plugins/microwave/__init__.py +++ b/tidy3d/plugins/microwave/__init__.py @@ -1,10 +1,15 @@ """ Imports from microwave plugin. """ -from .path_integrals import VoltageIntegralAA, CurrentIntegralAA +from .path_integrals import ( + AxisAlignedPathIntegral, + VoltageIntegralAxisAligned, + CurrentIntegralAxisAligned, +) from .impedance_calculator import ImpedanceCalculator __all__ = [ - "VoltageIntegralAA", - "CurrentIntegralAA", + "AxisAlignedPathIntegral", + "VoltageIntegralAxisAligned", + "CurrentIntegralAxisAligned", "ImpedanceCalculator", ] diff --git a/tidy3d/plugins/microwave/impedance_calculator.py b/tidy3d/plugins/microwave/impedance_calculator.py index 309db82a71..1a03da3269 100644 --- a/tidy3d/plugins/microwave/impedance_calculator.py +++ b/tidy3d/plugins/microwave/impedance_calculator.py @@ -6,40 +6,48 @@ import pydantic.v1 as pd import numpy as np from typing import Optional -from ...components.data.monitor_data import FieldData, FieldTimeData, ModeSolverData +from ...components.data.monitor_data import FieldTimeData from ...components.base import Tidy3dBaseModel -from ...exceptions import ValidationError +from ...exceptions import DataError, ValidationError -from .path_integrals import VoltageIntegralAA, CurrentIntegralAA +from .path_integrals import VoltageIntegralAxisAligned, CurrentIntegralAxisAligned +from .path_integrals import MonitorDataTypes, IntegralResultTypes class ImpedanceCalculator(Tidy3dBaseModel): """Tool for computing the characteristic impedance of a transmission line.""" - voltage_integral: Optional[VoltageIntegralAA] = pd.Field( + voltage_integral: Optional[VoltageIntegralAxisAligned] = pd.Field( ..., title="Voltage Integral", - description="Integral for computing voltage.", + description="Definition of path integral for computing voltage.", ) - current_integral: Optional[CurrentIntegralAA] = pd.Field( + current_integral: Optional[CurrentIntegralAxisAligned] = pd.Field( ..., - title="", - description="Integral for computing current.", + title="Current Integral", + description="Definition of contour integral for computing current.", ) - def compute_impedance(self, em_field: FieldData | ModeSolverData | FieldTimeData): + def compute_impedance(self, em_field: MonitorDataTypes) -> IntegralResultTypes: + """Compute impedance for the supplied ``em_field`` using ``voltage_integral`` and + ``current_integral``. If only a single integral has been defined, impedance is + computed using the total flux in ``em_field``.""" + if not isinstance(em_field, MonitorDataTypes): + raise DataError("'em_field' type not supported by impedance calculator.") + # If both voltage and current integrals have been defined then impedance is computed directly if self.voltage_integral: - voltage = self.voltage_integral.compute_voltage(em_field.field_components) + voltage = self.voltage_integral.compute_voltage(em_field) if self.current_integral: - current = self.current_integral.compute_current(em_field.field_components) + current = self.current_integral.compute_current(em_field) # If only one of the integrals has been provided then fall back to using total power (flux) - # with Ohm's law. The input field should cover an area large enough to render the flux computation accurate. - # If the input field is a time signal, then it is real and flux corresponds to the instantaneous power. - # Otherwise the input field is in frequency domain, where flux indicates the time-averaged power 0.5*Re(V*conj(I)) + # with Ohm's law. The input field should cover an area large enough to render the flux + # computation accurate. If the input field is a time signal, then it is real and flux + # corresponds to the instantaneous power. Otherwise the input field is in frequency domain, + # where flux indicates the time-averaged power 0.5*Re(V*conj(I)) if not self.voltage_integral: flux = em_field.flux if isinstance(em_field, FieldTimeData): @@ -53,11 +61,13 @@ def compute_impedance(self, em_field: FieldData | ModeSolverData | FieldTimeData else: current = np.conj(2 * flux / voltage) - return voltage / current + impedance = voltage / current + return impedance @pd.validator("current_integral", always=True) def check_voltage_or_current(cls, val, values): - """Ensure that 'voltage_integral' and/or 'current_integral' were provided.""" + """Raise validation error if both ``voltage_integral`` and ``current_integral`` + were not provided.""" if not values.get("voltage_integral") and not val: raise ValidationError( "Atleast one of 'voltage_integral' or 'current_integral' must be provided." diff --git a/tidy3d/plugins/microwave/path_integrals.py b/tidy3d/plugins/microwave/path_integrals.py index af75306115..b4f6a2462b 100644 --- a/tidy3d/plugins/microwave/path_integrals.py +++ b/tidy3d/plugins/microwave/path_integrals.py @@ -4,28 +4,39 @@ import pydantic.v1 as pd import numpy as np -from typing import Tuple -from ...components.data.dataset import FieldDataset, FieldTimeDataset, ModeSolverDataset -from ...components.data.data_array import AbstractSpatialDataArray -from ...components.base import cached_property +from abc import ABC, abstractmethod +from typing import List, Tuple, Union + +from ...components.data.monitor_data import FieldData, FieldTimeData, ModeSolverData +from ...components.data.data_array import ( + ScalarFieldDataArray, + ScalarFieldTimeDataArray, + ScalarModeFieldDataArray, +) +from ...components.data.data_array import FreqDataArray, TimeDataArray, FreqModeDataArray +from ...components.base import cached_property, Tidy3dBaseModel from ...components.types import Axis, Direction from ...components.geometry.base import Box from ...components.validators import assert_line, assert_plane -from ...exceptions import DataError +from ...exceptions import Tidy3dError, DataError +MonitorDataTypes = Union[FieldData, FieldTimeData, ModeSolverData] +EMScalarFieldType = Union[ScalarFieldDataArray, ScalarFieldTimeDataArray, ScalarModeFieldDataArray] +IntegralResultTypes = Union[FreqDataArray, FreqModeDataArray, TimeDataArray] -class AbstractAxesRH: + +class AbstractAxesRH(Tidy3dBaseModel, ABC): """Represents an axis-aligned right-handed coordinate system with one axis preferred.""" @cached_property + @abstractmethod def main_axis(self) -> Axis: - """Subclasses should implement this method.""" - raise NotImplementedError() + """Get the preferred axis.""" @cached_property def remaining_axes(self) -> Tuple[Axis, Axis]: - """Axes in plane, ordered to maintain a right-handed coordinate system""" - axes = [0, 1, 2] + """Get in-plane axes, ordered to maintain a right-handed coordinate system.""" + axes: List[Axis] = [0, 1, 2] axes.pop(self.main_axis) if self.main_axis == 1: return (axes[1], axes[0]) @@ -34,29 +45,30 @@ def remaining_axes(self) -> Tuple[Axis, Axis]: class AxisAlignedPathIntegral(AbstractAxesRH, Box): - """Class for defining the simplest type of path integral which is aligned with Cartesian axes.""" + """Class for defining the simplest type of path integral, which is aligned with Cartesian axes.""" _line_validator = assert_line() extrapolate_to_endpoints: bool = pd.Field( False, - title="Extrapolate to endpoints", + title="Extrapolate to Endpoints", description="If the endpoints of the path integral terminate at or near a material interface, " - "the field is likely discontinuous. This option ignores fields outside and on the bounds " - "of the integral. Should be turned on when computing voltage between two conductors. ", + "the field is likely discontinuous. When this field is ``True``, fields that are outside and on the bounds " + "of the integral are ignored. Should be enabled when computing voltage between two conductors.", ) snap_path_to_grid: bool = pd.Field( False, - title="Snap path to grid", + title="Snap Path to Grid", description="It might be desireable to integrate exactly along the Yee grid associated with " - "a field. If enabled, the integration path will be snapped to the grid.", + "a field. When this field is ``True``, the integration path will be snapped to the grid.", ) - def compute_integral(self, scalar_field: AbstractSpatialDataArray): - """Computes the defined integral given the input `scalar_field`.""" + def compute_integral(self, scalar_field: EMScalarFieldType) -> IntegralResultTypes: + """Computes the defined integral given the input ``scalar_field``.""" + if not scalar_field.does_cover(self.bounds): - raise DataError("scalar field does not cover the integration domain") + raise DataError("Scalar field does not cover the integration domain.") coord = "xyz"[self.main_axis] scalar_field = self._get_field_along_path(scalar_field) @@ -87,12 +99,17 @@ def compute_integral(self, scalar_field: AbstractSpatialDataArray): scalar_field = scalar_field.interp( coords_interp, method=method, kwargs={"fill_value": "extrapolate"} ) - return scalar_field.integrate(coord=coord) + result = scalar_field.integrate(coord=coord) + if isinstance(scalar_field, ScalarFieldDataArray): + return FreqDataArray(data=result.data, coords=result.coords) + elif isinstance(scalar_field, ScalarFieldTimeDataArray): + return TimeDataArray(data=result.data, coords=result.coords) + else: + assert isinstance(scalar_field, ScalarModeFieldDataArray) + return FreqModeDataArray(data=result.data, coords=result.coords) - def _get_field_along_path( - self, scalar_field: AbstractSpatialDataArray - ) -> AbstractSpatialDataArray: - """Returns a selection of the input `scalar_field` ready for integration.""" + def _get_field_along_path(self, scalar_field: EMScalarFieldType) -> EMScalarFieldType: + """Returns a selection of the input ``scalar_field`` ready for integration.""" axis1 = self.remaining_axes[0] axis2 = self.remaining_axes[1] coord1 = "xyz"[axis1] @@ -126,68 +143,74 @@ def _get_field_along_path( scalar_field = scalar_field.interp( coord2dict, method="linear", kwargs={"bounds_error": True} ) - + # Remove unneeded coordinates + scalar_field = scalar_field.reset_coords(drop=True) return scalar_field @cached_property def main_axis(self) -> Axis: """Axis for performing integration.""" - val = next((index for index, value in enumerate(self.size) if value != 0), None) - return val + for index, value in enumerate(self.size): + if value != 0: + return index + raise Tidy3dError("Failed to identify axis.") -class VoltageIntegralAA(AxisAlignedPathIntegral): +class VoltageIntegralAxisAligned(AxisAlignedPathIntegral): """Class for computing the voltage between two points defined by an axis-aligned line.""" sign: Direction = pd.Field( - "+", - title="Direction of path integral.", + ..., + title="Direction of Path Integral", description="Positive indicates V=Vb-Va where position b has a larger coordinate along the axis of integration.", ) - def compute_voltage(self, em_field: FieldDataset | ModeSolverDataset | FieldTimeDataset): + def compute_voltage(self, em_field: MonitorDataTypes) -> IntegralResultTypes: """Compute voltage along path defined by a line.""" + if not isinstance(em_field, MonitorDataTypes): + raise DataError("'em_field' type not supported.") e_component = "xyz"[self.main_axis] field_name = f"E{e_component}" # Validate that the field is present - if field_name not in em_field: - raise DataError(f"field_name '{field_name}' not found") - e_field = em_field[f"E{e_component}"] - integral = self.compute_integral(e_field) + if field_name not in em_field.field_components: + raise DataError(f"'field_name' '{field_name}' not found.") + e_field = em_field.field_components[field_name] + # V = -integral(E) + voltage = self.compute_integral(e_field) - voltage = -integral - if self.sign == "-": + if self.sign == "+": voltage *= -1 # Return data array of voltage while keeping coordinates of frequency|time|mode index return voltage -class CurrentIntegralAA(AbstractAxesRH, Box): - """Class for computing conduction current via Ampere's Circuital Law on an axis-aligned loop.""" +class CurrentIntegralAxisAligned(AbstractAxesRH, Box): + """Class for computing conduction current via Ampère's circuital law on an axis-aligned loop.""" _plane_validator = assert_plane() sign: Direction = pd.Field( - "+", - title="Direction of contour integral", + ..., + title="Direction of Contour Integral", description="Positive indicates current flowing in the positive normal axis direction.", ) extrapolate_to_endpoints: bool = pd.Field( False, - title="Extrapolate to endpoints", - description="This parameter is passed to `AxisAlignedPathIntegral` objects when computing the contour integral.", + title="Extrapolate to Endpoints", + description="This parameter is passed to ``AxisAlignedPathIntegral`` objects when computing the contour integral.", ) snap_contour_to_grid: bool = pd.Field( False, - title="Snap contour to grid", - description="This parameter is passed to `AxisAlignedPathIntegral` objects when computing the contour integral.", + title="Snap Contour to Grid", + description="This parameter is passed to ``AxisAlignedPathIntegral`` objects when computing the contour integral.", ) - def compute_current(self, em_field: FieldDataset | ModeSolverDataset | FieldTimeDataset): + def compute_current(self, em_field: MonitorDataTypes) -> IntegralResultTypes: """Compute current flowing in loop defined by the outer edge of a rectangle.""" - + if not isinstance(em_field, MonitorDataTypes): + raise DataError("'em_field' type not supported.") ax1 = self.remaining_axes[0] ax2 = self.remaining_axes[1] h_component = "xyz"[ax1] @@ -195,12 +218,12 @@ def compute_current(self, em_field: FieldDataset | ModeSolverDataset | FieldTime h_field_name = f"H{h_component}" v_field_name = f"H{v_component}" # Validate that fields are present - if h_field_name not in em_field: - raise DataError(f"field_name '{h_field_name}' not found") - if v_field_name not in em_field: - raise DataError(f"field_name '{v_field_name}' not found") - h_horizontal = em_field[f"H{h_component}"] - h_vertical = em_field[f"H{v_component}"] + if h_field_name not in em_field.field_components: + raise DataError(f"'field_name' '{h_field_name}' not found.") + if v_field_name not in em_field.field_components: + raise DataError(f"'field_name' '{v_field_name}' not found.") + h_horizontal = em_field.field_components[h_field_name] + h_vertical = em_field.field_components[v_field_name] # Decompose contour into path integrals (bottom, right, top, left) = self._to_path_integrals(h_horizontal, h_vertical) @@ -214,17 +237,20 @@ def compute_current(self, em_field: FieldDataset | ModeSolverDataset | FieldTime if self.sign == "-": current *= -1 - # Return data array of current while keeping coordinates of frequency|time|mode index - current = current.drop_vars((h_component, v_component)) + return current @cached_property def main_axis(self) -> Axis: """Axis normal to loop""" - val = next((index for index, value in enumerate(self.size) if value == 0), None) - return val + for index, value in enumerate(self.size): + if value == 0: + return index + raise Tidy3dError("Failed to identify axis.") def _to_path_integrals(self, h_horizontal, h_vertical) -> Tuple[AxisAlignedPathIntegral, ...]: + """Returns four ``AxisAlignedPathIntegral`` instances, which represent a contour + integral around the surface defined by ``self.size``.""" ax1 = self.remaining_axes[0] ax2 = self.remaining_axes[1] diff --git a/tidy3d/plugins/smatrix/component_modelers/terminal.py b/tidy3d/plugins/smatrix/component_modelers/terminal.py index e2faec7118..6d5aa88277 100644 --- a/tidy3d/plugins/smatrix/component_modelers/terminal.py +++ b/tidy3d/plugins/smatrix/component_modelers/terminal.py @@ -23,7 +23,7 @@ from .base import AbstractComponentModeler, FWIDTH_FRAC from ..ports.lumped import LumpedPortDataArray, LumpedPort -from ...microwave import VoltageIntegralAA, CurrentIntegralAA +from ...microwave import VoltageIntegralAxisAligned, CurrentIntegralAxisAligned class TerminalComponentModeler(AbstractComponentModeler): @@ -215,14 +215,14 @@ def port_voltage(port: LumpedPort, sim_data: SimulationData) -> xr.DataArray: size = list(port.size) size[port.current_axis] = 0 - voltage_integral = VoltageIntegralAA( + voltage_integral = VoltageIntegralAxisAligned( center=port.center, size=size, extrapolate_to_endpoints=True, snap_path_to_grid=True, sign="+", ) - voltage = voltage_integral.compute_voltage(field_data.field_components) + voltage = voltage_integral.compute_voltage(field_data) # Return data array of voltage with coordinates of frequency return voltage @@ -242,7 +242,8 @@ def port_current(port: LumpedPort, sim_data: SimulationData) -> xr.DataArray: # Get h field tangent to resistive sheet h_component = "xyz"[port.current_axis] inject_component = "xyz"[port.injection_axis] - field_data = sim_data[f"{port.name}_H{h_component}"].field_components + monitor_data = sim_data[f"{port.name}_H{h_component}"] + field_data = monitor_data.field_components h_field = field_data[f"H{h_component}"] # Coordinates as numpy array for h_field along curren and injection axis h_coords_along_current = h_field.coords[h_component].values @@ -290,14 +291,14 @@ def port_current(port: LumpedPort, sim_data: SimulationData) -> xr.DataArray: size[port.injection_axis] = dcap # H field is continuous at integral bounds, so extrapolation is turned off - I_integral = CurrentIntegralAA( + I_integral = CurrentIntegralAxisAligned( center=center, size=size, sign="+", extrapolate_to_endpoints=False, snap_contour_to_grid=True, ) - return I_integral.compute_current(field_data) + return I_integral.compute_current(monitor_data) def port_ab(port: LumpedPort, sim_data: SimulationData): """Helper to compute the port incident and reflected power waves."""