ler.image_properties.epl_shear_njit

Numba-JIT-compiled EPL + external shear lensing solver.

Provides analytical image-position finding, magnification computation, Fermat-potential evaluation, and a parallel batch solver for elliptical power-law (EPL) lens models with external shear.

The module uses the 1-D lens-equation approach, reducing the 2-D vector equation to a scalar root-finding problem parameterised by the image-plane angle. All computationally intensive functions are compiled with @njit (Numba) for performance.

The top-level entry points are:

Copyright (C) 2026 Author Name. Distributed under MIT License.

Module Contents

Functions

pol_to_ell(r, theta, q)

Convert polar coordinates to elliptical coordinates.

omega_scalar(phi, t, q[, niter_max, tol])

Scalar series expansion of the EPL deflection kernel Omega.

lensing_diagnostics_scalar(z, b, t, gamma1, gamma2, q, ...)

Compute magnification and image type at one image position.

fermat_potential_scalar(z, x, y, x_source, y_source, ...)

Compute the Fermat potential at one image position.

image_position_analytical_njit(x_src, y_src, q, phi, ...)

Standalone EPL + external shear analytical image finder.

create_epl_shear_solver([arrival_time_sort, max_img, ...])

Create a parallel EPL + shear solver for batched lens systems.

Attributes

EPS

MAX_ROOTS

MAX_IMGS

C_LIGHT

ler.image_properties.epl_shear_njit.EPS = '1e-16'[source]
ler.image_properties.epl_shear_njit.MAX_ROOTS = '16'[source]
ler.image_properties.epl_shear_njit.MAX_IMGS = '5'[source]
ler.image_properties.epl_shear_njit.C_LIGHT = '299792458.0'[source]
ler.image_properties.epl_shear_njit.pol_to_ell(r, theta, q)[source]

Convert polar coordinates to elliptical coordinates.

Parameters:
rfloat

Polar radial coordinate.

thetafloat

Polar angle in radians.

qfloat

Axis ratio of the ellipse (0 < q ≤ 1).

Returns:
rellfloat

Elliptical radial coordinate, rell = r * sqrt(q^2*cos^2(theta) + sin^2(theta)).

phifloat

Elliptical angle in radians, phi = arctan2(sin(theta), cos(theta)*q).

Examples

>>> import numpy as np
>>> from ler.image_properties.epl_shear_njit import pol_to_ell
>>> rell, phi = pol_to_ell(1.0, np.pi / 4, 0.7)
ler.image_properties.epl_shear_njit.omega_scalar(phi, t, q, niter_max=200, tol=1e-16)[source]

Scalar series expansion of the EPL deflection kernel Omega.

Evaluates the convergent series for the EPL Omega function at a single elliptical angle, avoiding temporary array allocation. The series converges geometrically with ratio f = (1-q)/(1+q).

Parameters:
phifloat

Elliptical angle in radians.

tfloat

EPL slope parameter (t = gamma - 1).

qfloat

Axis ratio (0 < q ≤ 1).

niter_maxint

Maximum number of series iterations. default: 200

tolfloat

Convergence tolerance for series truncation. default: 1e-16

Returns:
omega_sumcomplex

Complex EPL deflection kernel value at angle phi.

Examples

>>> import numpy as np
>>> from ler.image_properties.epl_shear_njit import omega_scalar
>>> omega = omega_scalar(np.pi / 3, t=1.0, q=0.8)
ler.image_properties.epl_shear_njit.lensing_diagnostics_scalar(z, b, t, gamma1, gamma2, q, phi, Omega)[source]

Compute magnification and image type at one image position.

Evaluates the total Hessian of the lensing potential (EPL + external shear), then derives the Jacobian determinant and trace to classify the image and compute its signed magnification. Inputs z, b, t, q, phi, and Omega are all scalars, not arrays.

Image type convention: - 1: Type I (minimum of the Fermat potential) - 2: Type II (saddle point) - 3: Type III (maximum of the Fermat potential) - 0: undefined / degenerate (on a critical curve)

Parameters:
zcomplex

Image position in the axis-aligned frame (z = exp(-i*phi) * (x + i*y)).

bfloat

EPL scale parameter (b = theta_E * sqrt(q)).

tfloat

EPL slope parameter (t = gamma - 1).

gamma1float

External shear component 1.

gamma2float

External shear component 2.

qfloat

Axis ratio (0 < q ≤ 1).

phifloat

Lens position angle in radians.

Omegacomplex

EPL deflection kernel at the image elliptical angle.

Returns:
mufloat

Signed magnification 1/det(A); np.inf on a critical curve.

image_typeint

Image-type code (0, 1, 2, or 3).

Examples

>>> import numpy as np
>>> from ler.image_properties.epl_shear_njit import (
...     lensing_diagnostics_scalar, omega_scalar
... )
>>> z = np.exp(-1j * 0.0) * (0.8 + 1j * 0.3)
>>> phi_ell = np.angle(z.real * 0.8 + 1j * z.imag)
>>> Omega = omega_scalar(phi_ell, t=1.0, q=0.8)
>>> mu, itype = lensing_diagnostics_scalar(
...     z, b=0.894, t=1.0, gamma1=0.0, gamma2=0.0, q=0.8, phi=0.0, Omega=Omega
... )
ler.image_properties.epl_shear_njit.fermat_potential_scalar(z, x, y, x_source, y_source, b, t, gamma1, gamma2, q, phi, Omega)[source]

