Materials Module

Materials Module

Provides optical material definitions with spatially-varying properties. Supports homogeneous materials with constant properties and dispersion models.

Module Structure

base/: Abstract base classes (MaterialField, HomogeneousMaterial, etc.)
implementations/: Concrete material implementations
utils/: Dispersion functions, factories, constants

Note: GPU propagation kernels are in lsurf.propagation.kernels.propagation

Classes

MaterialField: Abstract base class for material property fields.
HomogeneousMaterial: Uniform material with constant optical properties.
SimpleInhomogeneousModel: Radially-symmetric inhomogeneous materials (GPU via LUT).
GridInhomogeneousModel: 3D grid-based inhomogeneous materials (GPU via trilinear interp).
FullInhomogeneousModel: Arbitrary Python inhomogeneous materials (CPU-only).

Constants

VACUUM: Perfect vacuum with n = 1.0.
AIR_STP: Air at standard temperature and pressure, n = 1.000293.
WATER: Pure water with n = 1.333 and absorption.
BK7_GLASS: Schott N-BK7 optical glass, n = 1.5168.

Functions

sellmeier_refractive_index: Compute n using Sellmeier equation.
cauchy_refractive_index: Compute n using Cauchy equation.
create_sellmeier_material: Factory for Sellmeier-dispersion materials.
create_cauchy_material: Factory for Cauchy-dispersion materials.

Examples

>>> from lsurf.materials import VACUUM, WATER, HomogeneousMaterial
>>>
>>> # Use pre-defined material
>>> n = WATER.get_refractive_index(0, 0, 0, 550e-9)
>>>
>>> # Create custom material
>>> oil = HomogeneousMaterial("Oil", 1.47, absorption_coef=0.01)

class lsurf.materials.MaterialField(name='Material', kernel=None, propagator=None)[source]

Bases: ABC

Abstract base class for spatially-varying material properties.

All material properties are functions of position (x, y, z) and wavelength. Subclasses must implement methods that can be called from both CPU and GPU.

Parameters:: name (str, optional) – Descriptive name for this material. Default is “Material”.

name

Material identifier.

Type:: str

Notes

Material properties modeled: - Refractive index n: Real part of complex refractive index - Absorption coefficient α: Controls intensity decay - Scattering coefficient μ_s: Controls scattering mean free path - Anisotropy factor g: Henyey-Greenstein scattering parameter

References

__init__(name='Material', kernel=None, propagator=None)[source]

Initialize material field.

Parameters:

name (str, optional) – Descriptive name for this material. Default is “Material”.
kernel (PropagationKernelID, optional) – Override the default propagation kernel. Must be in supported_kernels(). If None, uses the class default.
propagator (PropagatorID, optional) – Override the default propagator. Must be in supported_propagators(). If None, uses the class default.

Raises:

ValueError – If kernel or propagator is not supported by this material type.

abstractmethod get_refractive_index(x, y, z, wavelength)[source]

Get refractive index at position (x, y, z) for given wavelength.

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

n – Real part of refractive index (dimensionless).

Return type:

float or ndarray

Notes

For absorbing materials, this returns only Re(n). Use get_absorption_coefficient() for the imaginary part.

abstractmethod get_refractive_index_gradient(x, y, z, wavelength)[source]

Get gradient of refractive index ∇n at position (x, y, z).

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

grad_n – (∂n/∂x, ∂n/∂y, ∂n/∂z) in units of m⁻¹.

Return type:

tuple of (float or ndarray)

Notes

For homogeneous materials, this returns (0, 0, 0). Gradient drives ray curvature via the eikonal equation: d²r/ds² = ∇n(r)/n(r)

References

abstractmethod get_absorption_coefficient(x, y, z, wavelength)[source]

Get absorption coefficient α at position (x, y, z).

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

alpha – Absorption coefficient in m⁻¹.

Return type:

float or ndarray

Notes

Intensity decays as I(d) = I₀ exp(-αd) (Beer-Lambert law).

References

abstractmethod get_scattering_coefficient(x, y, z, wavelength)[source]

Get scattering coefficient μ_s at position (x, y, z).

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

mu_s – Scattering coefficient in m⁻¹.

Return type:

float or ndarray

Notes

Mean free path between scattering events: ℓ_s = 1/μ_s.

get_extinction_coefficient(x, y, z, wavelength)[source]

Get total extinction coefficient (absorption + scattering).

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

mu_t – Total extinction coefficient in m⁻¹.

Return type:

float or ndarray

get_anisotropy_factor(x, y, z, wavelength)[source]

Get scattering anisotropy factor g (Henyey-Greenstein parameter).

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

g – Anisotropy factor, range [-1, 1]. - g = 0: isotropic scattering - g > 0: forward scattering - g < 0: backward scattering

Return type:

float or ndarray

Notes

Default implementation returns 0 (isotropic). Used in Henyey-Greenstein phase function: p(cos θ) = (1-g²) / (1 + g² - 2g cos θ)^(3/2)

References

is_homogeneous()[source]

Check if material has uniform properties (no spatial variation).

Returns:: True if material is homogeneous, False if inhomogeneous.
Return type:: bool

Notes

Homogeneous materials enable optimizations (straight-line propagation).

classmethod supported_kernels()[source]

Return list of propagation kernels supported by this material type.

Returns:: List of kernel IDs this material supports.
Return type:: list[PropagationKernelID]

Examples

>>> from lsurf.materials import ExponentialAtmosphere
>>> ExponentialAtmosphere.supported_kernels()
[<PropagationKernelID.SIMPLE_EULER: ...>, <PropagationKernelID.SIMPLE_RK4: ...>]

classmethod default_kernel()[source]

Return the default propagation kernel for this material type.

Returns:: The default kernel ID, or None if no GPU kernels are supported.
Return type:: PropagationKernelID or None

Examples

>>> from lsurf.materials import ExponentialAtmosphere
>>> ExponentialAtmosphere.default_kernel()
<PropagationKernelID.SIMPLE_RK4: ...>

classmethod supported_propagators()[source]

Return list of propagators supported by this material type.

Returns:: List of propagator IDs this material supports.
Return type:: list[PropagatorID]

Examples

>>> from lsurf.materials import ExponentialAtmosphere
>>> ExponentialAtmosphere.supported_propagators()
[<PropagatorID.GPU_GRADIENT: ...>, <PropagatorID.CPU_GRADIENT: ...>]

classmethod default_propagator()[source]

Return the default propagator for this material type.

Returns:: The default propagator ID, or None if not specified.
Return type:: PropagatorID or None

Examples

>>> from lsurf.materials import ExponentialAtmosphere
>>> ExponentialAtmosphere.default_propagator()
<PropagatorID.GPU_GRADIENT: ...>

property kernel_id: PropagationKernelID | None

Return the kernel ID configured for this material instance.

Returns:: The kernel ID selected for this instance, or None if not applicable.
Return type:: PropagationKernelID or None

property propagator_id: PropagatorID | None

Return the propagator ID configured for this material instance.

Returns:: The propagator ID selected for this instance, or None if not applicable.
Return type:: PropagatorID or None

get_refractive_index_gradient_magnitude(x, y, z, wavelength)[source]

Get magnitude of refractive index gradient |∇n|.

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

grad_mag – Magnitude of gradient |∇n| in m⁻¹.

