diff --git a/docs/index.rst b/docs/index.rst index 37ae8ad..49caa70 100644 --- a/docs/index.rst +++ b/docs/index.rst @@ -36,6 +36,7 @@ Features * Spherical, conical, and toroidal surface sag profiles * Circular, rectangular, and polygonal apertures * Support for mirrors and arbitrary multilayer coatings +* Refractive glass materials with Sellmeier dispersion (e.g. N-BK7, F2) * Diffraction grating support * Constant, polynomial and holographic ruling spacing @@ -54,10 +55,11 @@ Limitations system. * **Physical Optics**. Only geometric optics is supported right now, but adding a Fourier optics propagator is a longstanding goal of the project. -* **Glass Optical Constants**. :mod:`optika` has a wide array of optical +* **Glass Catalog**. :mod:`optika` has a wide array of optical constants from sources such as :cite:t:`Palik1997` and :cite:t:`Henke1993`, - but it does not yet have a database for different types of glass like Zemax - does. + and the :class:`~optika.materials.Glass` material provides Sellmeier dispersion + for a few common glasses (e.g. N-BK7, F2), but it does not yet have a + comprehensive glass database like Zemax does. Differences from Zemax ---------------------- diff --git a/optika/materials/__init__.py b/optika/materials/__init__.py index 1c4171d..2fdc2bb 100644 --- a/optika/materials/__init__.py +++ b/optika/materials/__init__.py @@ -21,6 +21,7 @@ AbstractMirror, Mirror, MeasuredMirror, + Glass, ) from ._multilayers import ( multilayer_coefficients, @@ -53,6 +54,7 @@ "AbstractMirror", "Mirror", "MeasuredMirror", + "Glass", "multilayer_coefficients", "multilayer_efficiency", "layer_absorbance", diff --git a/optika/materials/_materials.py b/optika/materials/_materials.py index f65d325..2c3a38a 100644 --- a/optika/materials/_materials.py +++ b/optika/materials/_materials.py @@ -1,6 +1,7 @@ from __future__ import annotations import abc import dataclasses +import numpy as np import astropy.units as u import named_arrays as na import optika @@ -12,6 +13,7 @@ "AbstractMirror", "Mirror", "MeasuredMirror", + "Glass", ] @@ -301,3 +303,153 @@ def efficiency( xp=wavelength, fp=efficiency, ) + + +@dataclasses.dataclass(eq=False, repr=False) +class Glass( + AbstractMaterial, +): + r""" + A transparent, refractive material whose index of refraction follows the + three-term Sellmeier dispersion equation, + + .. math:: + + n^2(\lambda) = 1 + + \frac{B_1 \lambda^2}{\lambda^2 - C_1} + + \frac{B_2 \lambda^2}{\lambda^2 - C_2} + + \frac{B_3 \lambda^2}{\lambda^2 - C_3}, + + where :math:`\lambda` is the vacuum wavelength of the light, and + :math:`B_i` (dimensionless) and :math:`C_i` (square length) are the + Sellmeier coefficients of the glass. + + Unlike :class:`Vacuum` and :class:`Mirror`, this material changes the index + of refraction of a transmitted ray, so a curved surface made of it has + optical power and bends light according to Snell's law. + + Examples + -------- + + Plot the index of refraction of N-BK7 and F2 glass across the visible + spectrum. + + .. jupyter-execute:: + + import matplotlib.pyplot as plt + import astropy.units as u + import named_arrays as na + import optika + + wavelength = na.linspace(380, 750, axis="wavelength", num=101) * u.nm + + glasses = { + "N-BK7": optika.materials.Glass.n_bk7(), + "F2": optika.materials.Glass.f2(), + } + + fig, ax = plt.subplots(constrained_layout=True) + for name, glass in glasses.items(): + rays = optika.rays.RayVectorArray(wavelength=wavelength) + na.plt.plot( + wavelength, + glass.index_refraction(rays), + ax=ax, + label=name, + ) + ax.set_xlabel(f"wavelength ({wavelength.unit:latex_inline})"); + ax.set_ylabel("index of refraction"); + ax.legend(); + """ + + b1: float | na.AbstractScalar = 0 + """The first dimensionless Sellmeier coefficient.""" + + b2: float | na.AbstractScalar = 0 + """The second dimensionless Sellmeier coefficient.""" + + b3: float | na.AbstractScalar = 0 + """The third dimensionless Sellmeier coefficient.""" + + c1: u.Quantity | na.AbstractScalar = 0 * u.um**2 + """The first Sellmeier resonance (units of square length).""" + + c2: u.Quantity | na.AbstractScalar = 0 * u.um**2 + """The second Sellmeier resonance (units of square length).""" + + c3: u.Quantity | na.AbstractScalar = 0 * u.um**2 + """The third Sellmeier resonance (units of square length).""" + + @classmethod + def n_bk7(cls) -> Glass: + """ + SCHOTT N-BK7 borosilicate crown glass + (:math:`n_d \\approx 1.5168`, :math:`V_d \\approx 64.2`). + """ + return cls( + b1=1.03961212, + b2=0.231792344, + b3=1.01046945, + c1=0.00600069867 * u.um**2, + c2=0.0200179144 * u.um**2, + c3=103.560653 * u.um**2, + ) + + @classmethod + def f2(cls) -> Glass: + """ + SCHOTT F2 flint glass + (:math:`n_d \\approx 1.6200`, :math:`V_d \\approx 36.4`). + """ + return cls( + b1=1.34533359, + b2=0.209073176, + b3=0.937357162, + c1=0.00997743871 * u.um**2, + c2=0.0470450767 * u.um**2, + c3=111.886764 * u.um**2, + ) + + @property + def shape(self) -> dict[str, int]: + return na.broadcast_shapes( + optika.shape(self.b1), + optika.shape(self.b2), + optika.shape(self.b3), + optika.shape(self.c1), + optika.shape(self.c2), + optika.shape(self.c3), + ) + + @property + def transformation(self) -> None: + return None + + def index_refraction( + self, + rays: optika.rays.RayVectorArray, + ) -> na.ScalarLike: + w2 = np.square(rays.wavelength) + n2 = 1 + ( + self.b1 * w2 / (w2 - self.c1) + + self.b2 * w2 / (w2 - self.c2) + + self.b3 * w2 / (w2 - self.c3) + ) + return np.sqrt(n2) + + def attenuation( + self, + rays: optika.rays.RayVectorArray, + ) -> na.ScalarLike: + return 0 / u.mm + + def efficiency( + self, + rays: optika.rays.RayVectorArray, + normal: na.AbstractCartesian3dVectorArray, + ) -> na.ScalarLike: + return 1 + + @property + def is_mirror(self) -> bool: + return False diff --git a/optika/materials/_tests/test_materials.py b/optika/materials/_tests/test_materials.py index bd746fb..dd0241e 100644 --- a/optika/materials/_tests/test_materials.py +++ b/optika/materials/_tests/test_materials.py @@ -107,3 +107,44 @@ class TestMeasuredMirror( AbstractTestAbstractMirror, ): pass + + +@pytest.mark.parametrize( + argnames="a", + argvalues=[ + optika.materials.Glass(), + optika.materials.Glass.n_bk7(), + optika.materials.Glass.f2(), + ], +) +class TestGlass( + AbstractTestAbstractMaterial, +): + pass + + +@pytest.mark.parametrize( + argnames="glass,n_d", + argvalues=[ + (optika.materials.Glass.n_bk7(), 1.5168), + (optika.materials.Glass.f2(), 1.6200), + ], +) +def test_glass_dispersion( + glass: optika.materials.Glass, + n_d: float, +): + # index of refraction at the helium d Fraunhofer line should match the + # published value of the glass. + rays_d = optika.rays.RayVectorArray(wavelength=587.5618 * u.nm) + n = glass.index_refraction(rays_d) + assert np.isclose(float(n), n_d, atol=1e-3) + + # the glass must be dispersive: a higher index toward the blue end of the + # spectrum (normal dispersion). + rays_F = optika.rays.RayVectorArray(wavelength=486.1327 * u.nm) + rays_C = optika.rays.RayVectorArray(wavelength=656.2725 * u.nm) + assert glass.index_refraction(rays_F) > glass.index_refraction(rays_C) + + # a glass transmits rather than reflects. + assert not glass.is_mirror