Propagation Module

The propagation module provides ray propagation engines for tracing rays through materials with varying refractive indices.

Overview

The propagation system includes:

  1. Propagators - Engines that advance rays through materials

  2. Kernels - GPU-accelerated computation kernels

  3. Factory - create_propagator() for automatic propagator selection

Propagators

GradientPropagator (CPU)

CPU-based propagator for gradient index materials using RK4 integration.

class lsurf.propagation.GradientPropagator(use_rk4=None, method=None, adaptive_stepping=True, min_step_size=1e-06, max_step_size=0.001, threads_per_block=256)[source]

Propagator for inhomogeneous media using ray equation integration.

This propagator handles materials with spatially-varying refractive indices by numerically integrating the ray equation (eikonal equation). Supports both Euler and RK4 integration methods.

use_rk4

If True, use RK4 integration; otherwise use Euler

Type:

bool

adaptive_stepping

If True, adapt step size to local gradient

Type:

bool

min_step_size

Minimum step size in meters

Type:

float

max_step_size

Maximum step size in meters

Type:

float

threads_per_block

CUDA threads per block (legacy, not used for CPU)

Type:

int

References

__init__(use_rk4=None, method=None, adaptive_stepping=True, min_step_size=1e-06, max_step_size=0.001, threads_per_block=256)[source]

Initialize gradient propagator.

Parameters:

  • use_rk4 (bool, optional) – Use RK4 instead of Euler (deprecated, use method instead)

  • method (str, optional) – Integration method: “euler” or “rk4”, default “rk4”

  • adaptive_stepping (bool, optional) – Enable adaptive step sizing, default True

  • min_step_size (float, optional) – Minimum step size in meters, default 1 μm

  • max_step_size (float, optional) – Maximum step size in meters, default 1 mm

  • threads_per_block (int, optional) – CUDA threads per block (legacy), default 256

compute_gradient_threshold(wavelength)[source]

Compute gradient threshold for switching to curved-path mode.

Parameters:

wavelength (float) – Wavelength in meters

Returns:

Gradient magnitude threshold in m^-1

Return type:

float

Notes

Threshold is chosen such that ray curvature becomes significant over a distance of ~wavelength: |∇n|/n > 1/λ

propagate_step_cpu(positions, directions, active, material, step_size, wavelength, geometric_path_lengths=None, optical_path_lengths=None, accumulated_time=None)[source]

Propagate rays through gradient medium using CPU (vectorized NumPy).

This method modifies arrays in-place using Euler or RK4 integration.

Parameters:

  • positions (ndarray) – Ray positions, shape (N, 3), modified in-place.

  • directions (ndarray) – Ray directions, shape (N, 3), modified in-place.

  • active (ndarray) – Active ray mask, shape (N,).

  • material (MaterialFieldProtocol) – Material with get_refractive_index and get_refractive_index_gradient.

  • step_size (float) – Step size in meters.

  • wavelength (float) – Wavelength in meters.

  • geometric_path_lengths (ndarray, optional) – Geometric path lengths, shape (N,), updated in-place.

  • optical_path_lengths (ndarray, optional) – Optical path lengths, shape (N,), updated in-place.

  • accumulated_time (ndarray, optional) – Accumulated time, shape (N,), updated in-place.

propagate_step(rays, step_size, wavelength=5.32e-07, material=None)[source]

Propagate rays by a single integration step.

This is the unified interface method compatible with GPU propagators. Delegates to propagate_step_cpu internally.

Parameters:

  • rays (RayBatch) – Ray batch containing positions, directions, and other ray properties. Modified in-place.

  • step_size (float) – Integration step size in meters.

  • wavelength (float, optional) – Wavelength in meters, default 532 nm.

  • material (MaterialFieldProtocol, optional) – Material to propagate through. If not provided, must be set via a separate mechanism (e.g., stored during initialization).

propagate(rays, total_distance, step_size, wavelength=5.32e-07, material=None)[source]