Return type:

float or ndarray

compute_phase_velocity(x, y, z, wavelength)[source]

Compute phase velocity v_p = c/n at position.

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

v_p – Phase velocity in m/s.

Return type:

float or ndarray

__repr__()[source]

Return string representation.

class lsurf.materials.HomogeneousMaterial(name, refractive_index, absorption_coef=0.0, scattering_coef=0.0, anisotropy=0.0, propagator=None)[source]

Bases: MaterialField

Homogeneous material with constant optical properties.

Properties are uniform throughout space but can depend on wavelength. Supports analytic dispersion models (Sellmeier, Cauchy) or custom functions.

Parameters:

name (str) – Descriptive name for this material.
refractive_index (float or callable) – Refractive index value or function f(wavelength) -> n. If float, uses constant n for all wavelengths. If callable, wavelength is in meters.
absorption_coef (float, optional) – Absorption coefficient α in m⁻¹. Default is 0.
scattering_coef (float, optional) – Scattering coefficient μ_s in m⁻¹. Default is 0.
anisotropy (float, optional) – Henyey-Greenstein anisotropy factor g ∈ [-1, 1]. Default is 0.

refractive_index

Refractive index specification.

Type:: float or callable

absorption_coef

Absorption coefficient in m⁻¹.

Type:: float

scattering_coef

Scattering coefficient in m⁻¹.

Type:: float

anisotropy

Scattering anisotropy factor.

Type:: float

Examples

>>> # Constant refractive index
>>> glass = HomogeneousMaterial("Glass", 1.5)

>>> # Wavelength-dependent (Cauchy dispersion)
>>> def cauchy_n(wavelength):
...     wl_um = wavelength * 1e6
...     return 1.5 + 0.01 / wl_um**2
>>> glass = HomogeneousMaterial("Glass", cauchy_n)

>>> # With absorption (colored glass)
>>> colored_glass = HomogeneousMaterial(
...     "Red Glass", 1.52, absorption_coef=0.1
... )

Notes

For homogeneous materials: - Rays propagate in straight lines (no gradient-driven bending) - Reflection/refraction occurs only at interfaces - Beer-Lambert absorption: I(d) = I₀ exp(-αd) - No GPU kernels are needed (straight-line propagation)

__init__(name, refractive_index, absorption_coef=0.0, scattering_coef=0.0, anisotropy=0.0, propagator=None)[source]

Initialize homogeneous material.

Parameters:

name (str) – Descriptive name for this material.
refractive_index (float or callable) – Refractive index value or function f(wavelength) -> n.
absorption_coef (float, optional) – Absorption coefficient α in m⁻¹. Default is 0.
scattering_coef (float, optional) – Scattering coefficient μ_s in m⁻¹. Default is 0.
anisotropy (float, optional) – Henyey-Greenstein anisotropy factor g ∈ [-1, 1]. Default is 0.
propagator (PropagatorID, optional) – Override the default propagator. Only CPU_GRADIENT is supported.

Raises:

ValueError – If anisotropy is outside [-1, 1].

get_refractive_index(x, y, z, wavelength)[source]

Get refractive index (position-independent, wavelength-dependent).

Parameters:

x (float or ndarray) – Position coordinates in meters (ignored for homogeneous).
y (float or ndarray) – Position coordinates in meters (ignored for homogeneous).
z (float or ndarray) – Position coordinates in meters (ignored for homogeneous).
wavelength (float or ndarray) – Wavelength in meters.

Returns:

n – Refractive index.

Return type:

float or ndarray

get_refractive_index_gradient(x, y, z, wavelength)[source]

Get refractive index gradient (always zero for homogeneous).

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters (unused).

Returns:

grad_n – (0, 0, 0) since ∇n = 0 for homogeneous materials.

Return type:

tuple of (float or ndarray)

get_absorption_coefficient(x, y, z, wavelength)[source]

Get absorption coefficient (constant for homogeneous).

Parameters:

x (float or ndarray) – Position coordinates in meters (ignored).
y (float or ndarray) – Position coordinates in meters (ignored).
z (float or ndarray) – Position coordinates in meters (ignored).
wavelength (float or ndarray) – Wavelength in meters (currently ignored, could extend).

Returns:

alpha – Absorption coefficient in m⁻¹.

Return type:

float or ndarray

get_scattering_coefficient(x, y, z, wavelength)[source]

Get scattering coefficient (constant for homogeneous).

Parameters:

x (float or ndarray) – Position coordinates in meters (ignored).
y (float or ndarray) – Position coordinates in meters (ignored).
z (float or ndarray) – Position coordinates in meters (ignored).
wavelength (float or ndarray) – Wavelength in meters (currently ignored, could extend).

Returns:

mu_s – Scattering coefficient in m⁻¹.

Return type:

float or ndarray

get_anisotropy_factor(x, y, z, wavelength)[source]

Get scattering anisotropy factor (constant for homogeneous).

Parameters:

x (float or ndarray) – Position coordinates in meters (ignored).
y (float or ndarray) – Position coordinates in meters (ignored).
z (float or ndarray) – Position coordinates in meters (ignored).
wavelength (float or ndarray) – Wavelength in meters (ignored).

Returns:

g – Anisotropy factor.

Return type:

float or ndarray

__repr__()[source]

Return string representation with key properties.

class lsurf.materials.SimpleInhomogeneousModel(name='Simple Inhomogeneous', center=(0.0, 0.0, 0.0), reference_radius=6371000.0, altitude_range=(0.0, 200000.0), lut_resolution=10000, kernel=None, propagator=None)[source]

Bases: MaterialField, ABC

Base class for radially-symmetric inhomogeneous materials.

User implements only n_at_altitude() - a simple 1D function that returns the refractive index at a given altitude. The system automatically provides: - get_refractive_index(x, y, z, wavelength) via altitude computation - get_refractive_index_gradient(x, y, z, wavelength) via radial direction - GPU support via precomputed lookup tables

Parameters:

name (str) – Name of the material model
center (tuple) – Center of spherical symmetry (x, y, z). Default: (0, 0, 0)
reference_radius (float) – Radius at which altitude = 0 (e.g., Earth radius). Default: 6,371,000 m
altitude_range (tuple) – (min_altitude, max_altitude) for lookup table. Default: (0, 200,000)
lut_resolution (int) – Number of samples in lookup table. Default: 10,000
kernel (PropagationKernelID, optional) – Override the default kernel. Supported: SIMPLE_EULER, SIMPLE_RK4
propagator (PropagatorID, optional) – Override the default propagator. Supported: GPU_GRADIENT, CPU_GRADIENT

Example

>>> class MyAtmosphere(SimpleInhomogeneousModel):
...     def n_at_altitude(self, altitude, wavelength=None):
...         return 1.0 + 0.000293 * math.exp(-altitude / 8500)
...
>>> atm = MyAtmosphere(name="My Atmosphere")
>>> propagator = GPUGradientPropagator(atm)  # Auto GPU support

__init__(name='Simple Inhomogeneous', center=(0.0, 0.0, 0.0), reference_radius=6371000.0, altitude_range=(0.0, 200000.0), lut_resolution=10000, kernel=None, propagator=None)[source]

Initialize material field.

Parameters:

name (str, optional) – Descriptive name for this material. Default is “Material”.
kernel (PropagationKernelID, optional) – Override the default propagation kernel. Must be in supported_kernels(). If None, uses the class default.
propagator (PropagatorID, optional) – Override the default propagator. Must be in supported_propagators(). If None, uses the class default.

Raises:

ValueError – If kernel or propagator is not supported by this material type.

abstractmethod n_at_altitude(altitude, wavelength=None)[source]

Return refractive index at given altitude.

This is the only method users must implement.

Parameters:

altitude (float) – Altitude above reference surface in meters (altitude >= 0)
wavelength (float, optional) – Wavelength in meters (for dispersion models)

Returns:

Refractive index n >= 1.0

Return type:

float

dn_dh_at_altitude(altitude, wavelength=None)[source]

Return dn/d(altitude) at given altitude.

Default implementation uses numerical differentiation. Override for analytical derivatives (more accurate).

Parameters:

altitude (float) – Altitude above reference surface in meters
wavelength (float, optional) – Wavelength in meters

Returns:

Derivative dn/dh in m^-1

Return type:

float

alpha_at_altitude(altitude, wavelength=None)[source]

Return absorption coefficient at given altitude.

Default implementation returns 0 (no absorption). Override in subclasses that model absorbing media.

Parameters:

altitude (float) – Altitude above reference surface in meters
wavelength (float, optional) – Wavelength in meters

Returns:

Absorption coefficient α in m^-1 (for Beer-Lambert: I = I₀·exp(-α·s))

Return type:

float

get_refractive_index(x, y, z, wavelength)[source]

Get refractive index at position (auto-generated from n_at_altitude).

get_refractive_index_gradient(x, y, z, wavelength)[source]

Get gradient of n at position (auto-generated).

get_absorption_coefficient(x, y, z, wavelength)[source]

Return absorption coefficient (default: 0).

get_scattering_coefficient(x, y, z, wavelength)[source]

Return scattering coefficient (default: 0).

property gpu_material_id: int: Return GPU material ID for kernel dispatch.

get_gpu_kernels()[source]

Return GPU kernels for propagation.

get_gpu_parameters()[source]

Return scalar parameters for GPU kernel.

get_gpu_arrays()[source]

Return device arrays for GPU kernel.

class lsurf.materials.GridInhomogeneousModel(name='Grid Inhomogeneous', n_grid=None, bounds=None, gradient_grid=None, alpha_grid=None, kernel=None, propagator=None)[source]

Bases: MaterialField

Base class for 3D inhomogeneous materials defined on a grid.

User provides a 3D array of refractive index values. The system provides: - Trilinear interpolation for arbitrary positions - Automatic gradient computation via finite differences - GPU support via 3D array interpolation

Parameters:

name (str) – Name of the material model
n_grid (ndarray (nx, ny, nz)) – Refractive index values on 3D grid
bounds (tuple) – ((x_min, x_max), (y_min, y_max), (z_min, z_max))
gradient_grid (ndarray (nx, ny, nz, 3), optional) – Pre-computed gradients. If None, computed via finite differences.
alpha_grid (ndarray (nx, ny, nz), optional) – Absorption coefficient values on 3D grid. If None, defaults to zeros.

Example

>>> n_grid = np.ones((100, 100, 50)) * 1.000293
>>> n_grid += generate_turbulence(100, 100, 50)  # Add perturbations
>>> material = GridInhomogeneousModel(
...     name="Turbulent Field",
...     n_grid=n_grid,
...     bounds=((-5000, 5000), (-5000, 5000), (0, 25000))
... )

__init__(name='Grid Inhomogeneous', n_grid=None, bounds=None, gradient_grid=None, alpha_grid=None, kernel=None, propagator=None)[source]

Initialize material field.

Parameters:

name (str, optional) – Descriptive name for this material. Default is “Material”.
kernel (PropagationKernelID, optional) – Override the default propagation kernel. Must be in supported_kernels(). If None, uses the class default.
propagator (PropagatorID, optional) – Override the default propagator. Must be in supported_propagators(). If None, uses the class default.

Raises:

ValueError – If kernel or propagator is not supported by this material type.

get_refractive_index(x, y, z, wavelength)[source]

Get refractive index via trilinear interpolation.

get_refractive_index_gradient(x, y, z, wavelength)[source]

Get gradient via trilinear interpolation in gradient grid.

get_absorption_coefficient(x, y, z, wavelength)[source]

Return absorption coefficient (default: 0).

get_scattering_coefficient(x, y, z, wavelength)[source]

Return scattering coefficient (default: 0).

property gpu_material_id: int: Return GPU material ID for kernel dispatch.

get_gpu_kernels()[source]

Return GPU kernels for propagation.

get_gpu_parameters()[source]

Return scalar parameters for GPU kernel.

get_gpu_arrays()[source]

Return device arrays for GPU kernel.

class lsurf.materials.FullInhomogeneousModel(name='Full Inhomogeneous')[source]

Bases: MaterialField, ABC

Base class for arbitrary 3D inhomogeneous materials (CPU-only).

Use this when you need maximum flexibility and cannot express your model using the other tiers. Any Python code is allowed, including external libraries, database lookups, or complex computations.

GPU acceleration is not available for this tier. Use GradientPropagator for CPU propagation.

User must implement: - get_refractive_index(x, y, z, wavelength)

User may optionally implement: - get_refractive_index_gradient(x, y, z, wavelength) (default: numerical)

Example

>>> class WeatherModel(FullInhomogeneousModel):
...     def __init__(self, weather_file):
...         super().__init__(name="Weather Model")
...         self._data = load_weather_data(weather_file)
...
...     def get_refractive_index(self, x, y, z, wavelength):
...         return self._data.interpolate(x, y, z)
...
>>> model = WeatherModel("weather.nc")
>>> propagator = GradientPropagator()  # CPU only

__init__(name='Full Inhomogeneous')[source]

Initialize material field.

Parameters:

name (str, optional) – Descriptive name for this material. Default is “Material”.
kernel (PropagationKernelID, optional) – Override the default propagation kernel. Must be in supported_kernels(). If None, uses the class default.
propagator (PropagatorID, optional) – Override the default propagator. Must be in supported_propagators(). If None, uses the class default.

Raises:

ValueError – If kernel or propagator is not supported by this material type.

abstractmethod get_refractive_index(x, y, z, wavelength)[source]

Get refractive index at position.

Any Python code is allowed here.

Parameters:

x (float or ndarray) – Position coordinates in meters
y (float or ndarray) – Position coordinates in meters
z (float or ndarray) – Position coordinates in meters
wavelength (float or ndarray) – Wavelength in meters

Returns:

Refractive index n >= 1.0

Return type:

float or ndarray

get_refractive_index_gradient(x, y, z, wavelength)[source]

Get gradient of n at position.

Default implementation uses numerical differentiation. Override for analytical gradients.

Parameters:

x (float or ndarray) – Position coordinates in meters
y (float or ndarray) – Position coordinates in meters
z (float or ndarray) – Position coordinates in meters
wavelength (float or ndarray) – Wavelength in meters

Returns:

(dn/dx, dn/dy, dn/dz) in m^-1

Return type:

tuple

get_absorption_coefficient(x, y, z, wavelength)[source]

Return absorption coefficient (default: 0).

get_scattering_coefficient(x, y, z, wavelength)[source]