Compute the Fermat potential at one image position.

Returns the geometric minus gravitational time-delay contribution: tau = 0.5*|theta - beta|^2 - psi_EPL(theta) - psi_shear(theta)

Parameters:
zcomplex

Image position in the axis-aligned frame (z = exp(-i*phi) * (x + i*y)).

xfloat

Image x-coordinate in sky frame.

yfloat

Image y-coordinate in sky frame.

x_sourcefloat

Source x-coordinate in sky frame.

y_sourcefloat

Source y-coordinate in sky frame.

bfloat

EPL scale parameter (b = theta_E * sqrt(q)).

tfloat

EPL slope parameter (t = gamma - 1).

gamma1float

External shear component 1.

gamma2float

External shear component 2.

qfloat

Axis ratio (0 < q ≤ 1).

phifloat

Lens position angle in radians.

Omegacomplex

EPL deflection kernel at the image elliptical angle.

Returns:
taufloat

Fermat potential (dimensionless; proportional to arrival-time delay).

ler.image_properties.epl_shear_njit.image_position_analytical_njit(x_src, y_src, q, phi, gamma, gamma1, gamma2, theta_E=1.0, alpha_scaling=1.0, magnification_limit=0.01, Nmeas=400, Nmeas_extra=80)[source]

Standalone EPL + external shear analytical image finder.

Locates all lensed images for a given source position, computes signed magnifications, Fermat potentials (arrival-time proxies), and image types. Results are sorted by ascending arrival time and filtered by a minimum magnification threshold.

Parameters:
x_srcfloat

Source x-coordinate (normalised to Einstein radius).

y_srcfloat

Source y-coordinate (normalised to Einstein radius).

qfloat

Lens axis ratio (0 < q ≤ 1).

phifloat

Lens position angle in radians.

gammafloat

EPL power-law slope (gamma = 2 for isothermal).

gamma1float

External shear component 1.

gamma2float

External shear component 2.

theta_Efloat

Einstein radius. default: 1.0

alpha_scalingfloat

Deflection scaling applied as theta_E_eff = theta_E * alpha_scaling^(1/(gamma-1)). default: 1.0

magnification_limitfloat

Minimum |mu| for an image to be retained. default: 0.01

Nmeasint

Number of uniformly-spaced angle samples for root finding. default: 400

Nmeas_extraint

Number of extra refinement samples near the source angle. default: 80

Returns:
x_imgnumpy.ndarray

Image x-coordinates.

y_imgnumpy.ndarray

Image y-coordinates.

fermat_potnumpy.ndarray

Fermat potentials at each image.

magnificationnumpy.ndarray

Signed magnifications at each image.

image_typenumpy.ndarray

Image-type codes (1 = min, 2 = saddle, 3 = max, 0 = degenerate).

nimgint

Number of valid images retained.

Examples

>>> from ler.image_properties.epl_shear_njit import image_position_analytical_njit
>>> x_img, y_img, fermat_pot, mu, itype, nimg = image_position_analytical_njit(
...     x_src=0.1, y_src=0.05,
...     q=0.8, phi=0.3, gamma=2.0,
...     gamma1=0.05, gamma2=0.02,
... )
>>> print(nimg)
ler.image_properties.epl_shear_njit.create_epl_shear_solver(arrival_time_sort=True, max_img=4, num_th=500, maginf=-100.0, max_source_sample_attempts=128, alpha_scaling=1.0, magnification_limit=0.01, Nmeas=400, Nmeas_extra=80)[source]

Create a parallel EPL + shear solver for batched lens systems.

Returns a JIT-compiled function that, for each system in a batch, samples a source from the double caustic and solves for image positions, magnifications, time delays, and image types.

Parameters:
arrival_time_sortbool

Whether to sort images by arrival time.

default: True

max_imgint

Maximum number of images to store per system.

default: 4

num_thint

Angular samples for caustic construction.

default: 500

maginffloat

Magnification cutoff for caustic boundary.

default: -100.0

max_source_sample_attemptsint

Maximum draws from sample_source_from_double_caustic per system when the returned (beta_x, beta_y) are non-finite (invalid caustic / geometry). Stops early on the first finite pair.

default: 128

alpha_scalingfloat

Deflection scaling factor.

default: 1.0

magnification_limitfloat

Minimum abs(mu) threshold for image retention.

default: 0.01

Nmeasint

Angular root-finding grid size.

default: 400

Nmeas_extraint

Extra refinement points.

default: 80

Returns:
solve_epl_shear_multithreadedcallable

Parallel solver function with signature

(theta_E, D_dt, q, phi, gamma, gamma1, gamma2)

returning a tuple of result arrays.

Examples

>>> solver = create_epl_shear_solver()
>>> results = solver(theta_E, D_dt, q, phi, gamma, gamma1, gamma2)