Propagate rays through a total distance.

This is the unified interface method compatible with GPU propagators.

Parameters:

  • rays (RayBatch) – Ray batch to propagate (modified in-place)

  • total_distance (float) – Total distance to propagate in meters

  • step_size (float) – Integration step size in meters

  • wavelength (float, optional) – Wavelength in meters, default 532 nm

  • material (MaterialFieldProtocol, optional) – Material to propagate through

GPUGradientPropagator (GPU)

GPU-accelerated propagator for gradient index materials.

class lsurf.propagation.GPUGradientPropagator(material, method='rk4', threads_per_block=256, enable_absorption=False)[source]

High-performance GPU-based propagator for inhomogeneous materials.

This propagator keeps all computation on the GPU for maximum performance. It works with any material that implements GPUMaterialProtocol, using the material’s gpu_material_id to dispatch to the appropriate GPU kernels.

Currently Supported Materials

  • ExponentialAtmosphere (GPUMaterialID.EXPONENTIAL_ATMOSPHERE)

For materials without GPU support, use GradientPropagator.propagate_step_cpu() which works with any MaterialField via the CPU.

param material:

The material to propagate through. Must implement GPUMaterialProtocol.

type material:

GPUMaterialProtocol

param method:

Integration method: ‘euler’ or ‘rk4’. Default: ‘rk4’

type method:

str

param threads_per_block:

CUDA threads per block. Default: 256

type threads_per_block:

int

Example

>>> from lsurf.materials import ExponentialAtmosphere
>>> atmosphere = ExponentialAtmosphere()
>>> propagator = GPUGradientPropagator(atmosphere, method='rk4')
>>> propagator.propagate(rays, total_distance=500e3, step_size=100.0)
raises ValueError:

If the material does not support GPU acceleration.

__init__(material, method='rk4', threads_per_block=256, enable_absorption=False)[source]
property material: GPUMaterialProtocol

The material being propagated through.

propagate_step(rays, step_size, wavelength=5.32e-07)[source]

Propagate rays by a single integration step.

This is the unified interface method compatible with CPU propagators. Performs a single step propagation on GPU.

Parameters:

  • rays (RayBatch) – Ray batch containing positions, directions, and other ray properties. Modified in-place.

  • step_size (float) – Integration step size in meters.

  • wavelength (float, optional) – Wavelength in meters, default 532 nm.

propagate(rays, total_distance, step_size, wavelength=5.32e-07, steps_per_kernel=100, target_dn_frac=None)[source]

Propagate rays through material on GPU.

Modifies the rays object in-place.

Parameters:

  • rays (RayBatch) – Ray batch to propagate (modified in-place)

  • total_distance (float) – Total distance to propagate in meters

  • step_size (float) – Integration step size in meters (max_step in adaptive mode)

  • wavelength (float, optional) – Wavelength in meters, default 532 nm (green light)

  • steps_per_kernel (int) – Number of steps per kernel launch. Higher values reduce kernel launch overhead but increase register pressure.

  • target_dn_frac (float or None) – If provided, use adaptive stepping with this fractional refractivity change target. The step_size becomes max_step and min_step_size from the material is used as min_step.

property last_adaptive_step_stats

Return adaptive step stats from the most recent propagate() call.

propagate_with_history(rays, total_distance, step_size, history_interval=100, wavelength=1e-06, target_dn_frac=None, history_distance_interval=None)[source]

Propagate rays and record position history at intervals.

Parameters:

  • rays (RayBatch) – Ray batch to propagate (modified in-place)

  • total_distance (float) – Total distance to propagate in meters

  • step_size (float) – Integration step size in meters (max_step in adaptive mode)

  • history_interval (int) – Record history every N steps (fixed-step mode only)

  • wavelength (float) – Wavelength in meters

  • target_dn_frac (float or None) – If provided, use adaptive stepping with this fractional refractivity change target.

  • history_distance_interval (float or None) – Distance between history records in adaptive mode. If None in adaptive mode, defaults to step_size * history_interval.