Return scattering coefficient (default: 0).

property gpu_material_id: int: GPU is not supported for FullInhomogeneousModel.

get_gpu_kernels()[source]

GPU is not supported for FullInhomogeneousModel.

get_gpu_parameters()[source]

GPU is not supported for FullInhomogeneousModel.

class lsurf.materials.ExponentialAtmosphere(name='Exponential Atmosphere', n_sea_level=1.000293, scale_height=8500.0, earth_radius=6371000.0, earth_center=(0.0, 0.0, 0.0), absorption_coef=0.0, absorption_scale_height=None, kernel=None, propagator=None)[source]

Bases: SimpleInhomogeneousModel

Exponential atmosphere with radially-dependent refractive index.

Models an atmosphere where air density decreases exponentially with altitude, causing the refractive index to approach 1 at high altitudes.

The model is spherically symmetric about Earth’s center, which is assumed to be at the origin (0, 0, 0) by default, or can be specified.

This class inherits from SimpleInhomogeneousModel, implementing the n_at_altitude() method for the exponential profile. GPU support is provided automatically via lookup table interpolation.

Parameters:

name (str, optional) – Descriptive name for this material. Default is “Exponential Atmosphere”.
n_sea_level (float, optional) – Refractive index at sea level (Earth’s surface). Default is 1.000293.
scale_height (float, optional) – Atmospheric scale height H in meters. Default is 8500.0 m.
earth_radius (float, optional) – Radius of Earth in meters. Default is 6,371,000 m.
earth_center (tuple of float, optional) – Position of Earth’s center in meters. Default is (0, 0, 0).
absorption_coef (float, optional) – Absorption coefficient at sea level in m⁻¹. Default is 0.0.
absorption_scale_height (float, optional) – Scale height for absorption (can differ from density). Default is same as scale_height.

Examples

>>> atmosphere = ExponentialAtmosphere()
>>> # Get refractive index at 10 km altitude
>>> n = atmosphere.get_refractive_index(0, 0, EARTH_RADIUS + 10000, 532e-9)
>>> print(f"n at 10 km: {n:.6f}")  # ~1.000089

>>> # Get gradient at sea level directly above Earth's center
>>> grad = atmosphere.get_refractive_index_gradient(0, 0, EARTH_RADIUS, 532e-9)
>>> print(f"dn/dz at sea level: {grad[2]:.2e}")  # ~ -3.4e-8 m^-1

__init__(name='Exponential Atmosphere', n_sea_level=1.000293, scale_height=8500.0, earth_radius=6371000.0, earth_center=(0.0, 0.0, 0.0), absorption_coef=0.0, absorption_scale_height=None, kernel=None, propagator=None)[source]

Initialize material field.

Parameters:

name (str, optional) – Descriptive name for this material. Default is “Material”.
kernel (PropagationKernelID, optional) – Override the default propagation kernel. Must be in supported_kernels(). If None, uses the class default.
propagator (PropagatorID, optional) – Override the default propagator. Must be in supported_propagators(). If None, uses the class default.

Raises:

ValueError – If kernel or propagator is not supported by this material type.

n_at_altitude(altitude, wavelength=None)[source]

Return refractive index at given altitude.

Implements the exponential profile:: n(h) = 1 + delta_n * exp(-h / H)

Parameters:

altitude (float) – Altitude above Earth’s surface in meters (clamped to >= 0)
wavelength (float, optional) – Wavelength in meters (not used - no dispersion)

Returns:

Refractive index at altitude

Return type:

float

dn_dh_at_altitude(altitude, wavelength=None)[source]

Return analytical dn/dh at given altitude.

Derivative of exponential profile:: dn/dh = -(delta_n / H) * exp(-h / H)

Parameters:

altitude (float) – Altitude above Earth’s surface in meters
wavelength (float, optional) – Wavelength in meters (not used)

Returns:

Derivative dn/dh in m^-1

Return type:

float

alpha_at_altitude(altitude, wavelength=None)[source]

Return absorption coefficient at given altitude.

Implements exponential absorption profile:: α(h) = α₀ * exp(-h / H_α)

Parameters:

altitude (float) – Altitude above Earth’s surface in meters
wavelength (float, optional) – Wavelength in meters (not used)

Returns:

Absorption coefficient α in m⁻¹

Return type:

float

get_absorption_coefficient(x, y, z, wavelength)[source]

Get absorption coefficient at position (x, y, z).

Absorption follows an exponential profile with altitude.

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

alpha – Absorption coefficient in m⁻¹.

Return type:

float or ndarray

get_density_ratio(x, y, z)[source]

Get atmospheric density ratio relative to sea level.

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.

Returns:

rho_ratio – Density relative to sea level (dimensionless). rho_ratio = rho(r) / rho_0 = exp(-altitude / H)

Return type:

float or ndarray

compute_curvature_radius(x, y, z, wavelength=5.32e-07)[source]

Compute the radius of curvature for a ray at given position.

The radius of curvature is R_c = n / |∇n|.

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float, optional) – Wavelength in meters. Default is 532 nm.

Returns:

R_c – Radius of curvature in meters. Very large values indicate nearly straight propagation.

Return type:

float or ndarray

__repr__()[source]

Return string representation.

lsurf.materials.create_exponential_atmosphere(n_sea_level=1.000293, scale_height=8500.0, earth_radius=6371000.0, earth_center=(0.0, 0.0, 0.0), name='Custom Exponential Atmosphere')[source]

Factory function to create an exponential atmosphere.

Parameters:

n_sea_level (float, optional) – Refractive index at sea level. Default is 1.000293.
scale_height (float, optional) – Atmospheric scale height in meters. Default is 8500.0 m.
earth_radius (float, optional) – Radius of Earth in meters. Default is 6,371,000 m.
earth_center (tuple of float, optional) – Position of Earth’s center. Default is (0, 0, 0).
name (str, optional) – Descriptive name for the material.

Returns:

Configured atmosphere material.

Return type:

ExponentialAtmosphere

Examples

>>> # Mars-like atmosphere (thinner, smaller planet)
>>> mars_atmo = create_exponential_atmosphere(
...     n_sea_level=1.0001,
...     scale_height=11_100.0,  # Mars scale height
...     earth_radius=3_389_500.0,  # Mars radius
... )

class lsurf.materials.DuctAtmosphere(name='Duct Atmosphere', n_sea_level=1.000293, scale_height=8500.0, earth_radius=6371000.0, earth_center=(0.0, 0.0, 0.0), duct_center=0.0, duct_width=100.0, duct_intensity=0.0, duct_sharpness=0.5)[source]

Bases: SimpleInhomogeneousModel

Exponential atmosphere with a refractive index duct layer.

Models an atmosphere where air density decreases exponentially with altitude, modified by a duct (inversion layer) at a specified altitude. The duct creates a localized region where the refractive index gradient is modified, which can trap or guide electromagnetic waves.

The duct is modeled using hyperbolic tangent functions to create smooth transitions at the duct boundaries.

Parameters:

name (str, optional) – Descriptive name for this material. Default is “Duct Atmosphere”.
n_sea_level (float, optional) – Refractive index at sea level. Default is 1.000293.
scale_height (float, optional) – Atmospheric scale height H in meters. Default is 8500.0 m.
earth_radius (float, optional) – Radius of Earth in meters. Default is 6,371,000 m.
earth_center (tuple of float, optional) – Position of Earth’s center in meters. Default is (0, 0, 0).
duct_center (float, optional) – Center altitude of the duct layer in meters. Default is 0.0.
duct_width (float, optional) – Width of the duct layer in meters. Default is 100.0.
duct_intensity (float, optional) – Intensity of the duct effect (0 = no duct, 1 = full strength). Default is 0.0.
duct_sharpness (float, optional) – Sharpness of duct boundaries (0 = gradual, 1 = sharp). Default is 0.5.

Examples

>>> # Surface duct at 100m altitude
>>> atmosphere = DuctAtmosphere(
...     duct_center=100.0,
...     duct_width=50.0,
...     duct_intensity=0.5,
...     duct_sharpness=0.8,
... )

>>> # Elevated duct for radio propagation study
>>> atmosphere = DuctAtmosphere(
...     duct_center=1000.0,
...     duct_width=200.0,
...     duct_intensity=0.8,
...     duct_sharpness=0.9,
...     earth_center=(0.0, 0.0, -EARTH_RADIUS),
... )

__init__(name='Duct Atmosphere', n_sea_level=1.000293, scale_height=8500.0, earth_radius=6371000.0, earth_center=(0.0, 0.0, 0.0), duct_center=0.0, duct_width=100.0, duct_intensity=0.0, duct_sharpness=0.5)[source]

Initialize material field.

Parameters:

name (str, optional) – Descriptive name for this material. Default is “Material”.
kernel (PropagationKernelID, optional) – Override the default propagation kernel. Must be in supported_kernels(). If None, uses the class default.
propagator (PropagatorID, optional) – Override the default propagator. Must be in supported_propagators(). If None, uses the class default.

Raises:

ValueError – If kernel or propagator is not supported by this material type.

n_at_altitude(altitude, wavelength=None)[source]

Return refractive index at given altitude.

Implements the exponential profile with duct modification:: n(h) = 1 + delta_n * exp(-h / H) * D(h)

Parameters:

altitude (float) – Altitude above Earth’s surface in meters (clamped to >= 0).
wavelength (float, optional) – Wavelength in meters (not used - no dispersion).

Returns:

Refractive index at altitude.

Return type:

float

dn_dh_at_altitude(altitude, wavelength=None)[source]

Return analytical dn/dh at given altitude.

Uses the product rule for the derivative:: dn/dh = delta_n * [exp’ * D + exp * D’]

Parameters:

altitude (float) – Altitude above Earth’s surface in meters.
wavelength (float, optional) – Wavelength in meters (not used).

Returns:

Derivative dn/dh in m^-1.

Return type:

float

property gpu_material_id: int: Return GPU material ID for analytical duct kernels.

get_gpu_kernels()[source]

Return GPU kernels for propagation (analytical, no LUT).

get_gpu_parameters()[source]

Return scalar parameters for GPU kernel.

get_duct_gpu_parameters()[source]

Return scalar parameters for analytical duct GPU kernels.

Equivalent to get_gpu_parameters() for DuctAtmosphere, but provides a consistent interface for propagators to detect analytical duct kernels.

get_gpu_arrays()[source]

Return device arrays for GPU kernel (none needed for analytical).

__repr__()[source]

Return string representation.

lsurf.materials.create_duct_atmosphere(n_sea_level=1.000293, scale_height=8500.0, earth_radius=6371000.0, earth_center=(0.0, 0.0, 0.0), duct_center=0.0, duct_width=100.0, duct_intensity=0.0, duct_sharpness=0.5, name='Custom Duct Atmosphere')[source]

Factory function to create a duct atmosphere.

Parameters:

n_sea_level (float, optional) – Refractive index at sea level. Default is 1.000293.
scale_height (float, optional) – Atmospheric scale height in meters. Default is 8500.0 m.
earth_radius (float, optional) – Radius of Earth in meters. Default is 6,371,000 m.
earth_center (tuple of float, optional) – Position of Earth’s center. Default is (0, 0, 0).
duct_center (float, optional) – Center altitude of the duct in meters. Default is 0.0.
duct_width (float, optional) – Width of the duct in meters. Default is 100.0.
duct_intensity (float, optional) – Intensity of the duct (0-1). Default is 0.0.
duct_sharpness (float, optional) – Sharpness of duct edges (0-1). Default is 0.5.
name (str, optional) – Descriptive name for the material.

Returns:

Configured atmosphere material with duct.

Return type:

DuctAtmosphere

Examples

>>> # Create a surface evaporation duct
>>> atmo = create_duct_atmosphere(
...     duct_center=50.0,
...     duct_width=30.0,
...     duct_intensity=0.7,
...     duct_sharpness=0.9,
... )

lsurf.materials.sellmeier_refractive_index(wavelength, B1=1.03961212, B2=0.231792344, B3=1.01046945, C1=0.00600069867, C2=0.0200179144, C3=103.560653)[source]

Compute refractive index using Sellmeier equation.

The Sellmeier equation is an empirical formula that accurately models the dispersion of optical glasses.

Parameters:

wavelength (float) – Wavelength in meters.
B1 (float) – Sellmeier B coefficients (dimensionless). Defaults are for N-BK7 glass.
B2 (float) – Sellmeier B coefficients (dimensionless). Defaults are for N-BK7 glass.
B3 (float) – Sellmeier B coefficients (dimensionless). Defaults are for N-BK7 glass.
C1 (float) – Sellmeier C coefficients in μm². Defaults are for N-BK7 glass.
C2 (float) – Sellmeier C coefficients in μm². Defaults are for N-BK7 glass.
C3 (float) – Sellmeier C coefficients in μm². Defaults are for N-BK7 glass.

Returns:

n – Refractive index at given wavelength.

Return type:

float

Notes

Sellmeier equation:

n² - 1 = B₁λ²/(λ² - C₁) + B₂λ²/(λ² - C₂) + B₃λ²/(λ² - C₃)

where λ is in micrometers.

Default coefficients are for Schott N-BK7 glass, valid 300-2500 nm.

References

Examples

>>> n_589nm = sellmeier_refractive_index(589e-9)
>>> print(f"n = {n_589nm:.6f}")  # ~1.5168

lsurf.materials.cauchy_refractive_index(wavelength, A=1.458, B=0.00354, C=0.0)[source]

Compute refractive index using Cauchy equation.

The Cauchy equation is a simple empirical dispersion formula, valid in regions of normal dispersion (away from absorption bands).

Parameters:

wavelength (float) – Wavelength in meters.
A (float) – Constant term (refractive index at infinite wavelength). Default is for typical crown glass.
B (float) – First-order dispersion coefficient in μm². Default is for typical crown glass.
C (float) – Second-order dispersion coefficient in μm⁴. Default is 0.

Returns:

n – Refractive index at given wavelength.

Return type:

float

Notes

Cauchy equation:

n(λ) = A + B/λ² + C/λ⁴

where λ is in micrometers.

Suitable for glasses and transparent materials in visible range. Less accurate than Sellmeier near absorption edges.

References

Examples

>>> n_550nm = cauchy_refractive_index(550e-9, A=1.5, B=0.004)
>>> print(f"n = {n_550nm:.4f}")

lsurf.materials.create_sellmeier_material(name, B1, B2, B3, C1, C2, C3, absorption_coef=0.0)[source]

