Enable bulk PMC material

CHANGELOG.md

@@ -11,6 +11,7 @@ and this project adheres to [Semantic Versioning](https://semver.org/spec/v2.0.0
 
 - Added `eps_lim` keyword argument to `Simulation.plot_eps()` for manual control over the permittivity color limits.
 - Added `thickness` parameter to `LossyMetalMedium` for computing surface impedance of a thin conductor.
+- Added material type `PMCMedium` for perfect magnetic conductor.
 
 ### Changed
 - Relaxed bounds checking of path integrals during `WavePort` validation.

docs/api/mediums.rst

@@ -16,6 +16,7 @@ Spatially uniform
    tidy3d.Medium
    tidy3d.LossyMetalMedium
    tidy3d.PECMedium
+   tidy3d.PMCMedium
    tidy3d.FullyAnisotropicMedium
 
 Spatially varying

tests/test_components/material/test_multi_physics.py

@@ -14,6 +14,7 @@ def test_delegated_attributes_work(dummy_optical):
 
     # delegated names resolve
     assert mp.is_pec is dummy_optical.is_pec
+    assert mp.is_pmc is dummy_optical.is_pmc
     assert mp._eps_plot == dummy_optical._eps_plot
     assert mp.viz_spec == dummy_optical.viz_spec
 
@@ -27,6 +28,9 @@ def test_delegated_attribute_without_optical_raises():
     with pytest.raises(AttributeError, match=r"optical medium is 'None'"):
         _ = mp_no_opt.is_pec
 
+    with pytest.raises(AttributeError, match=r"optical medium is 'None'"):
+        _ = mp_no_opt.is_pmc
+
 
 def test_has_cached_props(dummy_optical):
     mp = td.MultiPhysicsMedium(optical=dummy_optical)

tests/test_components/test_medium.py

@@ -14,12 +14,13 @@
 MEDIUM = td.Medium()
 ANIS_MEDIUM = td.AnisotropicMedium(xx=MEDIUM, yy=MEDIUM, zz=MEDIUM)
 PEC = td.PECMedium()
+PMC = td.PMCMedium()
 PR = td.PoleResidue(poles=[(-1 + 1j, 2 + 2j)])
 SM = td.Sellmeier(coeffs=[(1, 2)])
 LZ = td.Lorentz(coeffs=[(1, 2, 3)])
 DR = td.Drude(coeffs=[(1, 2)])
 DB = td.Debye(coeffs=[(1, 2)])
-MEDIUMS = [MEDIUM, ANIS_MEDIUM, PEC, PR, SM, LZ, DR, DB]
+MEDIUMS = [MEDIUM, ANIS_MEDIUM, PEC, PR, SM, LZ, DR, DB, PMC]
 
 f, AX = plt.subplots()
 
@@ -141,6 +142,10 @@ def test_PEC():
     _ = td.Structure(geometry=td.Box(size=(1, 1, 1)), medium=td.PEC)
 
 
+def test_PMC():
+    _ = td.Structure(geometry=td.Box(size=(1, 1, 1)), medium=td.PMC)
+
+
 def test_lossy_metal():
     # frequency_range shouldn't be None
     with pytest.raises(pydantic.ValidationError):

tests/test_components/test_time_modulation.py

@@ -212,6 +212,10 @@ def test_unsupported_modulated_medium_types():
     with pytest.raises(pydantic.ValidationError):
         td.PECMedium(modulation_spec=modulation_spec)
 
+    # PMC cannot be modulated
+    with pytest.raises(pydantic.ValidationError):
+        td.PMCMedium(modulation_spec=modulation_spec)
+
     # For Anisotropic medium, one should modulate the components, not the whole medium
     with pytest.raises(pydantic.ValidationError):
         td.AnisotropicMedium(

tests/utils.py

@@ -471,6 +471,10 @@ def make_custom_data(lims, unstructured):
             geometry=td.Box(size=(1, 1, 1), center=(-1, 0, 0)),
             medium=td.AnisotropicMedium(xx=td.PEC, yy=td.Medium(), zz=td.Medium()),
         ),
+        td.Structure(
+            geometry=td.Box(size=(1, 1, 1), center=(-1, 0, 0)),
+            medium=td.AnisotropicMedium(xx=td.PMC, yy=td.Medium(), zz=td.Medium()),
+        ),
         # Test a fully anistropic medium
         td.Structure(
             geometry=td.Box(size=(1, 1, 1), center=(-1, 0, 0)),

tidy3d/__init__.py

@@ -235,6 +235,7 @@
 from .components.medium import (
     PEC,
     PEC2D,
+    PMC,
     AbstractMedium,
     AnisotropicMedium,
     CustomAnisotropicMedium,
@@ -260,6 +261,7 @@
     PECMedium,
     PerturbationMedium,
     PerturbationPoleResidue,
+    PMCMedium,
     PoleResidue,
     Sellmeier,
     SurfaceImpedanceFitterParam,
@@ -419,7 +421,9 @@ def set_logging_level(level: str) -> None:
     "PoleResidue",
     "AnisotropicMedium",
     "PEC",
+    "PMC",
     "PECMedium",
+    "PMCMedium",
     "Medium2D",
     "PEC2D",
     "Sellmeier",

tidy3d/components/material/multi_physics.py

@@ -140,6 +140,7 @@ def __getattr__(self, name: str):
 
         DELEGATED_ATTRIBUTES = {
             "is_pec": self.optical,
+            "is_pmc": self.optical,
             "_eps_plot": self.optical,
             "viz_spec": self.optical,
         }

tidy3d/components/medium.py

@@ -1355,6 +1355,11 @@ def is_pec(self):
         """Whether the medium is a PEC."""
         return False
 
+    @cached_property
+    def is_pmc(self):
+        """Whether the medium is a PMC."""
+        return False
+
     def sel_inside(self, bounds: Bound) -> AbstractMedium:
         """Return a new medium that contains the minimal amount data necessary to cover
         a spatial region defined by ``bounds``.
@@ -1769,6 +1774,52 @@ def is_pec(self):
 PEC = PECMedium(name="PEC")
 
 
+# PMC keyword
+class PMCMedium(AbstractMedium):
+    """Perfect magnetic conductor class.
+
+    Note
+    ----
+
+        To avoid confusion from duplicate PMCs, must import ``tidy3d.PMC`` instance directly.
+
+
+
+    """
+
+    @pd.validator("modulation_spec", always=True)
+    def _validate_modulation_spec(cls, val):
+        """Check compatibility with modulation_spec."""
+        if val is not None:
+            raise ValidationError(
+                f"A 'modulation_spec' of class {type(val)} is not "
+                f"currently supported for medium class {cls}."
+            )
+        return val
+
+    @ensure_freq_in_range
+    def eps_model(self, frequency: float) -> complex:
+        # permittivity of a PMC.
+        return 1.0 + 0j
+
+    @cached_property
+    def n_cfl(self):
+        """This property computes the index of refraction related to CFL condition, so that
+        the FDTD with this medium is stable when the time step size that doesn't take
+        material factor into account is multiplied by ``n_cfl``.
+        """
+        return 1.0
+
+    @cached_property
+    def is_pmc(self):
+        """Whether the medium is a PMC."""
+        return True
+
+
+# PEC builtin instance
+PMC = PMCMedium(name="PMC")
+
+
 class Medium(AbstractMedium):
     """Dispersionless medium. Mediums define the optical properties of the materials within the simulation.
 
@@ -5681,6 +5732,9 @@ def plot(
 
 
 IsotropicUniformMediumType = Union[
+    Medium, LossyMetalMedium, PoleResidue, Sellmeier, Lorentz, Debye, Drude, PECMedium, PMCMedium
+]
+IsotropicUniformMediumFor2DType = Union[
     Medium, LossyMetalMedium, PoleResidue, Sellmeier, Lorentz, Debye, Drude, PECMedium
 ]
 IsotropicCustomMediumType = Union[
@@ -5890,10 +5944,19 @@ def is_pec(self):
         """Whether the medium is a PEC."""
         return any(self.is_comp_pec(i) for i in range(3))
 
+    @cached_property
+    def is_pmc(self):
+        """Whether the medium is a PMC."""
+        return any(self.is_comp_pmc(i) for i in range(3))
+
     def is_comp_pec(self, comp: Axis):
         """Whether the medium is a PEC."""
         return isinstance(self.components[["xx", "yy", "zz"][comp]], PECMedium)
 
+    def is_comp_pmc(self, comp: Axis):
+        """Whether the medium is a PMC."""
+        return isinstance(self.components[["xx", "yy", "zz"][comp]], PMCMedium)
+
     def sel_inside(self, bounds: Bound):
         """Return a new medium that contains the minimal amount data necessary to cover
         a spatial region defined by ``bounds``.
@@ -7007,6 +7070,7 @@ def perturbed_copy(
     Medium,
     AnisotropicMedium,
     PECMedium,
+    PMCMedium,
     PoleResidue,
     Sellmeier,
     Lorentz,
@@ -7041,7 +7105,7 @@ class Medium2D(AbstractMedium):
 
     """
 
-    ss: IsotropicUniformMediumType = pd.Field(
+    ss: IsotropicUniformMediumFor2DType = pd.Field(
         ...,
         title="SS Component",
         description="Medium describing the ss-component of the diagonal permittivity tensor. "
@@ -7052,7 +7116,7 @@ class Medium2D(AbstractMedium):
         discriminator=TYPE_TAG_STR,
     )
 
-    tt: IsotropicUniformMediumType = pd.Field(
+    tt: IsotropicUniformMediumFor2DType = pd.Field(
         ...,
         title="TT Component",
         description="Medium describing the tt-component of the diagonal permittivity tensor. "
@@ -7086,7 +7150,7 @@ def _validate_inplane_pec(cls, val, values):
 
     @classmethod
     def _weighted_avg(
-        cls, meds: List[IsotropicUniformMediumType], weights: List[float]
+        cls, meds: List[IsotropicUniformMediumFor2DType], weights: List[float]
     ) -> Union[PoleResidue, PECMedium]:
         """Average ``meds`` with weights ``weights``."""
         eps_inf = 1
@@ -7140,7 +7204,7 @@ def volumetric_equivalent(
             The 3D material corresponding to this 2D material.
         """
 
-        def get_component(med: MediumType3D, comp: Axis) -> IsotropicUniformMediumType:
+        def get_component(med: MediumType3D, comp: Axis) -> IsotropicUniformMediumFor2DType:
             """Extract the ``comp`` component of ``med``."""
             if isinstance(med, AnisotropicMedium):
                 dim = "xyz"[comp]
@@ -7402,7 +7466,7 @@ def sigma_model(self, freq: float) -> complex:
         return np.mean([self.ss.sigma_model(freq), self.tt.sigma_model(freq)], axis=0)
 
     @property
-    def elements(self) -> Dict[str, IsotropicUniformMediumType]:
+    def elements(self) -> Dict[str, IsotropicUniformMediumFor2DType]:
         """The diagonal elements of the 2D medium as a dictionary."""
         return dict(ss=self.ss, tt=self.tt)
 

tidy3d/components/mode/mode_solver.py

@@ -1738,6 +1738,8 @@ def _contain_good_conductor(self) -> bool:
         for medium in sim.scene.mediums:
             if medium.is_pec:
                 return True
+            if medium.is_pmc:
+                return True
             if apply_sibc and isinstance(medium, LossyMetalMedium):
                 return True
         return False

tidy3d/components/scene.py

@@ -495,9 +495,11 @@ def _get_structure_plot_params(
 
         if isinstance(medium, MultiPhysicsMedium):
             is_pec = medium.optical is not None and medium.optical.is_pec
+            is_pmc = medium.optical is not None and medium.optical.is_pmc
             is_time_modulated = medium.optical is not None and medium.optical.is_time_modulated
         else:
             is_pec = medium.is_pec
+            is_pmc = medium.is_pmc
             is_time_modulated = medium.is_time_modulated
 
         if mat_index == 0 or medium == self.medium:
@@ -508,6 +510,11 @@ def _get_structure_plot_params(
             plot_params = plot_params.copy(
                 update={"facecolor": "gold", "edgecolor": "k", "linewidth": 1}
             )
+        elif is_pmc:
+            # perfect magnetic conductor
+            plot_params = plot_params.copy(
+                update={"facecolor": "purple", "edgecolor": "k", "linewidth": 1}
+            )
         elif is_time_modulated:
             # time modulated medium
             plot_params = plot_params.copy(
@@ -1259,6 +1266,11 @@ def _get_structure_eps_plot_params(
             plot_params = plot_params.copy(
                 update={"facecolor": "gold", "edgecolor": "k", "linewidth": 1}
             )
+        elif medium.is_pmc:
+            # perfect magnetic conductor
+            plot_params = plot_params.copy(
+                update={"facecolor": "purple", "edgecolor": "k", "linewidth": 1}
+            )
         elif isinstance(medium, Medium2D):
             # 2d material
             plot_params = plot_params.copy(update={"edgecolor": "k", "linewidth": 1})

tidy3d/components/simulation.py

@@ -3487,8 +3487,8 @@ def _warn_grid_size_too_small(cls, val, values):
                 freq0 = source.source_time.freq0
 
                 for medium_index, medium in enumerate(mediums):
-                    # min wavelength in PEC is meaningless and we'll get divide by inf errors
-                    if medium.is_pec:
+                    # min wavelength in PEC/PMC is meaningless and we'll get divide by inf errors
+                    if medium.is_pec or medium.is_pmc:
                         continue
                     # min wavelength in Medium2D is meaningless
                     if isinstance(medium, Medium2D):
@@ -3500,8 +3500,10 @@ def _warn_grid_size_too_small(cls, val, values):
                     for comp, (key, grid_spec) in enumerate(
                         zip("xyz", (val.grid_x, val.grid_y, val.grid_z))
                     ):
-                        if medium.is_pec or (
-                            isinstance(medium, AnisotropicMedium) and medium.is_comp_pec(comp)
+                        if (
+                            medium.is_pec
+                            or medium.is_pmc
+                            or (isinstance(medium, AnisotropicMedium) and medium.is_comp_pec(comp))
                         ):
                             n_material = 1.0
                         lambda_min = C_0 / freq0 / n_material

tidy3d/components/subpixel_spec.py

@@ -129,6 +129,7 @@ def courant_ratio(self) -> float:
 
 
 PECSubpixelType = Union[Staircasing, HeuristicPECStaircasing, PECConformal]
+PMCSubpixelType = Union[Staircasing, HeuristicPECStaircasing]
 
 
 class SurfaceImpedance(PECConformal):
@@ -176,6 +177,13 @@ class SubpixelSpec(Tidy3dBaseModel):
         discriminator=TYPE_TAG_STR,
     )
 
+    pmc: PMCSubpixelType = pd.Field(
+        Staircasing(),
+        title="Subpixel Averaging Method For PMC Interfaces",
+        description="Subpixel averaging method applied to PMC structure interfaces.",
+        discriminator=TYPE_TAG_STR,
+    )
+
     lossy_metal: LossyMetalSubpixelType = pd.Field(
         SurfaceImpedance(),
         title="Subpixel Averaging Method for Lossy Metal Interfaces",
@@ -190,6 +198,7 @@ def staircasing(cls) -> SubpixelSpec:
             dielectric=Staircasing(),
             metal=Staircasing(),
             pec=Staircasing(),
+            pmc=Staircasing(),
             lossy_metal=Staircasing(),
         )