Returns:

Dictionary with ‘positions’, ‘directions’, ‘distances’ arrays

Return type:

dict

get_refractive_index(positions, wavelength=1e-06)[source]

Compute refractive index at positions (CPU, delegates to material).

Parameters:

  • positions (ndarray (N, 3)) – Position coordinates

  • wavelength (float) – Wavelength in meters

Returns:

Refractive index at each position

Return type:

ndarray (N,)

Propagator Factory

lsurf.propagation.create_propagator(material, propagator_id=None, kernel_id=None, method=None, prefer_gpu=True, **kwargs)[source]

Create the appropriate propagator for a material.

Selection priority (highest to lowest): 1. If propagator_id/kernel_id passed here → use those (override) 2. Else use material’s stored preference (material._propagator_id, material._kernel_id) 3. Else use material’s class defaults

Parameters:

  • material (MaterialField) – The material to propagate through.

  • propagator_id (PropagatorID, optional) – Override propagator selection. If None, uses material’s preference.

  • kernel_id (PropagationKernelID, optional) – Override kernel selection. If None, uses material’s preference. Note: Only used for GPU propagators.

  • method (str, optional) – Integration method (“euler” or “rk4”). If None, determined from kernel_id.

  • prefer_gpu (bool, default True) – If True and material supports GPU, use GPU propagator. If False, force CPU propagator.

  • **kwargs (Any) – Additional propagator configuration (e.g., threads_per_block).

Returns:

Configured propagator instance appropriate for the material.

Return type:

Propagator

Raises:

ValueError – If propagator_id or kernel_id is not supported by the material.

Examples

>>> from lsurf.materials import ExponentialAtmosphere
>>> from lsurf.propagation.propagators import create_propagator
>>> from lsurf.propagation.kernels import PropagatorID, PropagationKernelID

>>> # Use defaults (RK4 kernel, GPU propagator)
>>> atmo = ExponentialAtmosphere(n_sea_level=1.000293)
>>> propagator = create_propagator(atmo)

>>> # Override kernel at propagator creation time
>>> propagator_euler = create_propagator(
...     atmo,
...     kernel_id=PropagationKernelID.SIMPLE_EULER
... )

>>> # Force CPU propagator
>>> cpu_propagator = create_propagator(atmo, propagator_id=PropagatorID.CPU_GRADIENT)

>>> # Or use prefer_gpu=False
>>> cpu_propagator = create_propagator(atmo, prefer_gpu=False)

Usage Example

from lsurf.propagation import create_propagator, GPUMaterialID
from lsurf.materials import ExponentialAtmosphere

# Create material
atmosphere = ExponentialAtmosphere(n_sea_level=1.000293)

# Create propagator (auto-selects GPU if available)
propagator = create_propagator(atmosphere)

# Or explicitly request CPU
from lsurf.propagation.kernels.registry import PropagatorID
cpu_propagator = create_propagator(
    atmosphere,
    propagator_id=PropagatorID.CPU_GRADIENT
)

Material IDs

class lsurf.propagation.GPUMaterialID(*values)[source]

Registry of GPU-compatible material types.

Each material type has specialized device functions compiled into the CUDA kernels. This enum identifies which material’s functions to use.

Tiers: - EXPONENTIAL_ATMOSPHERE (1): Legacy specialized kernels - SIMPLE_INHOMOGENEOUS (2): 1D radial profile via LUT interpolation - GRID_INHOMOGENEOUS (3): 3D grid data via trilinear interpolation - SPECTRAL_INHOMOGENEOUS (4): 2D LUT-based (altitude × wavelength)

EXPONENTIAL_ATMOSPHERE = 1
SIMPLE_INHOMOGENEOUS = 2
GRID_INHOMOGENEOUS = 3
SPECTRAL_INHOMOGENEOUS = 4
DUCT_ATMOSPHERE = 5
US_DUCT_ATMOSPHERE = 6