Create a material with Sellmeier dispersion.

Parameters:

name (str) – Material name.
B1 (float) – Sellmeier B coefficients (dimensionless).
B2 (float) – Sellmeier B coefficients (dimensionless).
B3 (float) – Sellmeier B coefficients (dimensionless).
C1 (float) – Sellmeier C coefficients in μm².
C2 (float) – Sellmeier C coefficients in μm².
C3 (float) – Sellmeier C coefficients in μm².
absorption_coef (float, optional) – Absorption coefficient in m⁻¹. Default is 0.

Returns:

material – Material with wavelength-dependent refractive index.

Return type:

HomogeneousMaterial

Examples

>>> # Create N-BK7 glass with Sellmeier dispersion
>>> bk7 = create_sellmeier_material(
...     "N-BK7",
...     B1=1.03961212, B2=0.231792344, B3=1.01046945,
...     C1=6.00069867e-3, C2=2.00179144e-2, C3=1.03560653e2,
... )

lsurf.materials.create_cauchy_material(name, A, B, C=0.0, absorption_coef=0.0)[source]

Create a material with Cauchy dispersion.

Parameters:

name (str) – Material name.
A (float) – Constant term.
B (float) – First-order dispersion coefficient in μm².
C (float, optional) – Second-order dispersion coefficient in μm⁴. Default is 0.
absorption_coef (float, optional) – Absorption coefficient in m⁻¹. Default is 0.

Returns:

material – Material with wavelength-dependent refractive index.

Return type:

HomogeneousMaterial

Examples

>>> # Create simple glass with Cauchy dispersion
>>> glass = create_cauchy_material("Glass", A=1.52, B=0.004)

lsurf.materials.device_sellmeier_equation(wl_um, B1, B2, B3, C1, C2, C3)[source]

GPU-compatible Sellmeier equation.

For use in Numba CUDA kernels. Computes refractive index from wavelength and Sellmeier coefficients.

Parameters:

wl_um (float) – Wavelength in micrometers.
B1 (float) – Sellmeier B coefficients.
B2 (float) – Sellmeier B coefficients.
B3 (float) – Sellmeier B coefficients.
C1 (float) – Sellmeier C coefficients in μm².
C2 (float) – Sellmeier C coefficients in μm².
C3 (float) – Sellmeier C coefficients in μm².

Returns:

n – Refractive index.

Return type:

float

Notes

This function is designed to be compiled with @cuda.jit(device=True). Import and compile in your GPU module.

lsurf.materials.device_cauchy_equation(wl_um, A, B, C)[source]

GPU-compatible Cauchy equation.

For use in Numba CUDA kernels. Computes refractive index from wavelength and Cauchy coefficients.

Parameters:

wl_um (float) – Wavelength in micrometers.
A (float) – Constant term.
B (float) – First-order dispersion coefficient in μm².
C (float) – Second-order dispersion coefficient in μm⁴.

Returns:

n – Refractive index.

Return type:

float

Notes

This function is designed to be compiled with @cuda.jit(device=True). Import and compile in your GPU module.

Material Classes

`MaterialField`([name, kernel, propagator])	Abstract base class for spatially-varying material properties.
`HomogeneousMaterial`(name, refractive_index)	Homogeneous material with constant optical properties.

Inhomogeneous Material Base Classes

`SimpleInhomogeneousModel`([name, center, ...])	Base class for radially-symmetric inhomogeneous materials.
`GridInhomogeneousModel`([name, n_grid, ...])	Base class for 3D inhomogeneous materials defined on a grid.
`FullInhomogeneousModel`([name])	Base class for arbitrary 3D inhomogeneous materials (CPU-only).

Predefined Materials

`VACUUM`	Homogeneous material with constant optical properties.
`AIR_STP`	Homogeneous material with constant optical properties.
`WATER`	Homogeneous material with constant optical properties.
`BK7_GLASS`	Homogeneous material with constant optical properties.

Atmosphere Implementations

`ExponentialAtmosphere`([name, n_sea_level, ...])	Exponential atmosphere with radially-dependent refractive index.
`DuctAtmosphere`([name, n_sea_level, ...])	Exponential atmosphere with a refractive index duct layer.
`STANDARD_ATMOSPHERE`	Exponential atmosphere with radially-dependent refractive index.

Dispersion Functions

`sellmeier_refractive_index`(wavelength[, B1, ...])	Compute refractive index using Sellmeier equation.
`cauchy_refractive_index`(wavelength[, A, B, C])	Compute refractive index using Cauchy equation.
`create_sellmeier_material`(name, B1, B2, B3, ...)	Create a material with Sellmeier dispersion.
`create_cauchy_material`(name, A, B[, C, ...])	Create a material with Cauchy dispersion.

Material Classes

MaterialField

class lsurf.materials.MaterialField(name='Material', kernel=None, propagator=None)[source]

Bases: ABC

Abstract base class for spatially-varying material properties.

All material properties are functions of position (x, y, z) and wavelength. Subclasses must implement methods that can be called from both CPU and GPU.

Parameters:: name (str, optional) – Descriptive name for this material. Default is “Material”.

name

Material identifier.

Type:: str

Notes

References

__init__(name='Material', kernel=None, propagator=None)[source]

Initialize material field.

Parameters:

name (str, optional) – Descriptive name for this material. Default is “Material”.
kernel (PropagationKernelID, optional) – Override the default propagation kernel. Must be in supported_kernels(). If None, uses the class default.
propagator (PropagatorID, optional) – Override the default propagator. Must be in supported_propagators(). If None, uses the class default.

Raises:

ValueError – If kernel or propagator is not supported by this material type.

abstractmethod get_refractive_index(x, y, z, wavelength)[source]

Get refractive index at position (x, y, z) for given wavelength.

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

n – Real part of refractive index (dimensionless).

Return type:

float or ndarray

Notes

For absorbing materials, this returns only Re(n). Use get_absorption_coefficient() for the imaginary part.

abstractmethod get_refractive_index_gradient(x, y, z, wavelength)[source]

Get gradient of refractive index ∇n at position (x, y, z).

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

grad_n – (∂n/∂x, ∂n/∂y, ∂n/∂z) in units of m⁻¹.

Return type:

tuple of (float or ndarray)

Notes

For homogeneous materials, this returns (0, 0, 0). Gradient drives ray curvature via the eikonal equation: d²r/ds² = ∇n(r)/n(r)

References

abstractmethod get_absorption_coefficient(x, y, z, wavelength)[source]

Get absorption coefficient α at position (x, y, z).

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

alpha – Absorption coefficient in m⁻¹.

Return type:

float or ndarray

Notes

Intensity decays as I(d) = I₀ exp(-αd) (Beer-Lambert law).

References

abstractmethod get_scattering_coefficient(x, y, z, wavelength)[source]

Get scattering coefficient μ_s at position (x, y, z).

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

mu_s – Scattering coefficient in m⁻¹.

Return type:

float or ndarray

Notes

Mean free path between scattering events: ℓ_s = 1/μ_s.

get_extinction_coefficient(x, y, z, wavelength)[source]

Get total extinction coefficient (absorption + scattering).

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

mu_t – Total extinction coefficient in m⁻¹.

Return type:

float or ndarray

get_anisotropy_factor(x, y, z, wavelength)[source]

