From e53ab60f0e09414fefe361fb4cc432067c844b7a Mon Sep 17 00:00:00 2001 From: Weiliang Jin Date: Fri, 5 Nov 2021 12:04:11 -0700 Subject: [PATCH] add multi-pole Drude medium --- docs/api.rst | 1 + tests/test_components.py | 18 +++++++- tidy3d/__init__.py | 11 ++++- tidy3d/components/__init__.py | 2 +- tidy3d/components/medium.py | 82 ++++++++++++++++++++++++++++++++++- 5 files changed, 109 insertions(+), 5 deletions(-) diff --git a/docs/api.rst b/docs/api.rst index 883328767d..566ecc8f22 100644 --- a/docs/api.rst +++ b/docs/api.rst @@ -117,6 +117,7 @@ Mediums Sellmeier Debye Lorentz + Drude Methods ------- diff --git a/tests/test_components.py b/tests/test_components.py index dafb5e95e5..350114a8ba 100644 --- a/tests/test_components.py +++ b/tests/test_components.py @@ -218,10 +218,11 @@ def test_medium_dispersion(): m_PR = PoleResidue(eps_inf=1.0, poles=[(1 + 2j, 1 + 3j), (2 + 4j, 1 + 5j)]) m_SM = Sellmeier(coeffs=[(2, 3), (2, 4)]) m_LZ = Lorentz(eps_inf=1.0, coeffs=[(1, 3, 2), (2, 4, 1)]) + m_DR = Drude(eps_inf=1.0, coeffs=[(1, 3), (2, 4)]) m_DB = Debye(eps_inf=1.0, coeffs=[(1, 3), (2, 4)]) freqs = np.linspace(0.01, 1, 1001) - for medium in [m_PR, m_SM, m_LZ, m_DB]: + for medium in [m_PR, m_SM, m_LZ, m_DR, m_DB]: eps_c = medium.eps_model(freqs) @@ -230,15 +231,28 @@ def test_medium_dispersion_conversion(): m_PR = PoleResidue(eps_inf=1.0, poles=[((1 + 2j), (1 + 3j)), ((2 + 4j), (1 + 5j))]) m_SM = Sellmeier(coeffs=[(2, 3), (2, 4)]) m_LZ = Lorentz(eps_inf=1.0, coeffs=[(1, 3, 2), (2, 4, 1)]) + m_DR = Drude(eps_inf=1.0, coeffs=[(1, 3), (2, 4)]) m_DB = Debye(eps_inf=1.0, coeffs=[(1, 3), (2, 4)]) freqs = np.linspace(0.01, 1, 1001) - for medium in [m_PR, m_SM, m_DB, m_LZ]: # , m_DB]: + for medium in [m_PR, m_SM, m_DB, m_LZ, m_DR]: # , m_DB]: eps_model = medium.eps_model(freqs) eps_pr = medium.pole_residue.eps_model(freqs) np.testing.assert_allclose(eps_model, eps_pr) +def test_medium_dispersion_create(): + + m_PR = PoleResidue(eps_inf=1.0, poles=[((1 + 2j), (1 + 3j)), ((2 + 4j), (1 + 5j))]) + m_SM = Sellmeier(coeffs=[(2, 3), (2, 4)]) + m_LZ = Lorentz(eps_inf=1.0, coeffs=[(1, 3, 2), (2, 4, 1)]) + m_DR = Drude(eps_inf=1.0, coeffs=[(1, 3), (2, 4)]) + m_DB = Debye(eps_inf=1.0, coeffs=[(1, 3), (2, 4)]) + + for medium in [m_PR, m_SM, m_DB, m_LZ, m_DR]: + struct = Structure(geometry=Box(size=(1, 1, 1)), medium=medium) + + """ modes """ diff --git a/tidy3d/__init__.py b/tidy3d/__init__.py index 123b5da297..48e573334e 100644 --- a/tidy3d/__init__.py +++ b/tidy3d/__init__.py @@ -17,7 +17,16 @@ from .components import Box, Sphere, Cylinder, PolySlab # medium -from .components import Medium, PoleResidue, Sellmeier, Debye, Lorentz, AnisotropicMedium, PEC +from .components import ( + Medium, + PoleResidue, + Sellmeier, + Debye, + Drude, + Lorentz, + AnisotropicMedium, + PEC, +) from .components import nk_to_eps_complex, nk_to_eps_sigma, eps_complex_to_nk from .components import nk_to_medium, eps_sigma_to_eps_complex diff --git a/tidy3d/components/__init__.py b/tidy3d/components/__init__.py index 420c7d68b0..c65ffbde7f 100644 --- a/tidy3d/components/__init__.py +++ b/tidy3d/components/__init__.py @@ -13,7 +13,7 @@ from .geometry import Geometry # medium -from .medium import Medium, PoleResidue, Sellmeier, Debye, Lorentz, AnisotropicMedium, PEC +from .medium import Medium, PoleResidue, Sellmeier, Debye, Drude, Lorentz, AnisotropicMedium, PEC from .medium import nk_to_eps_complex, nk_to_eps_sigma, eps_complex_to_nk from .medium import nk_to_medium, eps_sigma_to_eps_complex diff --git a/tidy3d/components/medium.py b/tidy3d/components/medium.py index 3b5f27a2a1..e48817595f 100644 --- a/tidy3d/components/medium.py +++ b/tidy3d/components/medium.py @@ -511,6 +511,84 @@ def pole_residue(self): ) +class Drude(DispersiveMedium): + """A dispersive medium described by the Drude model. + The frequency-dependence of the complex-valued permittivity is described by: + + .. math:: + \\epsilon(f) = \\epsilon_\\infty - \\sum_i + \\frac{ f_i^2}{f^2 - jf\\delta_i} + + where :math:`f, f_i, \\delta_i` are in Hz. + + Parameters + ---------- + eps_inf : float = 1.0 + Relative permittivity at infinite frequency (:math:`\\epsilon_\\infty`). + coeffs : List[Tuple[float, float]] + List of (:math:`f_i, \\delta_i`) values for model. + frequeuncy_range : Tuple[float, float] = (-inf, inf) + Range of validity for the medium in Hz. + If simulation or plotting functions use frequency out of this range, a warning is thrown. + name : str = None + Optional name for the medium. + + Example + ------- + >>> drude_medium = Drude(eps_inf=2.0, coeffs=[(1,2), (3,4)]) + >>> eps = drude_medium.eps_model(200e12) + """ + + eps_inf: float = 1.0 + coeffs: List[Tuple[float, float]] + type: Literal["Drude"] = "Drude" + + @ensure_freq_in_range + def eps_model(self, frequency: float) -> complex: + """Complex-valued permittivity as a function of frequency. + + Parameters + ---------- + frequency : float + Frequency to evaluate permittivity at (Hz). + + Returns + ------- + complex + Complex-valued relative permittivity evaluated at the frequency. + """ + eps = self.eps_inf + 0.0j + for (f, delta) in self.coeffs: + eps -= (f ** 2) / (frequency ** 2 - 1j * frequency * delta) + return eps + + @property + def pole_residue(self): + """Representation of Medium as a pole-residue model.""" + + poles = [] + for (f, delta) in self.coeffs: + + w = 2 * np.pi * f + d = 2 * np.pi * delta + + c0 = -(w ** 2) / 2 / d + 0j + a0 = 0j + + c1 = w ** 2 / 2 / d + 0j + a1 = d + 0j + + poles.append((a0, c0)) + poles.append((a1, c1)) + + return PoleResidue( + eps_inf=self.eps_inf, + poles=poles, + frequency_range=self.frequency_range, + name=self.name, + ) + + class Debye(DispersiveMedium): """A dispersive medium described by the Debye model. The frequency-dependence of the complex-valued permittivity is described by: @@ -581,7 +659,9 @@ def pole_residue(self): # types of mediums that can be used in Simulation and Structures -MediumType = Union[Literal[PEC], Medium, AnisotropicMedium, PoleResidue, Sellmeier, Lorentz, Debye] +MediumType = Union[ + Literal[PEC], Medium, AnisotropicMedium, PoleResidue, Sellmeier, Lorentz, Debye, Drude +] """ Conversion helper functions """