Kernel Registry

The kernel registry manages available propagation kernels:

Kernel and Propagator Registry

Central registry for all kernel types in the ray tracing system: - PropagationKernelID: Material propagation kernels - IntersectionKernelID: Surface intersection kernels - DetectionKernelID: Ray detection kernels - FresnelKernelID: Fresnel reflection/refraction kernels - PropagatorID: Propagator implementations

This module provides: - Kernel ID enums for each category - Registration decorator and lookup functions - Unified registry for all kernel types

class lsurf.propagation.kernels.registry.PropagationKernelID(*values)[source]

Identifiers for material propagation kernels.

Each kernel corresponds to a specific integration method and material type. Materials declare which kernels they support, and propagators use these IDs to select the appropriate CUDA kernel.

Naming convention: {MATERIAL_TYPE}_{INTEGRATION_METHOD}

Examples

>>> from lsurf.propagation.kernels import PropagationKernelID
>>> PropagationKernelID.SIMPLE_RK4
<PropagationKernelID.SIMPLE_RK4: ...>
SIMPLE_EULER = 1
SIMPLE_RK4 = 2
SPECTRAL_EULER = 3
SPECTRAL_RK4 = 4
SPECTRAL_EULER_PERRAY = 5
SPECTRAL_RK4_PERRAY = 6
GRID_EULER = 7
GRID_RK4 = 8
DUCT_EULER = 9
DUCT_RK4 = 10
DUCT_EULER_ADAPTIVE = 11
DUCT_RK4_ADAPTIVE = 12
US_DUCT_EULER = 13
US_DUCT_RK4 = 14
US_DUCT_EULER_ADAPTIVE = 15
US_DUCT_RK4_ADAPTIVE = 16
class lsurf.propagation.kernels.registry.IntersectionKernelID(*values)[source]

Identifiers for surface intersection kernels.

Each kernel corresponds to a specific geometry type and intersection method. Surfaces declare which kernels they support for ray-surface intersection.

Note: Complex surfaces (waves, curved surfaces) use the generic signed distance kernels which dispatch to device functions based on geometry_id.

Naming convention: {GEOMETRY_TYPE}_{METHOD}

Examples

>>> from lsurf.propagation.kernels import IntersectionKernelID
>>> IntersectionKernelID.PLANE_ANALYTICAL
<IntersectionKernelID.PLANE_ANALYTICAL: ...>
PLANE_ANALYTICAL = 1
SPHERE_ANALYTICAL = 2
BOUNDED_PLANE_ANALYTICAL = 3
ANNULAR_PLANE_ANALYTICAL = 4
SIGNED_DISTANCE_GENERIC = 5
SURFACE_NORMAL_GENERIC = 6
class lsurf.propagation.kernels.registry.DetectionKernelID(*values)[source]

Identifiers for ray detection kernels.

Each kernel corresponds to a specific detector geometry.

Naming convention: {GEOMETRY_TYPE}_{VARIANT}

Examples

>>> from lsurf.propagation.kernels import DetectionKernelID
>>> DetectionKernelID.SPHERICAL_SINGLE
<DetectionKernelID.SPHERICAL_SINGLE: ...>
SPHERICAL_SINGLE = 1
SPHERICAL_MULTI = 2
PLANAR_SINGLE = 3
class lsurf.propagation.kernels.registry.FresnelKernelID(*values)[source]

Identifiers for Fresnel reflection/refraction kernels.

Each kernel implements a specific approach to computing Fresnel coefficients.

Examples

>>> from lsurf.propagation.kernels import FresnelKernelID
>>> FresnelKernelID.STANDARD
<FresnelKernelID.STANDARD: ...>
STANDARD = 1
POLARIZED = 2
class lsurf.propagation.kernels.registry.PropagatorID(*values)[source]