Get scattering anisotropy factor g (Henyey-Greenstein parameter).

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

g – Anisotropy factor, range [-1, 1]. - g = 0: isotropic scattering - g > 0: forward scattering - g < 0: backward scattering

Return type:

float or ndarray

Notes

Default implementation returns 0 (isotropic). Used in Henyey-Greenstein phase function: p(cos θ) = (1-g²) / (1 + g² - 2g cos θ)^(3/2)

References

is_homogeneous()[source]

Check if material has uniform properties (no spatial variation).

Returns:: True if material is homogeneous, False if inhomogeneous.
Return type:: bool

Notes

Homogeneous materials enable optimizations (straight-line propagation).

classmethod supported_kernels()[source]

Return list of propagation kernels supported by this material type.

Returns:: List of kernel IDs this material supports.
Return type:: list[PropagationKernelID]

Examples

>>> from lsurf.materials import ExponentialAtmosphere
>>> ExponentialAtmosphere.supported_kernels()
[<PropagationKernelID.SIMPLE_EULER: ...>, <PropagationKernelID.SIMPLE_RK4: ...>]

classmethod default_kernel()[source]

Return the default propagation kernel for this material type.

Returns:: The default kernel ID, or None if no GPU kernels are supported.
Return type:: PropagationKernelID or None

Examples

>>> from lsurf.materials import ExponentialAtmosphere
>>> ExponentialAtmosphere.default_kernel()
<PropagationKernelID.SIMPLE_RK4: ...>

classmethod supported_propagators()[source]

Return list of propagators supported by this material type.

Returns:: List of propagator IDs this material supports.
Return type:: list[PropagatorID]

Examples

>>> from lsurf.materials import ExponentialAtmosphere
>>> ExponentialAtmosphere.supported_propagators()
[<PropagatorID.GPU_GRADIENT: ...>, <PropagatorID.CPU_GRADIENT: ...>]

classmethod default_propagator()[source]

Return the default propagator for this material type.

Returns:: The default propagator ID, or None if not specified.
Return type:: PropagatorID or None

Examples

>>> from lsurf.materials import ExponentialAtmosphere
>>> ExponentialAtmosphere.default_propagator()
<PropagatorID.GPU_GRADIENT: ...>

property kernel_id: PropagationKernelID | None

Return the kernel ID configured for this material instance.

Returns:: The kernel ID selected for this instance, or None if not applicable.
Return type:: PropagationKernelID or None

property propagator_id: PropagatorID | None

Return the propagator ID configured for this material instance.

Returns:: The propagator ID selected for this instance, or None if not applicable.
Return type:: PropagatorID or None

get_refractive_index_gradient_magnitude(x, y, z, wavelength)[source]

Get magnitude of refractive index gradient |∇n|.

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

grad_mag – Magnitude of gradient |∇n| in m⁻¹.

Return type:

float or ndarray

compute_phase_velocity(x, y, z, wavelength)[source]

Compute phase velocity v_p = c/n at position.

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters.

Returns:

v_p – Phase velocity in m/s.

Return type:

float or ndarray

__repr__()[source]

Return string representation.

HomogeneousMaterial

class lsurf.materials.HomogeneousMaterial(name, refractive_index, absorption_coef=0.0, scattering_coef=0.0, anisotropy=0.0, propagator=None)[source]

Bases: MaterialField

Homogeneous material with constant optical properties.

Properties are uniform throughout space but can depend on wavelength. Supports analytic dispersion models (Sellmeier, Cauchy) or custom functions.

Parameters:

name (str) – Descriptive name for this material.
refractive_index (float or callable) – Refractive index value or function f(wavelength) -> n. If float, uses constant n for all wavelengths. If callable, wavelength is in meters.
absorption_coef (float, optional) – Absorption coefficient α in m⁻¹. Default is 0.
scattering_coef (float, optional) – Scattering coefficient μ_s in m⁻¹. Default is 0.
anisotropy (float, optional) – Henyey-Greenstein anisotropy factor g ∈ [-1, 1]. Default is 0.

refractive_index

Refractive index specification.

Type:: float or callable

absorption_coef

Absorption coefficient in m⁻¹.

Type:: float

scattering_coef

Scattering coefficient in m⁻¹.

Type:: float

anisotropy

Scattering anisotropy factor.

Type:: float

Examples

>>> # Constant refractive index
>>> glass = HomogeneousMaterial("Glass", 1.5)

>>> # Wavelength-dependent (Cauchy dispersion)
>>> def cauchy_n(wavelength):
...     wl_um = wavelength * 1e6
...     return 1.5 + 0.01 / wl_um**2
>>> glass = HomogeneousMaterial("Glass", cauchy_n)

>>> # With absorption (colored glass)
>>> colored_glass = HomogeneousMaterial(
...     "Red Glass", 1.52, absorption_coef=0.1
... )

Notes

__init__(name, refractive_index, absorption_coef=0.0, scattering_coef=0.0, anisotropy=0.0, propagator=None)[source]

Initialize homogeneous material.

Parameters:

name (str) – Descriptive name for this material.
refractive_index (float or callable) – Refractive index value or function f(wavelength) -> n.
absorption_coef (float, optional) – Absorption coefficient α in m⁻¹. Default is 0.
scattering_coef (float, optional) – Scattering coefficient μ_s in m⁻¹. Default is 0.
anisotropy (float, optional) – Henyey-Greenstein anisotropy factor g ∈ [-1, 1]. Default is 0.
propagator (PropagatorID, optional) – Override the default propagator. Only CPU_GRADIENT is supported.

Raises:

ValueError – If anisotropy is outside [-1, 1].

get_refractive_index(x, y, z, wavelength)[source]

Get refractive index (position-independent, wavelength-dependent).

Parameters:

x (float or ndarray) – Position coordinates in meters (ignored for homogeneous).
y (float or ndarray) – Position coordinates in meters (ignored for homogeneous).
z (float or ndarray) – Position coordinates in meters (ignored for homogeneous).
wavelength (float or ndarray) – Wavelength in meters.

Returns:

n – Refractive index.

Return type:

float or ndarray

get_refractive_index_gradient(x, y, z, wavelength)[source]

Get refractive index gradient (always zero for homogeneous).

Parameters:

x (float or ndarray) – Position coordinates in meters.
y (float or ndarray) – Position coordinates in meters.
z (float or ndarray) – Position coordinates in meters.
wavelength (float or ndarray) – Wavelength in meters (unused).

Returns:

grad_n – (0, 0, 0) since ∇n = 0 for homogeneous materials.

Return type:

tuple of (float or ndarray)

get_absorption_coefficient(x, y, z, wavelength)[source]

Get absorption coefficient (constant for homogeneous).

Parameters:

x (float or ndarray) – Position coordinates in meters (ignored).
y (float or ndarray) – Position coordinates in meters (ignored).
z (float or ndarray) – Position coordinates in meters (ignored).
wavelength (float or ndarray) – Wavelength in meters (currently ignored, could extend).

Returns:

alpha – Absorption coefficient in m⁻¹.

Return type:

float or ndarray

get_scattering_coefficient(x, y, z, wavelength)[source]

Get scattering coefficient (constant for homogeneous).

Parameters:

x (float or ndarray) – Position coordinates in meters (ignored).
y (float or ndarray) – Position coordinates in meters (ignored).
z (float or ndarray) – Position coordinates in meters (ignored).
wavelength (float or ndarray) – Wavelength in meters (currently ignored, could extend).

Returns:

mu_s – Scattering coefficient in m⁻¹.

Return type:

float or ndarray

get_anisotropy_factor(x, y, z, wavelength)[source]

Get scattering anisotropy factor (constant for homogeneous).

Parameters:

x (float or ndarray) – Position coordinates in meters (ignored).
y (float or ndarray) – Position coordinates in meters (ignored).
z (float or ndarray) – Position coordinates in meters (ignored).
wavelength (float or ndarray) – Wavelength in meters (ignored).

Returns:

g – Anisotropy factor.

Return type:

float or ndarray

__repr__()[source]

Return string representation with key properties.

Pre-defined Materials

lsurf.materials.VACUUM: Perfect vacuum with n = 1.0 (no dispersion).

lsurf.materials.AIR_STP: Air at standard temperature and pressure. Refractive index ~1.000293 at 589 nm.

lsurf.materials.WATER: Pure water with wavelength-dependent refractive index using empirical formula. n ≈ 1.333 at 589 nm.

lsurf.materials.BK7_GLASS: BK7 optical glass with Sellmeier dispersion formula. n ≈ 1.517 at 589 nm.

Dispersion Models

sellmeier_refractive_index

lsurf.materials.sellmeier_refractive_index(wavelength, B1=1.03961212, B2=0.231792344, B3=1.01046945, C1=0.00600069867, C2=0.0200179144, C3=103.560653)[source]

Compute refractive index using Sellmeier equation.

The Sellmeier equation is an empirical formula that accurately models the dispersion of optical glasses.

Parameters:

wavelength (float) – Wavelength in meters.
B1 (float) – Sellmeier B coefficients (dimensionless). Defaults are for N-BK7 glass.
B2 (float) – Sellmeier B coefficients (dimensionless). Defaults are for N-BK7 glass.
B3 (float) – Sellmeier B coefficients (dimensionless). Defaults are for N-BK7 glass.
C1 (float) – Sellmeier C coefficients in μm². Defaults are for N-BK7 glass.
C2 (float) – Sellmeier C coefficients in μm². Defaults are for N-BK7 glass.
C3 (float) – Sellmeier C coefficients in μm². Defaults are for N-BK7 glass.

Returns:

n – Refractive index at given wavelength.

Return type:

float

Notes

Sellmeier equation:

n² - 1 = B₁λ²/(λ² - C₁) + B₂λ²/(λ² - C₂) + B₃λ²/(λ² - C₃)

where λ is in micrometers.

Default coefficients are for Schott N-BK7 glass, valid 300-2500 nm.

References

Examples

>>> n_589nm = sellmeier_refractive_index(589e-9)
>>> print(f"n = {n_589nm:.6f}")  # ~1.5168

cauchy_refractive_index

lsurf.materials.cauchy_refractive_index(wavelength, A=1.458, B=0.00354, C=0.0)[source]

Compute refractive index using Cauchy equation.

The Cauchy equation is a simple empirical dispersion formula, valid in regions of normal dispersion (away from absorption bands).

Parameters:

wavelength (float) – Wavelength in meters.
A (float) – Constant term (refractive index at infinite wavelength). Default is for typical crown glass.
B (float) – First-order dispersion coefficient in μm². Default is for typical crown glass.
C (float) – Second-order dispersion coefficient in μm⁴. Default is 0.

Returns:

n – Refractive index at given wavelength.

Return type:

float

Notes

Cauchy equation:

n(λ) = A + B/λ² + C/λ⁴

where λ is in micrometers.

Suitable for glasses and transparent materials in visible range. Less accurate than Sellmeier near absorption edges.

References

Examples

>>> n_550nm = cauchy_refractive_index(550e-9, A=1.5, B=0.004)
>>> print(f"n = {n_550nm:.4f}")

Factory Functions

create_sellmeier_material

lsurf.materials.create_sellmeier_material(name, B1, B2, B3, C1, C2, C3, absorption_coef=0.0)[source]

Create a material with Sellmeier dispersion.

Parameters:

name (str) – Material name.
B1 (float) – Sellmeier B coefficients (dimensionless).
B2 (float) – Sellmeier B coefficients (dimensionless).
B3 (float) – Sellmeier B coefficients (dimensionless).
C1 (float) – Sellmeier C coefficients in μm².
C2 (float) – Sellmeier C coefficients in μm².
C3 (float) – Sellmeier C coefficients in μm².
absorption_coef (float, optional) – Absorption coefficient in m⁻¹. Default is 0.

Returns:

material – Material with wavelength-dependent refractive index.

Return type:

HomogeneousMaterial

Examples

>>> # Create N-BK7 glass with Sellmeier dispersion
>>> bk7 = create_sellmeier_material(
...     "N-BK7",
...     B1=1.03961212, B2=0.231792344, B3=1.01046945,
...     C1=6.00069867e-3, C2=2.00179144e-2, C3=1.03560653e2,
... )

create_cauchy_material

lsurf.materials.create_cauchy_material(name, A, B, C=0.0, absorption_coef=0.0)[source]

Create a material with Cauchy dispersion.

Parameters:

name (str) – Material name.
A (float) – Constant term.
B (float) – First-order dispersion coefficient in μm².
C (float, optional) – Second-order dispersion coefficient in μm⁴. Default is 0.
absorption_coef (float, optional) – Absorption coefficient in m⁻¹. Default is 0.

Returns:

material – Material with wavelength-dependent refractive index.

Return type:

HomogeneousMaterial

Examples

>>> # Create simple glass with Cauchy dispersion
>>> glass = create_cauchy_material("Glass", A=1.52, B=0.004)

Examples

Using Pre-defined Materials

from lsurf import AIR_STP, WATER, BK7_GLASS

# Get refractive index at 532 nm
n_air = AIR_STP.refractive_index(532e-9, (0, 0, 0))
n_water = WATER.refractive_index(532e-9, (0, 0, 0))
n_glass = BK7_GLASS.refractive_index(532e-9, (0, 0, 0))

print(f"n_air = {n_air:.6f}")     # ~1.000293
print(f"n_water = {n_water:.6f}") # ~1.333
print(f"n_glass = {n_glass:.6f}") # ~1.519

Creating Custom Materials

Using Sellmeier formula:

from lsurf import create_sellmeier_material

# Fused silica
silica = create_sellmeier_material(
    name="Fused Silica",
    B1=0.6961663, C1=0.0684043,
    B2=0.4079426, C2=0.1162414,
    B3=0.8974794, C3=9.896161,
)

Using Cauchy formula:

from lsurf import create_cauchy_material

# Simple glass
glass = create_cauchy_material(
    name="Simple Glass",
    A=1.5,  # Base refractive index
    B=0.01, # Dispersion coefficient
)

Wavelength Dependence

import numpy as np
import matplotlib.pyplot as plt

wavelengths = np.linspace(400e-9, 700e-9, 100)
n_values = [WATER.refractive_index(wl, (0,0,0)) for wl in wavelengths]

plt.plot(wavelengths * 1e9, n_values)
plt.xlabel('Wavelength (nm)')
plt.ylabel('Refractive Index')
plt.title('Water Dispersion')
plt.show()