Identifiers for available propagator implementations.

Each propagator ID corresponds to a specific propagator class that handles ray propagation through materials.

Examples

>>> from lsurf.propagation.kernels import PropagatorID
>>> PropagatorID.GPU_GRADIENT
<PropagatorID.GPU_GRADIENT: ...>
CPU_GRADIENT = 1
GPU_GRADIENT = 2
GPU_SPECTRAL = 3
GPU_SURFACE = 4
lsurf.propagation.kernels.registry.register_kernel(kernel_id)[source]

Decorator to register a kernel function with the global registry.

Parameters:

kernel_id (KernelID) – The unique identifier for this kernel. Can be any kernel ID type: PropagationKernelID, IntersectionKernelID, DetectionKernelID, or FresnelKernelID.

Returns:

A decorator that registers the function and returns it unchanged.

Return type:

decorator

Examples

>>> @register_kernel(PropagationKernelID.SIMPLE_EULER)
... def my_euler_kernel(positions, directions, ...):
...     ...

>>> @register_kernel(DetectionKernelID.SPHERICAL_SINGLE)
... def my_detection_kernel(positions, directions, ...):
...     ...

Notes

Registration happens at import time. The registered function can be retrieved later using get_kernel().

lsurf.propagation.kernels.registry.get_kernel(kernel_id)[source]

Retrieve a registered kernel function by ID.

Parameters:

kernel_id (KernelID) – The identifier of the kernel to retrieve. Can be any kernel ID type.

Returns:

The registered kernel function.

Return type:

Callable

Raises:

KeyError – If the kernel ID is not registered.

Examples

>>> kernel = get_kernel(PropagationKernelID.SIMPLE_RK4)
>>> kernel(positions, directions, ...)

>>> detect = get_kernel(DetectionKernelID.SPHERICAL_SINGLE)
>>> detect(positions, directions, ...)
lsurf.propagation.kernels.registry.get_registered_kernels()[source]

Return a copy of all registered kernels.

Returns:

Dictionary mapping kernel IDs to kernel functions.

Return type:

dict

lsurf.propagation.kernels.registry.list_kernels(kernel_type=None)[source]

List all registered kernels, optionally filtered by type.

Parameters:

kernel_type (type[Enum], optional) – If provided, only return kernels of this enum type. For example, PropagationKernelID to list only propagation kernels.

Returns:

List of registered kernel IDs.

Return type:

list[KernelID]

Examples

>>> # List all registered kernels
>>> all_kernels = list_kernels()

>>> # List only propagation kernels
>>> prop_kernels = list_kernels(PropagationKernelID)

>>> # List only detection kernels
>>> detect_kernels = list_kernels(DetectionKernelID)

Available Kernels

Propagation Kernels

  • SIMPLE_EULER - First-order Euler integration

  • SIMPLE_RK4 - Fourth-order Runge-Kutta integration

  • SPECTRAL_EULER - Wavelength-dependent Euler

  • SPECTRAL_RK4 - Wavelength-dependent RK4

  • GRID_EULER - Grid-based refractive index lookup (Euler)

  • GRID_RK4 - Grid-based refractive index lookup (RK4)

Intersection Kernels

  • PLANE_ANALYTICAL - Analytical plane intersection

  • SPHERE_ANALYTICAL - Analytical sphere intersection

  • SDF_BISECTION - Signed distance function with bisection refinement

Architecture

The propagation system follows a layered architecture:

  1. Simulation calls MaterialPropagator

  2. MaterialPropagator uses registered Kernels for:

    • Propagation (advancing rays)

    • Intersection (finding surface hits)

    • Detection (recording hits)

  3. Kernels are GPU-accelerated via Numba CUDA

Device Functions

Low-level device functions used by GPU kernels:

Shared GPU Device Functions

This module consolidates device functions used across GPU kernels: - Integration functions (Euler, RK4 steps) - Dispersion models (Sellmeier, Cauchy equations) - Common utilities (normalization, interpolation)

All functions are designed to be used with @cuda.jit(device=True).

lsurf.propagation.kernels.device_functions.device_adaptive_step_size(gradient_magnitude, refractive_index, wavelength, min_step, max_step)[source]

Compute adaptive step size based on local gradient on GPU device.

Parameters:

  • gradient_magnitude (float) – |∇n| at current position

  • refractive_index (float) – n at current position

  • wavelength (float) – Wavelength in meters

  • min_step (float) – Minimum allowed step size

  • max_step (float) – Maximum allowed step size

Returns:

Recommended step size

Return type:

float

lsurf.propagation.kernels.device_functions.device_euler_step(x, y, z, dx, dy, dz, n, grad_x, grad_y, grad_z, step_size)[source]

Single Euler integration step for ray equation on GPU.

The ray equation in gradient media is:

dr/ds = d̂ dd̂/ds = (∇n - (d̂·∇n)d̂) / n

Parameters:

  • x (float) – Current position

  • y (float) – Current position

  • z (float) – Current position

  • dx (float) – Current direction (unit vector)

  • dy (float) – Current direction (unit vector)

  • dz (float) – Current direction (unit vector)

  • n (float) – Refractive index at current position

  • grad_x (float) – Gradient of n at current position

  • grad_y (float) – Gradient of n at current position

  • grad_z (float) – Gradient of n at current position

  • step_size (float) – Step size ds

Returns:

(new_x, new_y, new_z, new_dx, new_dy, new_dz)

Return type:

tuple

lsurf.propagation.kernels.device_functions.device_rk4_step(x, y, z, dx, dy, dz, step_size, n_func, grad_func, *material_params)[source]

RK4 integration step for ray equation on GPU.

This is a generic RK4 step that takes material evaluation functions as parameters. For material-specific implementations with better performance, materials should define their own device functions.

Parameters:

  • x (float) – Current position

  • y (float) – Current position

  • z (float) – Current position

  • dx (float) – Current direction (unit vector)

  • dy (float) – Current direction (unit vector)

  • dz (float) – Current direction (unit vector)

  • step_size (float) – Integration step size

  • n_func (device function) – Function to evaluate n(x, y, z, *material_params)

  • grad_func (device function) – Function to evaluate ∇n(x, y, z, *material_params) -> (gx, gy, gz)

  • material_params (tuple) – Material-specific parameters passed to n_func and grad_func

Returns:

(new_x, new_y, new_z, new_dx, new_dy, new_dz, n_avg)

Return type:

tuple of float

Notes

For best performance, materials should implement their own specialized RK4 device functions that inline the n and gradient calculations. This generic version has function call overhead.

lsurf.propagation.kernels.device_functions.device_sellmeier_equation(wl_um, B1, B2, B3, C1, C2, C3)[source]

GPU-compatible Sellmeier equation.

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

lsurf.propagation.kernels.device_functions.device_cauchy_equation(wl_um, A, B, C)[source]

GPU-compatible Cauchy equation.

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

lsurf.propagation.kernels.device_functions.euler_step(x, y, z, dx, dy, dz, n, grad_x, grad_y, grad_z, step_size)[source]

Single Euler integration step for ray equation.

Parameters:

  • x (float) – Current position

  • y (float) – Current position

  • z (float) – Current position

  • dx (float) – Current direction (unit vector)

  • dy (float) – Current direction (unit vector)

  • dz (float) – Current direction (unit vector)

  • n (float) – Refractive index at current position

  • grad_x (float) – Gradient of n at current position

  • grad_y (float) – Gradient of n at current position

  • grad_z (float) – Gradient of n at current position

  • step_size (float) – Step size ds

Returns:

(new_x, new_y, new_z, new_dx, new_dy, new_dz)

Return type:

tuple

lsurf.propagation.kernels.device_functions.rk4_step(x, y, z, dx, dy, dz, n_func, grad_func, wavelength, step_size)[source]

RK4 integration step for ray equation in gradient medium.

Solves the coupled ODEs:

dr/ds = d̂ dd̂/ds = (∇n - (d̂·∇n)d̂) / n

Parameters:

  • x (float) – Current position

  • y (float) – Current position

  • z (float) – Current position

  • dx (float) – Current direction (unit vector)

  • dy (float) – Current direction (unit vector)

  • dz (float) – Current direction (unit vector)

  • n_func (callable) – Function to evaluate n(x, y, z, wavelength)

  • grad_func (callable) – Function to evaluate ∇n(x, y, z, wavelength) -> (gx, gy, gz)

  • wavelength (float) – Wavelength

  • step_size (float) – Integration step size

Returns:

(new_x, new_y, new_z, new_dx, new_dy, new_dz, optical_path_increment)

Return type:

tuple of float

lsurf.propagation.kernels.device_functions.compute_adaptive_step_size(gradient_magnitude, refractive_index, wavelength, min_step, max_step)[source]

Compute adaptive step size based on local gradient.

Parameters:

  • gradient_magnitude (float) – |∇n| at current position

  • refractive_index (float) – n at current position

  • wavelength (float) – Wavelength in meters

  • min_step (float) – Minimum allowed step size

  • max_step (float) – Maximum allowed step size

Returns:

Recommended step size

Return type:

float

Notes

Step size is chosen to resolve curvature: Δs ≤ R_c / 10 where R_c = n / |∇n| is the radius of curvature.

lsurf.propagation.kernels.device_functions.normalize_directions(directions)[source]

Normalize direction vectors to unit length.

lsurf.propagation.kernels.device_functions.euler_step_batch(positions, directions, active_mask, n, grad_n, step_size)[source]

Vectorized Euler step for batch of rays.

Parameters:

  • positions (ndarray of shape (N, 3)) – Ray positions

  • directions (ndarray of shape (N, 3)) – Ray directions (unit vectors)

  • active_mask (ndarray of shape (N,)) – Boolean mask for active rays

  • n (ndarray of shape (N,)) – Refractive index at each position

  • grad_n (ndarray of shape (N, 3)) – Gradient of n at each position

  • step_size (float) – Step size

Returns:

  • new_positions (ndarray of shape (N, 3))

  • new_directions (ndarray of shape (N, 3))

Return type:

tuple[ndarray[tuple[Any, …], dtype[float32]], ndarray[tuple[Any, …], dtype[float32]]]

lsurf.propagation.kernels.device_functions.rk4_step_batch(positions, directions, active_mask, material, step_size, wavelength)[source]

Vectorized RK4 step for batch of rays.

Parameters:

  • positions (ndarray of shape (N, 3)) – Ray positions

  • directions (ndarray of shape (N, 3)) – Ray directions (unit vectors)

  • active_mask (ndarray of shape (N,)) – Boolean mask for active rays

  • material (MaterialFieldProtocol) – Material providing n and ∇n evaluation

  • step_size (float) – Step size

  • wavelength (float) – Wavelength

Returns:

  • new_positions (ndarray of shape (N, 3))

  • new_directions (ndarray of shape (N, 3))

  • n_avg (ndarray of shape (N,)) – Average refractive index over step (for optical path)

Return type:

tuple[ndarray[tuple[Any, …], dtype[float32]], ndarray[tuple[Any, …], dtype[float32]], ndarray[tuple[Any, …], dtype[float32]]]

Integration with Simulation

The Simulation class automatically creates and manages propagators:

from lsurf.simulation import Simulation, SimulationConfig

# Simulation creates propagator internally
config = SimulationConfig(use_gpu=True)  # Use GPU if available
sim = Simulation(geometry, config)

# Or force CPU-only
config = SimulationConfig(use_gpu=False)
sim = Simulation(geometry, config)